songtianlun/openfoam
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Merge branch 'master' into molecularDynamics

Browse Source
This commit is contained in:
graham 2009-06-29 13:40:29 +01:00
commit c73fa61a97
762 changed files with 16609 additions and 10749 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

4
applications/solvers/combustion/PDRFoam/Make/options
View File

@ -8,7 +8,7 @@ EXE_INC = \
-I$(LIB_SRC)/sampling/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/combustion/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
-I$(LIB_SRC)/turbulenceModels \
-I$(LIB_SRC)/turbulenceModels/compressible/RAS/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/laminarFlameSpeed/lnInclude \
@ -23,7 +23,7 @@ EXE_LIBS = \
-lmeshTools \
-lcompressibleRASModels \
-lbasicThermophysicalModels \
-lcombustionThermophysicalModels \
-lreactionThermophysicalModels \
-lspecie \
-llaminarFlameSpeedModels \
-lfiniteVolume \

50
applications/solvers/combustion/PDRFoam/PDRFoam.C
View File

@ -36,7 +36,7 @@ Description
to be appropriate by comparison with the results from the
spectral model.
Strain effects are encorporated directly into the Xi equation
Strain effects are incorporated directly into the Xi equation
but not in the algebraic approximation. Further work need to be
done on this issue, particularly regarding the enhanced removal rate
caused by flame compression. Analysis using results of the spectral
@ -70,53 +70,53 @@ Description
int main(int argc, char *argv[])
{
# include "setRootCase.H"
#include "setRootCase.H"
# include "createTime.H"
# include "createMesh.H"
# include "readCombustionProperties.H"
# include "readEnvironmentalProperties.H"
# include "createFields.H"
# include "readPISOControls.H"
# include "initContinuityErrs.H"
# include "readTimeControls.H"
# include "CourantNo.H"
# include "setInitialDeltaT.H"
#include "createTime.H"
#include "createMesh.H"
#include "readCombustionProperties.H"
#include "readEnvironmentalProperties.H"
#include "createFields.H"
#include "readPISOControls.H"
#include "initContinuityErrs.H"
#include "readTimeControls.H"
#include "CourantNo.H"
#include "setInitialDeltaT.H"
scalar StCoNum = 0.0;
scalar StCoNum = 0.0;
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
while (runTime.run())
{
# include "readTimeControls.H"
# include "readPISOControls.H"
# include "CourantNo.H"
# include "setDeltaT.H"
#include "readTimeControls.H"
#include "readPISOControls.H"
#include "CourantNo.H"
#include "setDeltaT.H"
runTime++;
Info<< "\n\nTime = " << runTime.timeName() << endl;
# include "rhoEqn.H"
# include "UEqn.H"
#include "rhoEqn.H"
#include "UEqn.H"
// --- PISO loop
for (int corr=1; corr<=nCorr; corr++)
{
# include "bEqn.H"
# include "ftEqn.H"
# include "huEqn.H"
# include "hEqn.H"
#include "bEqn.H"
#include "ftEqn.H"
#include "huEqn.H"
#include "hEqn.H"
if (!ign.ignited())
{
hu == h;
}
# include "pEqn.H"
#include "pEqn.H"
}
turbulence->correct();

28
applications/solvers/combustion/PDRFoam/StCourantNo.H
View File

@ -31,23 +31,25 @@ Description
\*---------------------------------------------------------------------------*/
{
scalar meanStCoNum = 0.0;
scalar meanStCoNum = 0.0;
if (mesh.nInternalFaces())
{
surfaceScalarField SfUfbyDelta =
mesh.surfaceInterpolation::deltaCoeffs()
*mag(phiSt/fvc::interpolate(rho));
if (mesh.nInternalFaces())
{
surfaceScalarField SfUfbyDelta =
mesh.surfaceInterpolation::deltaCoeffs()
*mag(phiSt/fvc::interpolate(rho));
StCoNum = max(SfUfbyDelta/mesh.magSf())
.value()*runTime.deltaT().value();
StCoNum =
max(SfUfbyDelta/mesh.magSf()).value()
*runTime.deltaT().value();
meanStCoNum = (sum(SfUfbyDelta)/sum(mesh.magSf()))
.value()*runTime.deltaT().value();
}
meanStCoNum =
(sum(SfUfbyDelta)/sum(mesh.magSf())).value()
*runTime.deltaT().value();
}
Info<< "St courant Number mean: " << meanStCoNum
<< " max: " << StCoNum << endl;
Info<< "St courant Number mean: " << meanStCoNum
<< " max: " << StCoNum << endl;
}
// ************************************************************************* //

4
applications/solvers/combustion/PDRFoam/bEqn.H
View File

@ -1,7 +1,7 @@
tmp<fv::convectionScheme<scalar> > mvConvection
(
fv::convectionScheme<scalar>::New
(
(
mesh,
fields,
phi,
@ -25,7 +25,7 @@ if (ign.ignited())
// Unburnt gas density
// ~~~~~~~~~~~~~~~~~~~
volScalarField rhou = thermo->rhou();
volScalarField rhou = thermo.rhou();
// Calculate flame normal etc.
// ~~~~~~~~~~~~~~~~~~~~~~~~~~~

17
applications/solvers/combustion/PDRFoam/createFields.H
View File

@ -1,10 +1,11 @@
Info<< "Reading thermophysical properties\n" << endl;
autoPtr<hhuCombustionThermo> thermo
autoPtr<hhuCombustionThermo> pThermo
(
hhuCombustionThermo::New(mesh)
);
combustionMixture& composition = thermo->composition();
hhuCombustionThermo& thermo = pThermo();
basicMultiComponentMixture& composition = thermo.composition();
volScalarField rho
(
@ -16,13 +17,13 @@
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
thermo->rho()
thermo.rho()
);
volScalarField& p = thermo->p();
const volScalarField& psi = thermo->psi();
volScalarField& h = thermo->h();
volScalarField& hu = thermo->hu();
volScalarField& p = thermo.p();
const volScalarField& psi = thermo.psi();
volScalarField& h = thermo.h();
volScalarField& hu = thermo.hu();
volScalarField& b = composition.Y("b");
Info<< "min(b) = " << min(b).value() << endl;
@ -54,7 +55,7 @@
rho,
U,
phi,
thermo()
thermo
)
);

2
applications/solvers/combustion/PDRFoam/hEqn.H
View File

@ -8,5 +8,5 @@
betav*DpDt
);
thermo->correct();
thermo.correct();
}

2
applications/solvers/combustion/PDRFoam/huEqn.H
View File

@ -13,6 +13,6 @@ if (ign.ignited())
//+ fvm::Sp(fvc::div(muEff*fvc::grad(b)/(b + 0.001)), hu)
==
betav*DpDt*rho/thermo->rhou()
betav*DpDt*rho/thermo.rhou()
);
}

3
applications/solvers/combustion/PDRFoam/laminarFlameSpeed/SCOPE/SCOPELaminarFlameSpeed.H
View File

@ -196,8 +196,7 @@ public:
// Destructor
~SCOPE();
~SCOPE();
// Member functions

6
applications/solvers/combustion/PDRFoam/pEqn.H
View File

@ -1,4 +1,4 @@
rho = thermo->rho();
rho = thermo.rho();
volScalarField rUA = 1.0/UEqn.A();
U = invA & UEqn.H();
@ -8,7 +8,7 @@ if (transonic)
surfaceScalarField phid
(
"phid",
fvc::interpolate(thermo->psi())
fvc::interpolate(psi)
*(
(fvc::interpolate(U) & mesh.Sf())
+ fvc::ddtPhiCorr(rUA, rho, U, phi)
@ -34,7 +34,7 @@ if (transonic)
}
else
{
phi =
phi =
fvc::interpolate(rho)*
(
(fvc::interpolate(U) & mesh.Sf())

4
applications/solvers/combustion/XiFoam/Make/options
View File

@ -2,7 +2,7 @@ EXE_INC = \
-I$(LIB_SRC)/engine/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/combustion/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
-I$(LIB_SRC)/turbulenceModels/compressible/turbulenceModel \
-I$(LIB_SRC)/thermophysicalModels/laminarFlameSpeed/lnInclude \
-I$(LIB_SRC)/sampling/lnInclude \
@ -13,7 +13,7 @@ EXE_LIBS = \
-lcompressibleRASModels \
-lcompressibleLESModels \
-lbasicThermophysicalModels \
-lcombustionThermophysicalModels \
-lreactionThermophysicalModels \
-lspecie \
-llaminarFlameSpeedModels \
-lfiniteVolume \

2
applications/solvers/combustion/XiFoam/XiFoam.C
View File

@ -109,7 +109,7 @@ int main(int argc, char *argv[])
turbulence->correct();
rho = thermo->rho();
rho = thermo.rho();
runTime.write();

6
applications/solvers/combustion/XiFoam/bEqn.H
View File

@ -6,7 +6,7 @@ if (ign.ignited())
// Unburnt gas density
// ~~~~~~~~~~~~~~~~~~~
volScalarField rhou = thermo->rhou();
volScalarField rhou = thermo.rhou();
// Calculate flame normal etc.
// ~~~~~~~~~~~~~~~~~~~~~~~~~~~
@ -76,7 +76,7 @@ if (ign.ignited())
volScalarField epsilon = pow(uPrimeCoef, 3)*turbulence->epsilon();
volScalarField tauEta = sqrt(thermo->muu()/(rhou*epsilon));
volScalarField tauEta = sqrt(thermo.muu()/(rhou*epsilon));
volScalarField Reta = up/
(
@ -180,7 +180,7 @@ if (ign.ignited())
// with a linear correction function to give a plausible profile for Xi
// ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
volScalarField XiEqStar =
volScalarField XiEqStar =
scalar(1.001) + XiCoef*sqrt(up/(Su + SuMin))*Reta;
volScalarField XiEq =

19
applications/solvers/combustion/XiFoam/createFields.H
View File

@ -1,10 +1,11 @@
Info<< "Reading thermophysical properties\n" << endl;
autoPtr<hhuCombustionThermo> thermo
autoPtr<hhuCombustionThermo> pThermo
(
hhuCombustionThermo::New(mesh)
);
combustionMixture& composition = thermo->composition();
hhuCombustionThermo& thermo = pThermo();
basicMultiComponentMixture& composition = thermo.composition();
volScalarField rho
(
@ -16,18 +17,18 @@
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
thermo->rho()
thermo.rho()
);
volScalarField& p = thermo->p();
const volScalarField& psi = thermo->psi();
volScalarField& h = thermo->h();
volScalarField& hu = thermo->hu();
volScalarField& p = thermo.p();
const volScalarField& psi = thermo.psi();
volScalarField& h = thermo.h();
volScalarField& hu = thermo.hu();
volScalarField& b = composition.Y("b");
Info<< "min(b) = " << min(b).value() << endl;
const volScalarField& T = thermo->T();
const volScalarField& T = thermo.T();
Info<< "\nReading field U\n" << endl;
@ -55,7 +56,7 @@
rho,
U,
phi,
thermo()
thermo
)
);

2
applications/solvers/combustion/XiFoam/hEqn.H
View File

@ -8,5 +8,5 @@
DpDt
);
thermo->correct();
thermo.correct();
}

2
applications/solvers/combustion/XiFoam/huEqn.H
View File

@ -13,6 +13,6 @@ if (ign.ignited())
//+ fvm::Sp(fvc::div(muEff*fvc::grad(b)/(b + 0.001)), hu)
==
DpDt*rho/thermo->rhou()
DpDt*rho/thermo.rhou()
);
}

10
applications/solvers/combustion/XiFoam/pEqn.H
View File

@ -1,4 +1,4 @@
rho = thermo->rho();
rho = thermo.rho();
volScalarField rUA = 1.0/UEqn.A();
U = rUA*UEqn.H();
@ -8,7 +8,7 @@ if (transonic)
surfaceScalarField phid
(
"phid",
fvc::interpolate(thermo->psi())
fvc::interpolate(psi)
*(
(fvc::interpolate(U) & mesh.Sf())
+ fvc::ddtPhiCorr(rUA, rho, U, phi)
@ -34,9 +34,9 @@ if (transonic)
}
else
{
phi =
fvc::interpolate(rho)*
(
phi =
fvc::interpolate(rho)
*(
(fvc::interpolate(U) & mesh.Sf())
+ fvc::ddtPhiCorr(rUA, rho, U, phi)
);

48
applications/solvers/combustion/coalChemistryFoam/additionalOutput.H
View File

@ -1,48 +0,0 @@
{
tmp<volScalarField> tdQ
(
new volScalarField
(
IOobject
(
"dQ",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar
(
"zero",
dimensionSet(1, -3, -1, 0, 0, 0, 0),
0.0
)
)
);
scalarField& dQ = tdQ();
scalarField cp(dQ.size(), 0.0);
forAll(Y, i)
{
volScalarField RRi = chemistry.RR(i);
forAll(h, celli)
{
scalar Ti = T[celli];
cp[celli] += Y[i][celli]*chemistry.specieThermo()[i].Cp(Ti);
scalar hi = chemistry.specieThermo()[i].h(Ti);
scalar RR = RRi[celli];
dQ[celli] -= hi*RR;
}
}
forAll(dQ, celli)
{
dQ[celli] /= cp[celli];
}
tdQ().write();
}

2
applications/solvers/combustion/coldEngineFoam/Make/options
View File

@ -4,7 +4,7 @@ EXE_INC = \
-I$(LIB_SRC)/engine/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/combustion/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
-I$(LIB_SRC)/turbulenceModels/compressible/turbulenceModel \
-I$(LIB_SRC)/finiteVolume/lnInclude

40
applications/solvers/combustion/coldEngineFoam/coldEngineFoam.C
View File

@ -33,7 +33,7 @@ Description
#include "fvCFD.H"
#include "engineTime.H"
#include "engineMesh.H"
#include "basicThermo.H"
#include "basicPsiThermo.H"
#include "turbulenceModel.H"
#include "OFstream.H"
@ -41,27 +41,27 @@ Description
int main(int argc, char *argv[])
{
# include "setRootCase.H"
#include "setRootCase.H"
# include "createEngineTime.H"
# include "createEngineMesh.H"
# include "createFields.H"
# include "initContinuityErrs.H"
# include "readEngineTimeControls.H"
# include "compressibleCourantNo.H"
# include "setInitialDeltaT.H"
# include "startSummary.H"
#include "createEngineTime.H"
#include "createEngineMesh.H"
#include "createFields.H"
#include "initContinuityErrs.H"
#include "readEngineTimeControls.H"
#include "compressibleCourantNo.H"
#include "setInitialDeltaT.H"
#include "startSummary.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
while (runTime.run())
{
# include "readPISOControls.H"
# include "readEngineTimeControls.H"
# include "compressibleCourantNo.H"
# include "setDeltaT.H"
#include "readPISOControls.H"
#include "readEngineTimeControls.H"
#include "compressibleCourantNo.H"
#include "setDeltaT.H"
runTime++;
@ -70,22 +70,22 @@ int main(int argc, char *argv[])
mesh.move();
# include "rhoEqn.H"
#include "rhoEqn.H"
# include "UEqn.H"
#include "UEqn.H"
// --- PISO loop
for (int corr=1; corr<=nCorr; corr++)
{
# include "hEqn.H"
# include "pEqn.H"
#include "hEqn.H"
#include "pEqn.H"
}
turbulence->correct();
runTime.write();
# include "logSummary.H"
#include "logSummary.H"
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"

19
applications/solvers/combustion/coldEngineFoam/createFields.H
View File

@ -1,9 +1,10 @@
Info<< "Reading thermophysical properties\n" << endl;
autoPtr<basicThermo> thermo
autoPtr<basicPsiThermo> pThermo
(
basicThermo::New(mesh)
basicPsiThermo::New(mesh)
);
basicPsiThermo& thermo = pThermo();
volScalarField rho
(
@ -15,13 +16,13 @@
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
thermo->rho()
thermo.rho()
);
volScalarField& p = thermo->p();
const volScalarField& psi = thermo->psi();
volScalarField& h = thermo->h();
const volScalarField& T = thermo->T();
volScalarField& p = thermo.p();
const volScalarField& psi = thermo.psi();
volScalarField& h = thermo.h();
const volScalarField& T = thermo.T();
Info<< "\nReading field U\n" << endl;
@ -38,7 +39,7 @@
mesh
);
# include "compressibleCreatePhi.H"
#include "compressibleCreatePhi.H"
Info<< "Creating turbulence model\n" << endl;
@ -49,7 +50,7 @@
rho,
U,
phi,
thermo()
thermo
)
);

2
applications/solvers/combustion/coldEngineFoam/hEqn.H
View File

@ -8,5 +8,5 @@
DpDt
);
thermo->correct();
thermo.correct();
}

6
applications/solvers/combustion/dieselEngineFoam/Make/options
View File

@ -7,10 +7,10 @@ EXE_INC = \
-I$(LIB_SRC)/thermophysicalModels/liquidMixture/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/thermophysicalFunctions/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/combustion/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/laminarFlameSpeed/lnInclude \
-I$(LIB_SRC)/../applications/solvers/combustion/XiFoam \
-I$(LIB_SRC)/../applications/solvers/reactionThermo/XiFoam \
-I$(LIB_SRC)/thermophysicalModels/chemistryModel/lnInclude \
-I$(LIB_SRC)/ODE/lnInclude \
-I$(LIB_SRC)/engine/lnInclude \
@ -20,7 +20,7 @@ EXE_LIBS = \
-lengine \
-lcompressibleRASModels \
-lcompressibleLESModels \
-lcombustionThermophysicalModels \
-lreactionThermophysicalModels \
-lfiniteVolume \
-llagrangian \
-ldieselSpray \

5
applications/solvers/combustion/dieselEngineFoam/YEqn.H
View File

@ -14,7 +14,7 @@ tmp<fv::convectionScheme<scalar> > mvConvection
label inertIndex = -1;
volScalarField Yt = 0.0*Y[0];
for(label i=0; i<Y.size(); i++)
for (label i=0; i<Y.size(); i++)
{
if (Y[i].name() != inertSpecie)
{
@ -39,8 +39,7 @@ tmp<fv::convectionScheme<scalar> > mvConvection
inertIndex = i;
}
}
Y[inertIndex] = scalar(1) - Yt;
Y[inertIndex].max(0.0);
}

50
applications/solvers/combustion/dieselEngineFoam/createFields.H
View File

@ -1,13 +1,17 @@
Info<< nl << "Reading thermophysicalProperties" << endl;
autoPtr<hCombustionThermo> thermo
(
hCombustionThermo::New(mesh)
);
combustionMixture& composition = thermo->composition();
autoPtr<psiChemistryModel> pChemistry
(
psiChemistryModel::New(mesh)
);
psiChemistryModel& chemistry = pChemistry();
hCombustionThermo& thermo = chemistry.thermo();
basicMultiComponentMixture& composition = thermo.composition();
PtrList<volScalarField>& Y = composition.Y();
word inertSpecie(thermo->lookup("inertSpecie"));
word inertSpecie(thermo.lookup("inertSpecie"));
volScalarField rho
(
@ -17,7 +21,7 @@ volScalarField rho
runTime.timeName(),
mesh
),
thermo->rho()
thermo.rho()
);
Info<< "Reading field U\n" << endl;
@ -35,10 +39,10 @@ volVectorField U
);
volScalarField& p = thermo->p();
const volScalarField& psi = thermo->psi();
const volScalarField& T = thermo->T();
volScalarField& h = thermo->h();
volScalarField& p = thermo.p();
const volScalarField& psi = thermo.psi();
const volScalarField& T = thermo.T();
volScalarField& h = thermo.h();
#include "compressibleCreatePhi.H"
@ -65,7 +69,7 @@ autoPtr<compressible::turbulenceModel> turbulence
rho,
U,
phi,
thermo()
thermo
)
);
@ -73,31 +77,11 @@ Info<< "Creating field DpDt\n" << endl;
volScalarField DpDt =
fvc::DDt(surfaceScalarField("phiU", phi/fvc::interpolate(rho)), p);
Info << "Constructing chemical mechanism" << endl;
chemistryModel chemistry
(
thermo(),
rho
);
multivariateSurfaceInterpolationScheme<scalar>::fieldTable fields;
for(label i=0; i<Y.size(); i++)
forAll (Y, i)
{
fields.add(Y[i]);
}
fields.add(h);
volScalarField dQ
(
IOobject
(
"dQ",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("zero", dimensionSet(1,-3,-1,0,0,0,0), 0.0)
);

7
applications/solvers/combustion/dieselEngineFoam/createSpray.H
View File

@ -1,14 +1,15 @@
Info << "Constructing Spray" << endl;
PtrList<specieProperties> gasProperties(Y.size());
PtrList<gasThermoPhysics> gasProperties(Y.size());
forAll(gasProperties, i)
{
gasProperties.set
(
i,
new specieProperties
new gasThermoPhysics
(
dynamic_cast<const reactingMixture&>(thermo()).speciesData()[i]
dynamic_cast<const reactingMixture<gasThermoPhysics>&>
(thermo).speciesData()[i]
)
);
}

53
applications/solvers/combustion/dieselEngineFoam/dieselEngineFoam.C
View File

@ -36,40 +36,41 @@ Description
#include "hCombustionThermo.H"
#include "turbulenceModel.H"
#include "spray.H"
#include "chemistryModel.H"
#include "psiChemistryModel.H"
#include "chemistrySolver.H"
#include "multivariateScheme.H"
#include "Switch.H"
#include "OFstream.H"
#include "volPointInterpolation.H"
#include "thermoPhysicsTypes.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
# include "setRootCase.H"
# include "createEngineTime.H"
# include "createEngineMesh.H"
# include "createFields.H"
# include "readEnvironmentalProperties.H"
# include "readCombustionProperties.H"
# include "createSpray.H"
# include "initContinuityErrs.H"
# include "readEngineTimeControls.H"
# include "compressibleCourantNo.H"
# include "setInitialDeltaT.H"
# include "startSummary.H"
#include "setRootCase.H"
#include "createEngineTime.H"
#include "createEngineMesh.H"
#include "createFields.H"
#include "readEnvironmentalProperties.H"
#include "readCombustionProperties.H"
#include "createSpray.H"
#include "initContinuityErrs.H"
#include "readEngineTimeControls.H"
#include "compressibleCourantNo.H"
#include "setInitialDeltaT.H"
#include "startSummary.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info << "\nStarting time loop\n" << endl;
while (runTime.run())
{
# include "readPISOControls.H"
# include "readEngineTimeControls.H"
# include "compressibleCourantNo.H"
# include "setDeltaT.H"
#include "readPISOControls.H"
#include "readEngineTimeControls.H"
#include "compressibleCourantNo.H"
#include "setDeltaT.H"
runTime++;
@ -101,27 +102,27 @@ int main(int argc, char *argv[])
kappa = (runTime.deltaT() + tc)/(runTime.deltaT() + tc + tk);
}
# include "rhoEqn.H"
# include "UEqn.H"
#include "rhoEqn.H"
#include "UEqn.H"
for (label ocorr=1; ocorr <= nOuterCorr; ocorr++)
{
# include "YEqn.H"
# include "hEqn.H"
#include "YEqn.H"
#include "hEqn.H"
// --- PISO loop
for (int corr=1; corr<=nCorr; corr++)
{
# include "pEqn.H"
#include "pEqn.H"
}
}
turbulence->correct();
# include "logSummary.H"
# include "spraySummary.H"
#include "logSummary.H"
#include "spraySummary.H"
rho = thermo->rho();
rho = thermo.rho();
runTime.write();

29
applications/solvers/combustion/dieselEngineFoam/hEqn.H
View File

@ -9,32 +9,5 @@
+ dieselSpray.heatTransferSource()
);
thermo->correct();
forAll(dQ, i)
{
dQ[i] = 0.0;
}
scalarField cp(dQ.size(), 0.0);
forAll(Y, i)
{
volScalarField RRi = chemistry.RR(i);
forAll(h, celli)
{
scalar Ti = T[celli];
cp[celli] += Y[i][celli]*chemistry.specieThermo()[i].Cp(Ti);
scalar hi = chemistry.specieThermo()[i].h(Ti);
scalar RR = RRi[celli];
dQ[celli] -= hi*RR;
}
}
forAll(dQ, celli)
{
dQ[celli] /= cp[celli];
}
thermo.correct();
}

4
applications/solvers/combustion/dieselEngineFoam/pEqn.H
View File

@ -1,4 +1,4 @@
rho = thermo->rho();
rho = thermo.rho();
volScalarField A = UEqn.A();
U = UEqn.H()/A;
@ -8,7 +8,7 @@ if (transonic)
surfaceScalarField phid
(
"phid",
fvc::interpolate(thermo->psi())
fvc::interpolate(psi)
*((fvc::interpolate(U) & mesh.Sf()) - fvc::meshPhi(rho, U))
);

2
applications/solvers/combustion/dieselEngineFoam/rhoEqn.H
View File

@ -44,7 +44,7 @@ volScalarField Sevap
dimensionedScalar("zero", dimensionSet(1, -3, -1, 0, 0), 0.0)
);
for(label i=0; i<Y.size(); i++)
for (label i=0; i<Y.size(); i++)
{
if (dieselSpray.isLiquidFuel()[i])
{

48
applications/solvers/combustion/dieselEngineFoam/spraySummary.H
View File

@ -1,30 +1,30 @@
label Nparcels = dieselSpray.size();
reduce(Nparcels, sumOp<label>());
label Nparcels = dieselSpray.size();
reduce(Nparcels, sumOp<label>());
Info<< "\nNumber of parcels in system.... | "
<< Nparcels << endl
<< "Injected liquid mass........... | "
<< 1e6*dieselSpray.injectedMass(runTime.value()) << " mg" << endl
<< "Liquid Mass in system.......... | "
<< 1e6*dieselSpray.liquidMass() << " mg" << endl
<< "SMD, Dmax...................... | "
<< dieselSpray.smd()*1e6 << " mu, "
<< dieselSpray.maxD()*1e6 << " mu"
<< endl;
Info<< "\nNumber of parcels in system.... | "
<< Nparcels << endl
<< "Injected liquid mass........... | "
<< 1e6*dieselSpray.injectedMass(runTime.value()) << " mg" << endl
<< "Liquid Mass in system.......... | "
<< 1e6*dieselSpray.liquidMass() << " mg" << endl
<< "SMD, Dmax...................... | "
<< dieselSpray.smd()*1e6 << " mu, "
<< dieselSpray.maxD()*1e6 << " mu"
<< endl;
scalar evapMass =
dieselSpray.injectedMass(runTime.value())
- dieselSpray.liquidMass();
scalar evapMass =
dieselSpray.injectedMass(runTime.value())
- dieselSpray.liquidMass();
scalar gasMass = fvc::domainIntegrate(rho).value();
scalar gasMass = fvc::domainIntegrate(rho).value();
if (dieselSpray.twoD())
{
gasMass *= 2.0*mathematicalConstant::pi/dieselSpray.angleOfWedge();
}
if (dieselSpray.twoD())
{
gasMass *= 2.0*mathematicalConstant::pi/dieselSpray.angleOfWedge();
}
scalar addedMass = gasMass - gasMass0;
scalar addedMass = gasMass - gasMass0;
Info<< "Added gas mass................. | " << 1e6*addedMass << " mg"
<< nl << "Evaporation Continuity Error... | "
<< 1e6*(addedMass - evapMass) << " mg" << endl;
Info<< "Added gas mass................. | " << 1e6*addedMass << " mg"
<< nl << "Evaporation Continuity Error... | "
<< 1e6*(addedMass - evapMass) << " mg" << endl;

6
applications/solvers/combustion/dieselFoam/Make/options
View File

@ -8,17 +8,17 @@ EXE_INC = \
-I$(LIB_SRC)/thermophysicalModels/liquidMixture/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/thermophysicalFunctions/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/combustion/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/laminarFlameSpeed/lnInclude \
-I$(LIB_SRC)/../applications/solvers/combustion/XiFoam \
-I$(LIB_SRC)/../applications/solvers/reactionThermo/XiFoam \
-I$(LIB_SRC)/thermophysicalModels/chemistryModel/lnInclude \
-I$(LIB_SRC)/ODE/lnInclude
EXE_LIBS = \
-lcompressibleRASModels \
-lcompressibleLESModels \
-lcombustionThermophysicalModels \
-lreactionThermophysicalModels \
-llagrangian \
-ldieselSpray \
-lliquids \

47
applications/solvers/combustion/dieselFoam/dieselFoam.C
View File

@ -34,7 +34,7 @@ Description
#include "hCombustionThermo.H"
#include "turbulenceModel.H"
#include "spray.H"
#include "chemistryModel.H"
#include "psiChemistryModel.H"
#include "chemistrySolver.H"
#include "multivariateScheme.H"
@ -46,28 +46,27 @@ Description
int main(int argc, char *argv[])
{
#include "setRootCase.H"
#include "createTime.H"
#include "createMesh.H"
#include "createFields.H"
#include "readEnvironmentalProperties.H"
#include "readCombustionProperties.H"
#include "createSpray.H"
#include "initContinuityErrs.H"
#include "readTimeControls.H"
#include "compressibleCourantNo.H"
#include "setInitialDeltaT.H"
# include "setRootCase.H"
# include "createTime.H"
# include "createMesh.H"
# include "createFields.H"
# include "readEnvironmentalProperties.H"
# include "readCombustionProperties.H"
# include "createSpray.H"
# include "initContinuityErrs.H"
# include "readTimeControls.H"
# include "compressibleCourantNo.H"
# include "setInitialDeltaT.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info << "\nStarting time loop\n" << endl;
while (runTime.run())
{
# include "readPISOControls.H"
# include "compressibleCourantNo.H"
# include "setDeltaT.H"
#include "readPISOControls.H"
#include "compressibleCourantNo.H"
#include "setDeltaT.H"
runTime++;
Info<< "Time = " << runTime.timeName() << nl << endl;
@ -94,26 +93,26 @@ int main(int argc, char *argv[])
kappa = (runTime.deltaT() + tc)/(runTime.deltaT()+tc+tk);
}
# include "rhoEqn.H"
# include "UEqn.H"
#include "rhoEqn.H"
#include "UEqn.H"
for (label ocorr=1; ocorr <= nOuterCorr; ocorr++)
{
# include "YEqn.H"
# include "hEqn.H"
#include "YEqn.H"
#include "hEqn.H"
// --- PISO loop
for (int corr=1; corr<=nCorr; corr++)
{
# include "pEqn.H"
#include "pEqn.H"
}
}
turbulence->correct();
# include "spraySummary.H"
#include "spraySummary.H"
rho = thermo->rho();
rho = thermo.rho();
runTime.write();

10
applications/solvers/combustion/dieselFoam/pEqn.H
View File

@ -1,4 +1,4 @@
rho = thermo->rho();
rho = thermo.rho();
volScalarField rUA = 1.0/UEqn.A();
U = rUA*UEqn.H();
@ -8,7 +8,7 @@ if (transonic)
surfaceScalarField phid
(
"phid",
fvc::interpolate(thermo->psi())
fvc::interpolate(psi)
*(
(fvc::interpolate(U) & mesh.Sf())
+ fvc::ddtPhiCorr(rUA, rho, U, phi)
@ -36,9 +36,9 @@ if (transonic)
}
else
{
phi =
fvc::interpolate(rho)*
(
phi =
fvc::interpolate(rho)
*(
(fvc::interpolate(U) & mesh.Sf())
+ fvc::ddtPhiCorr(rUA, rho, U, phi)
);

4
applications/solvers/combustion/engineFoam/Make/options
View File

@ -3,7 +3,7 @@ EXE_INC = \
-I$(LIB_SRC)/engine/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/combustion/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
-I$(LIB_SRC)/turbulenceModels/compressible/turbulenceModel \
-I$(LIB_SRC)/thermophysicalModels/laminarFlameSpeed/lnInclude \
-I$(LIB_SRC)/finiteVolume/lnInclude
@ -13,7 +13,7 @@ EXE_LIBS = \
-lcompressibleRASModels \
-lcompressibleLESModels \
-lbasicThermophysicalModels \
-lcombustionThermophysicalModels \
-lreactionThermophysicalModels \
-lspecie \
-llaminarFlameSpeedModels \
-lfiniteVolume

50
applications/solvers/combustion/engineFoam/engineFoam.C
View File

@ -63,29 +63,29 @@ Description
int main(int argc, char *argv[])
{
# include "setRootCase.H"
#include "setRootCase.H"
# include "createEngineTime.H"
# include "createEngineMesh.H"
# include "readPISOControls.H"
# include "readCombustionProperties.H"
# include "createFields.H"
# include "initContinuityErrs.H"
# include "readEngineTimeControls.H"
# include "compressibleCourantNo.H"
# include "setInitialDeltaT.H"
# include "startSummary.H"
#include "createEngineTime.H"
#include "createEngineMesh.H"
#include "readPISOControls.H"
#include "readCombustionProperties.H"
#include "createFields.H"
#include "initContinuityErrs.H"
#include "readEngineTimeControls.H"
#include "compressibleCourantNo.H"
#include "setInitialDeltaT.H"
#include "startSummary.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info << "\nStarting time loop\n" << endl;
while (runTime.run())
{
# include "readPISOControls.H"
# include "readEngineTimeControls.H"
# include "compressibleCourantNo.H"
# include "setDeltaT.H"
#include "readPISOControls.H"
#include "readEngineTimeControls.H"
#include "compressibleCourantNo.H"
#include "setDeltaT.H"
runTime++;
@ -93,31 +93,31 @@ int main(int argc, char *argv[])
mesh.move();
# include "rhoEqn.H"
#include "rhoEqn.H"
# include "UEqn.H"
#include "UEqn.H"
// --- PISO loop
for (int corr=1; corr<=nCorr; corr++)
{
# include "ftEqn.H"
# include "bEqn.H"
# include "huEqn.H"
# include "hEqn.H"
#include "ftEqn.H"
#include "bEqn.H"
#include "huEqn.H"
#include "hEqn.H"
if (!ign.ignited())
{
hu == h;
}
# include "pEqn.H"
#include "pEqn.H"
}
turbulence->correct();
# include "logSummary.H"
#include "logSummary.H"
rho = thermo->rho();
rho = thermo.rho();
runTime.write();

6
applications/solvers/combustion/engineFoam/pEqn.H
View File

@ -1,4 +1,4 @@
rho = thermo->rho();
rho = thermo.rho();
volScalarField rUA = 1.0/UEqn.A();
U = rUA*UEqn.H();
@ -8,8 +8,8 @@ if (transonic)
surfaceScalarField phid
(
"phid",
fvc::interpolate(thermo->psi())
*((fvc::interpolate(U) & mesh.Sf()) - fvc::meshPhi(rho, U))
fvc::interpolate(psi)
*((fvc::interpolate(U) & mesh.Sf()) - fvc::meshPhi(rho, U))
);
for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)

2
applications/solvers/combustion/engineFoam/readEngineTimeControls.H
View File

@ -1,3 +1,3 @@
#include "readTimeControls.H"
#include "readTimeControls.H"
maxDeltaT = runTime.userTimeToTime(maxDeltaT);

4
applications/solvers/combustion/reactingFoam/Make/options
View File

@ -2,7 +2,7 @@ EXE_INC = \
-I../XiFoam \
-I$(LIB_SRC)/turbulenceModels/compressible/turbulenceModel \
-I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/combustion/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/chemistryModel/lnInclude \
-I$(LIB_SRC)/ODE/lnInclude \
@ -11,7 +11,7 @@ EXE_INC = \
EXE_LIBS = \
-lcompressibleRASModels \
-lcompressibleLESModels \
-lcombustionThermophysicalModels \
-lreactionThermophysicalModels \
-lspecie \
-lbasicThermophysicalModels \
-lchemistryModel \

4
applications/solvers/combustion/reactingFoam/YEqn.H
View File

@ -13,7 +13,7 @@ tmp<fv::convectionScheme<scalar> > mvConvection
label inertIndex = -1;
volScalarField Yt = 0.0*Y[0];
for(label i=0; i<Y.size(); i++)
for (label i=0; i<Y.size(); i++)
{
if (Y[i].name() != inertSpecie)
{
@ -37,7 +37,7 @@ tmp<fv::convectionScheme<scalar> > mvConvection
inertIndex = i;
}
}
Y[inertIndex] = scalar(1) - Yt;
Y[inertIndex].max(0.0);
}

43
applications/solvers/combustion/reactingFoam/createFields.H
View File

@ -1,13 +1,16 @@
Info<< nl << "Reading thermophysicalProperties" << endl;
autoPtr<hCombustionThermo> thermo
autoPtr<psiChemistryModel> pChemistry
(
hCombustionThermo::New(mesh)
psiChemistryModel::New(mesh)
);
psiChemistryModel& chemistry = pChemistry();
combustionMixture& composition = thermo->composition();
hCombustionThermo& thermo = chemistry.thermo();
basicMultiComponentMixture& composition = thermo.composition();
PtrList<volScalarField>& Y = composition.Y();
word inertSpecie(thermo->lookup("inertSpecie"));
word inertSpecie(thermo.lookup("inertSpecie"));
volScalarField rho
(
@ -17,7 +20,7 @@ volScalarField rho
runTime.timeName(),
mesh
),
thermo->rho()
thermo.rho()
);
Info<< "Reading field U\n" << endl;
@ -35,9 +38,9 @@ volVectorField U
);
volScalarField& p = thermo->p();
const volScalarField& psi = thermo->psi();
volScalarField& h = thermo->h();
volScalarField& p = thermo.p();
const volScalarField& psi = thermo.psi();
volScalarField& h = thermo.h();
#include "compressibleCreatePhi.H"
@ -64,7 +67,7 @@ autoPtr<compressible::turbulenceModel> turbulence
rho,
U,
phi,
thermo()
thermo
)
);
@ -72,31 +75,11 @@ Info<< "Creating field DpDt\n" << endl;
volScalarField DpDt =
fvc::DDt(surfaceScalarField("phiU", phi/fvc::interpolate(rho)), p);
Info << "Constructing chemical mechanism" << endl;
chemistryModel chemistry
(
thermo(),
rho
);
multivariateSurfaceInterpolationScheme<scalar>::fieldTable fields;
for(label i=0; i<Y.size(); i++)
forAll (Y, i)
{
fields.add(Y[i]);
}
fields.add(h);
volScalarField dQ
(
IOobject
(
"dQ",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("zero", dimensionSet(1, -3, -1, 0, 0, 0, 0), 0.0)
);

48
applications/solvers/combustion/reactingFoam/reactingFoam.C
View File

@ -33,7 +33,7 @@ Description
#include "fvCFD.H"
#include "hCombustionThermo.H"
#include "turbulenceModel.H"
#include "chemistryModel.H"
#include "psiChemistryModel.H"
#include "chemistrySolver.H"
#include "multivariateScheme.H"
@ -41,52 +41,52 @@ Description
int main(int argc, char *argv[])
{
# include "setRootCase.H"
# include "createTime.H"
# include "createMesh.H"
# include "readChemistryProperties.H"
# include "readEnvironmentalProperties.H"
# include "createFields.H"
# include "initContinuityErrs.H"
# include "readTimeControls.H"
# include "compressibleCourantNo.H"
# include "setInitialDeltaT.H"
#include "setRootCase.H"
#include "createTime.H"
#include "createMesh.H"
#include "readChemistryProperties.H"
#include "readEnvironmentalProperties.H"
#include "createFields.H"
#include "initContinuityErrs.H"
#include "readTimeControls.H"
#include "compressibleCourantNo.H"
#include "setInitialDeltaT.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info << "\nStarting time loop\n" << endl;
while (runTime.run())
{
# include "readTimeControls.H"
# include "readPISOControls.H"
# include "compressibleCourantNo.H"
# include "setDeltaT.H"
#include "readTimeControls.H"
#include "readPISOControls.H"
#include "compressibleCourantNo.H"
#include "setDeltaT.H"
runTime++;
Info<< "Time = " << runTime.timeName() << nl << endl;
# include "chemistry.H"
# include "rhoEqn.H"
# include "UEqn.H"
#include "chemistry.H"
#include "rhoEqn.H"
#include "UEqn.H"
for (label ocorr=1; ocorr <= nOuterCorr; ocorr++)
{
# include "YEqn.H"
#include "YEqn.H"
# define Db turbulence->alphaEff()
# include "hEqn.H"
#define Db turbulence->alphaEff()
#include "hEqn.H"
// --- PISO loop
for (int corr=1; corr<=nCorr; corr++)
{
# include "pEqn.H"
#include "pEqn.H"
}
}
turbulence->correct();
rho = thermo->rho();
rho = thermo.rho();
runTime.write();

3
applications/solvers/combustion/reactingFoam/readChemistryProperties.H
View File

@ -8,7 +8,8 @@ IOdictionary chemistryProperties
runTime.constant(),
mesh,
IOobject::MUST_READ,
IOobject::NO_WRITE
IOobject::NO_WRITE,
false
)
);

8
applications/solvers/compressible/rhoCentralFoam/compressibleCourantNo.H
View File

@ -38,11 +38,11 @@ if (mesh.nInternalFaces())
surfaceScalarField amaxSfbyDelta =
mesh.surfaceInterpolation::deltaCoeffs()*amaxSf;
CoNum = max(amaxSfbyDelta/mesh.magSf())
.value()*runTime.deltaT().value();
CoNum = max(amaxSfbyDelta/mesh.magSf()).value()*runTime.deltaT().value();
meanCoNum = (sum(amaxSfbyDelta)/sum(mesh.magSf()))
.value()*runTime.deltaT().value();
meanCoNum =
(sum(amaxSfbyDelta)/sum(mesh.magSf())).value()
*runTime.deltaT().value();
}
Info<< "Mean and max Courant Numbers = "

17
applications/solvers/compressible/rhoCentralFoam/createFields.H
View File

@ -1,15 +1,16 @@
Info<< "Reading thermophysical properties\n" << endl;
autoPtr<basicThermo> thermo
autoPtr<basicPsiThermo> pThermo
(
basicThermo::New(mesh)
basicPsiThermo::New(mesh)
);
basicPsiThermo& thermo = pThermo();
volScalarField& p = thermo->p();
volScalarField& h = thermo->h();
const volScalarField& T = thermo->T();
const volScalarField& psi = thermo->psi();
const volScalarField& mu = thermo->mu();
volScalarField& p = thermo.p();
volScalarField& h = thermo.h();
const volScalarField& T = thermo.T();
const volScalarField& psi = thermo.psi();
const volScalarField& mu = thermo.mu();
bool inviscid(true);
if (max(mu.internalField()) > 0.0)
@ -42,7 +43,7 @@ volScalarField rho
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
thermo->rho(),
thermo.rho(),
rhoBoundaryTypes
);

5
applications/solvers/compressible/rhoCentralFoam/rhoBoundaryTypes.H
View File

@ -3,10 +3,7 @@ wordList rhoBoundaryTypes = pbf.types();
forAll(rhoBoundaryTypes, patchi)
{
if
(
rhoBoundaryTypes[patchi] == "waveTransmissive"
)
if (rhoBoundaryTypes[patchi] == "waveTransmissive")
{
rhoBoundaryTypes[patchi] = zeroGradientFvPatchScalarField::typeName;
}

37
applications/solvers/compressible/rhoCentralFoam/rhoCentralFoam.C
View File

@ -32,7 +32,7 @@ Description
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "basicThermo.H"
#include "basicPsiThermo.H"
#include "zeroGradientFvPatchFields.H"
#include "fixedRhoFvPatchScalarField.H"
@ -40,18 +40,17 @@ Description
int main(int argc, char *argv[])
{
#include "setRootCase.H"
# include "setRootCase.H"
#include "createTime.H"
#include "createMesh.H"
#include "createFields.H"
#include "readThermophysicalProperties.H"
#include "readTimeControls.H"
# include "createTime.H"
# include "createMesh.H"
# include "createFields.H"
# include "readThermophysicalProperties.H"
# include "readTimeControls.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
# include "readFluxScheme.H"
#include "readFluxScheme.H"
dimensionedScalar v_zero("v_zero",dimVolume/dimTime, 0.0);
@ -91,7 +90,7 @@ int main(int argc, char *argv[])
surfaceScalarField phiv_pos = U_pos & mesh.Sf();
surfaceScalarField phiv_neg = U_neg & mesh.Sf();
volScalarField c = sqrt(thermo->Cp()/thermo->Cv()*rPsi);
volScalarField c = sqrt(thermo.Cp()/thermo.Cv()*rPsi);
surfaceScalarField cSf_pos = fvc::interpolate(c, pos, "reconstruct(T)")*mesh.magSf();
surfaceScalarField cSf_neg = fvc::interpolate(c, neg, "reconstruct(T)")*mesh.magSf();
@ -102,9 +101,9 @@ int main(int argc, char *argv[])
surfaceScalarField amaxSf("amaxSf", max(mag(am), mag(ap)));
# include "compressibleCourantNo.H"
# include "readTimeControls.H"
# include "setDeltaT.H"
#include "compressibleCourantNo.H"
#include "readTimeControls.H"
#include "setDeltaT.H"
runTime++;
@ -183,7 +182,7 @@ int main(int argc, char *argv[])
h = (rhoE + p)/rho - 0.5*magSqr(U);
h.correctBoundaryConditions();
thermo->correct();
thermo.correct();
rhoE.boundaryField() =
rho.boundaryField()*
(
@ -193,15 +192,15 @@ int main(int argc, char *argv[])
if (!inviscid)
{
volScalarField k("k", thermo->Cp()*mu/Pr);
volScalarField k("k", thermo.Cp()*mu/Pr);
solve
(
fvm::ddt(rho, h) - fvc::ddt(rho, h)
- fvm::laplacian(thermo->alpha(), h)
+ fvc::laplacian(thermo->alpha(), h)
- fvm::laplacian(thermo.alpha(), h)
+ fvc::laplacian(thermo.alpha(), h)
- fvc::laplacian(k, T)
);
thermo->correct();
thermo.correct();
rhoE = rho*(h + 0.5*magSqr(U)) - p;
}

15
applications/solvers/compressible/rhoPimpleFoam/createFields.H
View File

@ -1,13 +1,14 @@
Info<< "Reading thermophysical properties\n" << endl;
autoPtr<basicThermo> thermo
autoPtr<basicPsiThermo> pThermo
(
basicThermo::New(mesh)
basicPsiThermo::New(mesh)
);
basicPsiThermo& thermo = pThermo();
volScalarField& p = thermo->p();
volScalarField& h = thermo->h();
const volScalarField& psi = thermo->psi();
volScalarField& p = thermo.p();
volScalarField& h = thermo.h();
const volScalarField& psi = thermo.psi();
volScalarField rho
(
@ -19,7 +20,7 @@
IOobject::READ_IF_PRESENT,
IOobject::AUTO_WRITE
),
thermo->rho()
thermo.rho()
);
Info<< "Reading field U\n" << endl;
@ -51,7 +52,7 @@
rho,
U,
phi,
thermo()
thermo
)
);

2
applications/solvers/compressible/rhoPimpleFoam/hEqn.H
View File

@ -19,5 +19,5 @@
hEqn.solve();
}
thermo->correct();
thermo.correct();
}

10
applications/solvers/compressible/rhoPimpleFoam/pEqn.H
View File

@ -1,4 +1,4 @@
rho = thermo->rho();
rho = thermo.rho();
volScalarField rUA = 1.0/UEqn().A();
U = rUA*UEqn().H();
@ -13,7 +13,7 @@ if (transonic)
surfaceScalarField phid
(
"phid",
fvc::interpolate(thermo->psi())
fvc::interpolate(psi)
*(
(fvc::interpolate(U) & mesh.Sf())
+ fvc::ddtPhiCorr(rUA, rho, U, phi)
@ -99,7 +99,7 @@ else
// Explicitly relax pressure for momentum corrector
p.relax();
rho = thermo->rho();
rho = thermo.rho();
rho.relax();
Info<< "rho max/min : " << max(rho).value()
<< " " << min(rho).value() << endl;
@ -117,7 +117,7 @@ bound(p, pMin);
/*
if (closedVolume)
{
p += (initialMass - fvc::domainIntegrate(thermo->psi()*p))
/fvc::domainIntegrate(thermo->psi());
p += (initialMass - fvc::domainIntegrate(psi*p))
/fvc::domainIntegrate(psi);
}
*/

2
applications/solvers/compressible/rhoPimpleFoam/rhoPimpleFoam.C
View File

@ -35,7 +35,7 @@ Description
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "basicThermo.H"
#include "basicPsiThermo.H"
#include "turbulenceModel.H"
#include "bound.H"

15
applications/solvers/compressible/rhoPisoFoam/createFields.H
View File

@ -1,13 +1,14 @@
Info<< "Reading thermophysical properties\n" << endl;
autoPtr<basicThermo> thermo
autoPtr<basicPsiThermo> pThermo
(
basicThermo::New(mesh)
basicPsiThermo::New(mesh)
);
basicPsiThermo& thermo = pThermo();
volScalarField& p = thermo->p();
volScalarField& h = thermo->h();
const volScalarField& psi = thermo->psi();
volScalarField& p = thermo.p();
volScalarField& h = thermo.h();
const volScalarField& psi = thermo.psi();
volScalarField rho
(
@ -19,7 +20,7 @@
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
thermo->rho()
thermo.rho()
);
Info<< "\nReading field U\n" << endl;
@ -47,7 +48,7 @@
rho,
U,
phi,
thermo()
thermo
)
);

2
applications/solvers/compressible/rhoPisoFoam/hEqn.H
View File

@ -8,5 +8,5 @@
DpDt
);
thermo->correct();
thermo.correct();
}

6
applications/solvers/compressible/rhoPisoFoam/pEqn.H
View File

@ -1,4 +1,4 @@
rho = thermo->rho();
rho = thermo.rho();
volScalarField rUA = 1.0/UEqn.A();
U = rUA*UEqn.H();
@ -8,7 +8,7 @@ if (transonic)
surfaceScalarField phid
(
"phid",
fvc::interpolate(thermo->psi())
fvc::interpolate(psi)
*(
(fvc::interpolate(U) & mesh.Sf())
+ fvc::ddtPhiCorr(rUA, rho, U, phi)
@ -34,7 +34,7 @@ if (transonic)
}
else
{
phi =
phi =
fvc::interpolate(rho)*
(
(fvc::interpolate(U) & mesh.Sf())

4
applications/solvers/compressible/rhoPisoFoam/rhoPisoFoam.C
View File

@ -31,7 +31,7 @@ Description
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "basicThermo.H"
#include "basicPsiThermo.H"
#include "turbulenceModel.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -77,7 +77,7 @@ int main(int argc, char *argv[])
turbulence->correct();
rho = thermo->rho();
rho = thermo.rho();
runTime.write();

16
applications/solvers/compressible/rhoPorousSimpleFoam/createFields.H
View File

@ -1,9 +1,10 @@
Info<< "Reading thermophysical properties\n" << endl;
autoPtr<basicThermo> thermo
autoPtr<basicPsiThermo> pThermo
(
basicThermo::New(mesh)
basicPsiThermo::New(mesh)
);
basicPsiThermo& thermo = pThermo();
volScalarField rho
(
@ -15,11 +16,12 @@
IOobject::READ_IF_PRESENT,
IOobject::AUTO_WRITE
),
thermo->rho()
thermo.rho()
);
volScalarField& p = thermo->p();
volScalarField& h = thermo->h();
volScalarField& p = thermo.p();
volScalarField& h = thermo.h();
const volScalarField& psi = thermo.psi();
Info<< "Reading field U\n" << endl;
@ -36,7 +38,7 @@
mesh
);
# include "compressibleCreatePhi.H"
#include "compressibleCreatePhi.H"
label pRefCell = 0;
@ -56,7 +58,7 @@
rho,
U,
phi,
thermo()
thermo
)
);

2
applications/solvers/compressible/rhoPorousSimpleFoam/hEqn.H
View File

@ -14,5 +14,5 @@
eqnResidual = hEqn.solve().initialResidual();
maxResidual = max(eqnResidual, maxResidual);
thermo->correct();
thermo.correct();
}

6
applications/solvers/compressible/rhoPorousSimpleFoam/pEqn.H
View File

@ -65,10 +65,10 @@ bound(p, pMin);
// to obey overall mass continuity
if (closedVolume)
{
p += (initialMass - fvc::domainIntegrate(thermo->psi()*p))
/fvc::domainIntegrate(thermo->psi());
p += (initialMass - fvc::domainIntegrate(psi*p))
/fvc::domainIntegrate(psi);
}
rho = thermo->rho();
rho = thermo.rho();
rho.relax();
Info<< "rho max/min : " << max(rho).value() << " " << min(rho).value() << endl;

4
applications/solvers/compressible/rhoPorousSimpleFoam/rhoPorousSimpleFoam.C
View File

@ -27,12 +27,12 @@ Application
Description
Steady-state solver for turbulent flow of compressible fluids with
implicit or explicit porosity treatment
RANS turbulence modelling, and implicit or explicit porosity treatment
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "basicThermo.H"
#include "basicPsiThermo.H"
#include "RASModel.H"
#include "porousZones.H"

15
applications/solvers/compressible/rhoSimpleFoam/createFields.H
View File

@ -1,9 +1,10 @@
Info<< "Reading thermophysical properties\n" << endl;
autoPtr<basicThermo> thermo
autoPtr<basicPsiThermo> pThermo
(
basicThermo::New(mesh)
basicPsiThermo::New(mesh)
);
basicPsiThermo& thermo = pThermo();
volScalarField rho
(
@ -15,12 +16,12 @@
IOobject::READ_IF_PRESENT,
IOobject::AUTO_WRITE
),
thermo->rho()
thermo.rho()
);
volScalarField& p = thermo->p();
volScalarField& h = thermo->h();
volScalarField& p = thermo.p();
volScalarField& h = thermo.h();
const volScalarField& psi = thermo.psi();
Info<< "Reading field U\n" << endl;
volVectorField U
@ -56,7 +57,7 @@
rho,
U,
phi,
thermo()
thermo
)
);

2
applications/solvers/compressible/rhoSimpleFoam/hEqn.H
View File

@ -14,5 +14,5 @@
eqnResidual = hEqn.solve().initialResidual();
maxResidual = max(eqnResidual, maxResidual);
thermo->correct();
thermo.correct();
}

10
applications/solvers/compressible/rhoSimpleFoam/pEqn.H
View File

@ -1,4 +1,4 @@
rho = thermo->rho();
rho = thermo.rho();
volScalarField rUA = 1.0/UEqn().A();
U = rUA*UEqn().H();
@ -11,7 +11,7 @@ if (transonic)
surfaceScalarField phid
(
"phid",
fvc::interpolate(thermo->psi())*(fvc::interpolate(U) & mesh.Sf())
fvc::interpolate(psi)*(fvc::interpolate(U) & mesh.Sf())
);
for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
@ -82,7 +82,7 @@ else
// Explicitly relax pressure for momentum corrector
p.relax();
rho = thermo->rho();
rho = thermo.rho();
rho.relax();
Info<< "rho max/min : " << max(rho).value() << " " << min(rho).value() << endl;
@ -95,6 +95,6 @@ bound(p, pMin);
// to obey overall mass continuity
if (closedVolume)
{
p += (initialMass - fvc::domainIntegrate(thermo->psi()*p))
/fvc::domainIntegrate(thermo->psi());
p += (initialMass - fvc::domainIntegrate(psi*p))
/fvc::domainIntegrate(psi);
}

4
applications/solvers/compressible/rhoSimpleFoam/rhoSimpleFoam.C
View File

@ -26,13 +26,13 @@ Application
rhoSimpleFoam
Description
Steady-state SIMPLE solver for laminar or turbulent flow of
Steady-state SIMPLE solver for laminar or turbulent RANS flow of
compressible fluids.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "basicThermo.H"
#include "basicPsiThermo.H"
#include "RASModel.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

2
applications/solvers/compressible/sonicDyMFoam/compressibleContinuityErrs.H
View File

@ -1,5 +1,5 @@
{
# include "rhoEqn.H"
#include "rhoEqn.H"
}
{
scalar sumLocalContErr =

15
applications/solvers/compressible/sonicDyMFoam/createFields.H
View File

@ -1,13 +1,14 @@
Info<< "Reading thermophysical properties\n" << endl;
autoPtr<basicThermo> thermo
autoPtr<basicPsiThermo> pThermo
(
basicThermo::New(mesh)
basicPsiThermo::New(mesh)
);
basicPsiThermo& thermo = pThermo();
volScalarField& p = thermo->p();
volScalarField& h = thermo->h();
const volScalarField& psi = thermo->psi();
volScalarField& p = thermo.p();
volScalarField& h = thermo.h();
const volScalarField& psi = thermo.psi();
volScalarField rho
(
@ -17,7 +18,7 @@
runTime.timeName(),
mesh
),
thermo->rho()
thermo.rho()
);
Info<< "Reading field U\n" << endl;
@ -45,7 +46,7 @@
rho,
U,
phi,
thermo()
thermo
)
);

6
applications/solvers/compressible/sonicDyMFoam/sonicDyMFoam.C
View File

@ -32,7 +32,7 @@ Description
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "basicThermo.H"
#include "basicPsiThermo.H"
#include "turbulenceModel.H"
#include "motionSolver.H"
@ -84,8 +84,8 @@ int main(int argc, char *argv[])
surfaceScalarField phid
(
"phid",
fvc::interpolate(psi)*
(
fvc::interpolate(psi)
*(
(fvc::interpolate(U) & mesh.Sf()) - fvc::meshPhi(rho, U)
)
);

8
applications/solvers/compressible/sonicFoam/UEqn.H Normal file
View File

@ -0,0 +1,8 @@
fvVectorMatrix UEqn
(
fvm::ddt(rho, U)
+ fvm::div(phi, U)
+ turbulence->divDevRhoReff(U)
);
solve(UEqn == -fvc::grad(p));

12
applications/solvers/compressible/sonicFoam/compressibleContinuityErrs.H
View File

@ -1,12 +0,0 @@
{
# include "rhoEqn.H"
}
{
scalar sumLocalContErr = (sum(mag(rho - psi*p))/sum(rho)).value();
scalar globalContErr = (sum(rho - psi*p)/sum(rho)).value();
cumulativeContErr += globalContErr;
Info<< "time step continuity errors : sum local = " << sumLocalContErr
<< ", global = " << globalContErr
<< ", cumulative = " << cumulativeContErr << endl;
}

21
applications/solvers/compressible/sonicFoam/createFields.H
View File

@ -1,13 +1,14 @@
Info<< "Reading thermophysical properties\n" << endl;
autoPtr<basicThermo> thermo
autoPtr<basicPsiThermo> pThermo
(
basicThermo::New(mesh)
basicPsiThermo::New(mesh)
);
basicPsiThermo& thermo = pThermo();
volScalarField& p = thermo->p();
volScalarField& h = thermo->h();
const volScalarField& psi = thermo->psi();
volScalarField& p = thermo.p();
volScalarField& e = thermo.e();
const volScalarField& psi = thermo.psi();
volScalarField rho
(
@ -17,7 +18,7 @@
runTime.timeName(),
mesh
),
thermo->rho()
thermo.rho()
);
Info<< "Reading field U\n" << endl;
@ -34,7 +35,7 @@
mesh
);
# include "compressibleCreatePhi.H"
#include "compressibleCreatePhi.H"
Info<< "Creating turbulence model\n" << endl;
@ -45,10 +46,6 @@
rho,
U,
phi,
thermo()
thermo
)
);
Info<< "Creating field DpDt\n" << endl;
volScalarField DpDt =
fvc::DDt(surfaceScalarField("phiU", phi/fvc::interpolate(rho)), p);

12
applications/solvers/compressible/sonicFoam/eEqn.H Normal file
View File

@ -0,0 +1,12 @@
{
solve
(
fvm::ddt(rho, e)
+ fvm::div(phi, e)
- fvm::laplacian(turbulence->alphaEff(), e)
==
- p*fvc::div(phi/fvc::interpolate(rho))
);
thermo.correct();
}

12
applications/solvers/compressible/sonicFoam/hEqn.H
View File

@ -1,12 +0,0 @@
{
solve
(
fvm::ddt(rho, h)
+ fvm::div(phi, h)
- fvm::laplacian(turbulence->alphaEff(), h)
==
DpDt
);
thermo->correct();
}

37
applications/solvers/compressible/sonicFoam/pEqn.H Normal file
View File

@ -0,0 +1,37 @@
rho = thermo.rho();
volScalarField rUA = 1.0/UEqn.A();
U = rUA*UEqn.H();
surfaceScalarField phid
(
"phid",
fvc::interpolate(psi)
*(
(fvc::interpolate(U) & mesh.Sf())
+ fvc::ddtPhiCorr(rUA, rho, U, phi)
)
);
for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
{
fvScalarMatrix pEqn
(
fvm::ddt(psi, p)
+ fvm::div(phid, p)
- fvm::laplacian(rho*rUA, p)
);
pEqn.solve();
if (nonOrth == nNonOrthCorr)
{
phi = pEqn.flux();
}
}
#include "rhoEqn.H"
#include "compressibleContinuityErrs.H"
U -= rUA*fvc::grad(p);
U.correctBoundaryConditions();

23
applications/solvers/compressible/sonicFoam/readThermodynamicProperties.H
View File

@ -1,23 +0,0 @@
Info<< "Reading thermodynamicProperties\n" << endl;
IOdictionary thermodynamicProperties
(
IOobject
(
"thermodynamicProperties",
runTime.constant(),
mesh,
IOobject::MUST_READ,
IOobject::NO_WRITE
)
);
dimensionedScalar R
(
thermodynamicProperties.lookup("R")
);
dimensionedScalar Cv
(
thermodynamicProperties.lookup("Cv")
);

18
applications/solvers/compressible/sonicFoam/readTransportProperties.H
View File

@ -1,18 +0,0 @@
Info<< "Reading transportProperties\n" << endl;
IOdictionary transportProperties
(
IOobject
(
"transportProperties",
runTime.constant(),
mesh,
IOobject::MUST_READ,
IOobject::NO_WRITE
)
);
dimensionedScalar mu
(
transportProperties.lookup("mu")
);

53
applications/solvers/compressible/sonicFoam/sonicFoam.C
View File

@ -32,7 +32,7 @@ Description
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "basicThermo.H"
#include "basicPsiThermo.H"
#include "turbulenceModel.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -58,64 +58,21 @@ int main(int argc, char *argv[])
#include "rhoEqn.H"
fvVectorMatrix UEqn
(
fvm::ddt(rho, U)
+ fvm::div(phi, U)
+ turbulence->divDevRhoReff(U)
);
#include "UEqn.H"
solve(UEqn == -fvc::grad(p));
#include "hEqn.H"
#include "eEqn.H"
// --- PISO loop
for (int corr=0; corr<nCorr; corr++)
{
volScalarField rUA = 1.0/UEqn.A();
U = rUA*UEqn.H();
surfaceScalarField phid
(
"phid",
fvc::interpolate(thermo->psi())
*(
(fvc::interpolate(U) & mesh.Sf())
+ fvc::ddtPhiCorr(rUA, rho, U, phi)
)
);
for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
{
fvScalarMatrix pEqn
(
fvm::ddt(psi, p)
+ fvm::div(phid, p)
- fvm::laplacian(rho*rUA, p)
);
pEqn.solve();
if (nonOrth == nNonOrthCorr)
{
phi = pEqn.flux();
}
}
#include "compressibleContinuityErrs.H"
U -= rUA*fvc::grad(p);
U.correctBoundaryConditions();
#include "pEqn.H"
}
DpDt =
fvc::DDt(surfaceScalarField("phiU", phi/fvc::interpolate(rho)), p);
turbulence->correct();
rho = psi*p;
rho = thermo.rho();
runTime.write();

2
applications/solvers/compressible/sonicLiquidFoam/createFields.H
View File

@ -41,4 +41,4 @@
);
# include "compressibleCreatePhi.H"
#include "compressibleCreatePhi.H"

29
applications/solvers/compressible/sonicLiquidFoam/sonicLiquidFoam.C
View File

@ -37,16 +37,15 @@ Description
int main(int argc, char *argv[])
{
#include "setRootCase.H"
#include "createTime.H"
#include "createMesh.H"
#include "readThermodynamicProperties.H"
#include "readTransportProperties.H"
#include "createFields.H"
#include "initContinuityErrs.H"
# include "setRootCase.H"
# include "createTime.H"
# include "createMesh.H"
# include "readThermodynamicProperties.H"
# include "readTransportProperties.H"
# include "createFields.H"
# include "initContinuityErrs.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
@ -54,10 +53,10 @@ int main(int argc, char *argv[])
{
Info<< "Time = " << runTime.timeName() << nl << endl;
# include "readPISOControls.H"
# include "compressibleCourantNo.H"
#include "readPISOControls.H"
#include "compressibleCourantNo.H"
# include "rhoEqn.H"
#include "rhoEqn.H"
fvVectorMatrix UEqn
(
@ -79,8 +78,8 @@ int main(int argc, char *argv[])
surfaceScalarField phid
(
"phid",
psi*
(
psi
*(
(fvc::interpolate(U) & mesh.Sf())
+ fvc::ddtPhiCorr(rUA, rho, U, phi)
)
@ -100,7 +99,7 @@ int main(int argc, char *argv[])
phi += pEqn.flux();
# include "compressibleContinuityErrs.H"
#include "compressibleContinuityErrs.H"
U -= rUA*fvc::grad(p);
U.correctBoundaryConditions();

1
applications/solvers/heatTransfer/buoyantBoussinesqPisoFoam/Make/options
View File

@ -1,4 +1,5 @@
EXE_INC = \
-I../buoyantBoussinesqSimpleFoam \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/turbulenceModels \
-I$(LIB_SRC)/turbulenceModels/incompressible/RAS/lnInclude \

4
applications/solvers/heatTransfer/buoyantBoussinesqPisoFoam/TEqn.H
View File

@ -2,7 +2,7 @@
volScalarField kappaEff
(
"kappaEff",
turbulence->nu() + turbulence->nut()/Prt
turbulence->nu()/Pr + turbulence->nut()/Prt
);
fvScalarMatrix TEqn
@ -15,4 +15,6 @@
TEqn.relax();
TEqn.solve();
rhok = 1.0 - beta*(T - TRef);
}

27
applications/solvers/heatTransfer/buoyantBoussinesqPisoFoam/UEqn.H
View File

@ -1,23 +1,26 @@
// Solve the momentum equation
tmp<fvVectorMatrix> UEqn
fvVectorMatrix UEqn
(
fvm::ddt(U)
+ fvm::div(phi, U)
+ turbulence->divDevReff(U)
);
UEqn().relax();
UEqn.relax();
solve
(
UEqn()
==
-fvc::reconstruct
if (momentumPredictor)
{
solve
(
UEqn
==
fvc::reconstruct
(
fvc::snGrad(pd)
- betaghf*fvc::snGrad(T)
) * mesh.magSf()
)
);
(
fvc::interpolate(rhok)*(g & mesh.Sf())
- fvc::snGrad(p)*mesh.magSf()
)
)
);
}

40
applications/solvers/heatTransfer/buoyantBoussinesqPisoFoam/buoyantBoussinesqPisoFoam.C
View File

@ -54,18 +54,17 @@ Description
int main(int argc, char *argv[])
{
#include "setRootCase.H"
#include "createTime.H"
#include "createMesh.H"
#include "readEnvironmentalProperties.H"
#include "createFields.H"
#include "initContinuityErrs.H"
#include "readTimeControls.H"
#include "CourantNo.H"
#include "setInitialDeltaT.H"
# include "setRootCase.H"
# include "createTime.H"
# include "createMesh.H"
# include "readEnvironmentalProperties.H"
# include "createFields.H"
# include "initContinuityErrs.H"
# include "readTimeControls.H"
# include "CourantNo.H"
# include "setInitialDeltaT.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
@ -73,26 +72,23 @@ int main(int argc, char *argv[])
{
Info<< "Time = " << runTime.timeName() << nl << endl;
# include "readTimeControls.H"
# include "readPISOControls.H"
# include "CourantNo.H"
# include "setDeltaT.H"
#include "readTimeControls.H"
#include "readPISOControls.H"
#include "CourantNo.H"
#include "setDeltaT.H"
# include "UEqn.H"
#include "UEqn.H"
#include "TEqn.H"
// --- PISO loop
for (int corr=0; corr<nCorr; corr++)
{
# include "TEqn.H"
# include "pdEqn.H"
#include "pEqn.H"
}
turbulence->correct();
if (runTime.write())
{
# include "writeAdditionalFields.H"
}
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"

35
applications/solvers/heatTransfer/buoyantBoussinesqPisoFoam/createFields.H
View File

@ -14,13 +14,12 @@
mesh
);
// kinematic pd
Info<< "Reading field pd\n" << endl;
volScalarField pd
Info<< "Reading field p\n" << endl;
volScalarField p
(
IOobject
(
"pd",
"p",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
@ -53,15 +52,25 @@
incompressible::RASModel::New(U, phi, laminarTransport)
);
Info<< "Calculating field beta*(g.h)\n" << endl;
surfaceScalarField betaghf("betagh", beta*(g & mesh.Cf()));
label pdRefCell = 0;
scalar pdRefValue = 0.0;
label pRefCell = 0;
scalar pRefValue = 0.0;
setRefCell
(
pd,
mesh.solutionDict().subDict("SIMPLE"),
pdRefCell,
pdRefValue
p,
mesh.solutionDict().subDict("PISO"),
pRefCell,
pRefValue
);
// Kinematic density for buoyancy force
volScalarField rhok
(
IOobject
(
"rhok",
runTime.timeName(),
mesh
),
1.0 - beta*(T - TRef)
);

19
applications/solvers/heatTransfer/buoyantBoussinesqPisoFoam/pdEqn.H → applications/solvers/heatTransfer/buoyantBoussinesqPisoFoam/pEqn.H
View File

@ -1,9 +1,8 @@
{
volScalarField rUA("rUA", 1.0/UEqn().A());
volScalarField rUA("rUA", 1.0/UEqn.A());
surfaceScalarField rUAf("(1|A(U))", fvc::interpolate(rUA));
U = rUA*UEqn().H();
UEqn.clear();
U = rUA*UEqn.H();
surfaceScalarField phiU
(
@ -11,31 +10,31 @@
+ fvc::ddtPhiCorr(rUA, U, phi)
);
phi = phiU + betaghf*fvc::snGrad(T)*rUAf*mesh.magSf();
phi = phiU + rUAf*fvc::interpolate(rhok)*(g & mesh.Sf());
for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
{
fvScalarMatrix pdEqn
fvScalarMatrix pEqn
(
fvm::laplacian(rUAf, pd) == fvc::div(phi)
fvm::laplacian(rUAf, p) == fvc::div(phi)
);
if (corr == nCorr-1 && nonOrth == nNonOrthCorr)
{
pdEqn.solve(mesh.solver(pd.name() + "Final"));
pEqn.solve(mesh.solver(p.name() + "Final"));
}
else
{
pdEqn.solve(mesh.solver(pd.name()));
pEqn.solve(mesh.solver(p.name()));
}
if (nonOrth == nNonOrthCorr)
{
phi += pdEqn.flux();
phi -= pEqn.flux();
}
}
U -= rUA*fvc::reconstruct((phi - phiU)/rUAf);
U += rUA*fvc::reconstruct((phi - phiU)/rUAf);
U.correctBoundaryConditions();
#include "continuityErrs.H"

13
applications/solvers/heatTransfer/buoyantBoussinesqPisoFoam/readTransportProperties.H
View File

@ -1,13 +0,0 @@
singlePhaseTransportModel laminarTransport(U, phi);
// thermal expansion coefficient [1/K]
dimensionedScalar beta(laminarTransport.lookup("beta"));
// reference temperature [K]
dimensionedScalar TRef(laminarTransport.lookup("TRef"));
// reference kinematic pressure [m2/s2]
dimensionedScalar pRef(laminarTransport.lookup("pRef"));
// turbulent Prandtl number
dimensionedScalar Prt(laminarTransport.lookup("Prt"));

29
applications/solvers/heatTransfer/buoyantBoussinesqPisoFoam/writeAdditionalFields.H
View File

@ -1,29 +0,0 @@
{
volScalarField rhoEff
(
IOobject
(
"rhoEff",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
1.0 - beta*(T - TRef)
);
rhoEff.write();
volScalarField p
(
IOobject
(
"p",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
pd + rhoEff*(g & mesh.C()) + pRef
);
p.write();
}

4
applications/solvers/heatTransfer/buoyantBoussinesqSimpleFoam/TEqn.H
View File

@ -2,7 +2,7 @@
volScalarField kappaEff
(
"kappaEff",
turbulence->nu() + turbulence->nut()/Prt
turbulence->nu()/Pr + turbulence->nut()/Prt
);
fvScalarMatrix TEqn
@ -16,4 +16,6 @@
eqnResidual = TEqn.solve().initialResidual();
maxResidual = max(eqnResidual, maxResidual);
rhok = 1.0 - beta*(T - TRef);
}

8
applications/solvers/heatTransfer/buoyantBoussinesqSimpleFoam/UEqn.H
View File

@ -13,12 +13,12 @@
(
UEqn()
==
-fvc::reconstruct
fvc::reconstruct
(
(
fvc::snGrad(pd)
- betaghf*fvc::snGrad(T)
) * mesh.magSf()
fvc::interpolate(rhok)*(g & mesh.Sf())
- fvc::snGrad(p)*mesh.magSf()
)
)
).initialResidual();

34
applications/solvers/heatTransfer/buoyantBoussinesqSimpleFoam/buoyantBoussinesqSimpleFoam.C
View File

@ -54,15 +54,14 @@ Description
int main(int argc, char *argv[])
{
#include "setRootCase.H"
#include "createTime.H"
#include "createMesh.H"
#include "readEnvironmentalProperties.H"
#include "createFields.H"
#include "initContinuityErrs.H"
# include "setRootCase.H"
# include "createTime.H"
# include "createMesh.H"
# include "readEnvironmentalProperties.H"
# include "createFields.H"
# include "initContinuityErrs.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
@ -70,30 +69,27 @@ int main(int argc, char *argv[])
{
Info<< "Time = " << runTime.timeName() << nl << endl;
# include "readSIMPLEControls.H"
# include "initConvergenceCheck.H"
#include "readSIMPLEControls.H"
#include "initConvergenceCheck.H"
pd.storePrevIter();
p.storePrevIter();
// Pressure-velocity SIMPLE corrector
{
# include "UEqn.H"
# include "TEqn.H"
# include "pdEqn.H"
#include "UEqn.H"
#include "TEqn.H"
#include "pEqn.H"
}
turbulence->correct();
if (runTime.write())
{
# include "writeAdditionalFields.H"
}
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
# include "convergenceCheck.H"
#include "convergenceCheck.H"
}
Info<< "End\n" << endl;

30
applications/solvers/heatTransfer/buoyantBoussinesqSimpleFoam/createFields.H
View File

@ -14,13 +14,12 @@
mesh
);
// kinematic pd
Info<< "Reading field pd\n" << endl;
volScalarField pd
Info<< "Reading field p\n" << endl;
volScalarField p
(
IOobject
(
"pd",
"p",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
@ -56,12 +55,25 @@
Info<< "Calculating field beta*(g.h)\n" << endl;
surfaceScalarField betaghf("betagh", beta*(g & mesh.Cf()));
label pdRefCell = 0;
scalar pdRefValue = 0.0;
label pRefCell = 0;
scalar pRefValue = 0.0;
setRefCell
(
pd,
p,
mesh.solutionDict().subDict("SIMPLE"),
pdRefCell,
pdRefValue
pRefCell,
pRefValue
);
// Kinematic density for buoyancy force
volScalarField rhok
(
IOobject
(
"rhok",
runTime.timeName(),
mesh
),
1.0 - beta*(T - TRef)
);

23
applications/solvers/heatTransfer/buoyantBoussinesqSimpleFoam/pdEqn.H → applications/solvers/heatTransfer/buoyantBoussinesqSimpleFoam/pEqn.H
View File

@ -6,41 +6,42 @@
UEqn.clear();
phi = fvc::interpolate(U) & mesh.Sf();
adjustPhi(phi, U, pd);
surfaceScalarField buoyancyPhi = -betaghf*fvc::snGrad(T)*rUAf*mesh.magSf();
phi -= buoyancyPhi;
adjustPhi(phi, U, p);
surfaceScalarField buoyancyPhi =
rUAf*fvc::interpolate(rhok)*(g & mesh.Sf());
phi += buoyancyPhi;
for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
{
fvScalarMatrix pdEqn
fvScalarMatrix pEqn
(
fvm::laplacian(rUAf, pd) == fvc::div(phi)
fvm::laplacian(rUAf, p) == fvc::div(phi)
);
pdEqn.setReference(pdRefCell, pdRefValue);
pEqn.setReference(pRefCell, pRefValue);
// retain the residual from the first iteration
if (nonOrth == 0)
{
eqnResidual = pdEqn.solve().initialResidual();
eqnResidual = pEqn.solve().initialResidual();
maxResidual = max(eqnResidual, maxResidual);
}
else
{
pdEqn.solve();
pEqn.solve();
}
if (nonOrth == nNonOrthCorr)
{
// Calculate the conservative fluxes
phi -= pdEqn.flux();
phi -= pEqn.flux();
// Explicitly relax pressure for momentum corrector
pd.relax();
p.relax();
// Correct the momentum source with the pressure gradient flux
// calculated from the relaxed pressure
U -= rUA*fvc::reconstruct((buoyancyPhi + pdEqn.flux())/rUAf);
U += rUA*fvc::reconstruct((buoyancyPhi - pEqn.flux())/rUAf);
U.correctBoundaryConditions();
}
}

10
applications/solvers/heatTransfer/buoyantBoussinesqSimpleFoam/readTransportProperties.H
View File

@ -1,13 +1,13 @@
singlePhaseTransportModel laminarTransport(U, phi);
// thermal expansion coefficient [1/K]
// Thermal expansion coefficient [1/K]
dimensionedScalar beta(laminarTransport.lookup("beta"));
// reference temperature [K]
// Reference temperature [K]
dimensionedScalar TRef(laminarTransport.lookup("TRef"));
// reference kinematic pressure [m2/s2]
dimensionedScalar pRef(laminarTransport.lookup("pRef"));
// Laminar Prandtl number
dimensionedScalar Pr(laminarTransport.lookup("Pr"));
// turbulent Prandtl number
// Turbulent Prandtl number
dimensionedScalar Prt(laminarTransport.lookup("Prt"));

29
applications/solvers/heatTransfer/buoyantBoussinesqSimpleFoam/writeAdditionalFields.H
View File

@ -1,29 +0,0 @@
{
volScalarField rhoEff
(
IOobject
(
"rhoEff",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
1.0 - beta*(T - TRef)
);
rhoEff.write();
volScalarField p
(
IOobject
(
"p",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
pd + rhoEff*(g & mesh.C()) + pRef
);
p.write();
}

8
applications/solvers/heatTransfer/buoyantPisoFoam/buoyantPisoFoam.C
View File

@ -27,13 +27,13 @@ Application
Description
Transient Solver for buoyant, turbulent flow of compressible fluids for
ventilation and heat-transfer. Turbulence is modelled using a run-time
selectable compressible RAS model.
ventilation and heat-transfer. Turbulence is modelled using a run-time
selectable compressible RAS or LES model.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "basicThermo.H"
#include "basicRhoThermo.H"
#include "turbulenceModel.H"
#include "fixedGradientFvPatchFields.H"
@ -81,6 +81,8 @@ int main(int argc, char *argv[])
turbulence->correct();
rho = thermo.rho();
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"

21
applications/solvers/heatTransfer/buoyantPisoFoam/createFields.H
View File

@ -1,9 +1,10 @@
Info<< "Reading thermophysical properties\n" << endl;
autoPtr<basicThermo> thermo
autoPtr<basicRhoThermo> pThermo
(
basicThermo::New(mesh)
basicRhoThermo::New(mesh)
);
basicRhoThermo& thermo = pThermo();
volScalarField rho
(
@ -15,12 +16,12 @@
IOobject::NO_READ,
IOobject::NO_WRITE
),
thermo->rho()
thermo.rho()
);
volScalarField& p = thermo->p();
volScalarField& h = thermo->h();
const volScalarField& psi = thermo->psi();
volScalarField& p = thermo.p();
volScalarField& h = thermo.h();
const volScalarField& psi = thermo.psi();
Info<< "Reading field U\n" << endl;
@ -37,7 +38,7 @@
mesh
);
# include "compressibleCreatePhi.H"
#include "compressibleCreatePhi.H"
Info<< "Creating turbulence model\n" << endl;
@ -48,7 +49,7 @@
rho,
U,
phi,
thermo()
thermo
)
);
@ -59,6 +60,8 @@
fvc::DDt(surfaceScalarField("phiU", phi/fvc::interpolate(rho)), p)
);
thermo->correct();
thermo.correct();
dimensionedScalar initialMass = fvc::domainIntegrate(rho);
dimensionedScalar totalVolume = sum(mesh.V());

2
applications/solvers/heatTransfer/buoyantPisoFoam/hEqn.H
View File

@ -11,5 +11,5 @@
hEqn.relax();
hEqn.solve();
thermo->correct();
thermo.correct();
}

19
applications/solvers/heatTransfer/buoyantPisoFoam/pEqn.H
View File

@ -1,7 +1,11 @@
{
bool closedVolume = p.needReference();
rho = thermo->rho();
rho = thermo.rho();
// Thermodynamic density needs to be updated by psi*d(p) after the
// pressure solution - done in 2 parts. Part 1:
thermo.rho() -= psi*p;
volScalarField rUA = 1.0/UEqn.A();
surfaceScalarField rhorUAf("(rho*(1|A(U)))", fvc::interpolate(rho*rUA));
@ -23,7 +27,7 @@
{
fvScalarMatrix pEqn
(
fvm::ddt(psi, p)
fvc::ddt(rho) + psi*correction(fvm::ddt(p))
+ fvc::div(phi)
- fvm::laplacian(rhorUAf, p)
);
@ -43,6 +47,9 @@
}
}
// Second part of thermodynamic density update
thermo.rho() += psi*p;
U += rUA*fvc::reconstruct((phi - phiU)/rhorUAf);
U.correctBoundaryConditions();
@ -55,8 +62,10 @@
// to obey overall mass continuity
if (closedVolume)
{
p += (initialMass - fvc::domainIntegrate(thermo->psi()*p))
/fvc::domainIntegrate(thermo->psi());
rho = thermo->rho();
p +=
(initialMass - fvc::domainIntegrate(psi*p))
/fvc::domainIntegrate(psi);
thermo.rho() = psi*p;
rho += (initialMass - fvc::domainIntegrate(rho))/totalVolume;
}
}

2
applications/solvers/heatTransfer/buoyantSimpleFoam/buoyantSimpleFoam.C
View File

@ -31,7 +31,7 @@ Description
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "basicThermo.H"
#include "basicPsiThermo.H"
#include "RASModel.H"
#include "fixedGradientFvPatchFields.H"

18
applications/solvers/heatTransfer/buoyantSimpleFoam/createFields.H
View File

@ -1,9 +1,10 @@
Info<< "Reading thermophysical properties\n" << endl;
autoPtr<basicThermo> thermo
autoPtr<basicPsiThermo> pThermo
(
basicThermo::New(mesh)
basicPsiThermo::New(mesh)
);
basicPsiThermo& thermo = pThermo();
volScalarField rho
(
@ -15,11 +16,12 @@
IOobject::NO_READ,
IOobject::NO_WRITE
),
thermo->rho()
thermo.rho()
);
volScalarField& p = thermo->p();
volScalarField& h = thermo->h();
volScalarField& p = thermo.p();
volScalarField& h = thermo.h();
const volScalarField& psi = thermo.psi();
Info<< "Reading field U\n" << endl;
@ -36,7 +38,7 @@
mesh
);
# include "compressibleCreatePhi.H"
#include "compressibleCreatePhi.H"
Info<< "Creating turbulence model\n" << endl;
@ -47,11 +49,11 @@
rho,
U,
phi,
thermo()
thermo
)
);
thermo->correct();
thermo.correct();
label pRefCell = 0;
scalar pRefValue = 0.0;

2
applications/solvers/heatTransfer/buoyantSimpleFoam/hEqn.H
View File

@ -14,5 +14,5 @@
eqnResidual = hEqn.solve().initialResidual();
maxResidual = max(eqnResidual, maxResidual);
thermo->correct();
thermo.correct();
}

2
applications/solvers/heatTransfer/buoyantSimpleFoam/initConvergenceCheck.H
View File

@ -1,4 +1,4 @@
// initialize values for convergence checks
// initialize values for convergence checks
scalar eqnResidual = 1, maxResidual = 0;
scalar convergenceCriterion = 0;

8
applications/solvers/heatTransfer/buoyantSimpleFoam/pEqn.H
View File

@ -1,5 +1,5 @@
{
rho = thermo->rho();
rho = thermo.rho();
volScalarField rUA = 1.0/UEqn().A();
surfaceScalarField rhorUAf("(rho*(1|A(U)))", fvc::interpolate(rho*rUA));
@ -39,8 +39,8 @@
// to obey overall mass continuity
if (closedVolume)
{
p += (initialMass - fvc::domainIntegrate(thermo->psi()*p))
/fvc::domainIntegrate(thermo->psi());
p += (initialMass - fvc::domainIntegrate(psi*p))
/fvc::domainIntegrate(psi);
}
// Calculate the conservative fluxes
@ -58,7 +58,7 @@
#include "continuityErrs.H"
rho = thermo->rho();
rho = thermo.rho();
rho.relax();
Info<< "rho max/min : " << max(rho).value() << " " << min(rho).value()
<< endl;

Some files were not shown because too many files have changed in this diff Show More