Commit 80393339 authored by Henry's avatar Henry
Browse files

Upgrade of compressible solvers merging in developments from OpenFOAM-1.7.x

parent 06f2882b
EXE_INC = \
-IBCs/lnInclude \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
-IBCs/lnInclude \
-I$(LIB_SRC)/sampling/lnInclude
-I$(LIB_SRC)/turbulenceModels/compressible/turbulenceModel \
-I$(LIB_SRC)/dynamicMesh/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude
EXE_LIBS = \
-lfiniteVolume \
-lbasicThermophysicalModels \
-lspecie \
-lrhoCentralFoam
-lrhoCentralFoam \
-lcompressibleTurbulenceModel \
-lcompressibleRASModels \
-lcompressibleLESModels \
-lmeshTools
......@@ -96,3 +96,18 @@ surfaceScalarField neg
mesh,
dimensionedScalar("neg", dimless, -1.0)
);
surfaceScalarField phi("phi", mesh.Sf() & fvc::interpolate(rhoU));
Info<< "Creating turbulence model\n" << endl;
autoPtr<compressible::turbulenceModel> turbulence
(
compressible::turbulenceModel::New
(
rho,
U,
phi,
thermo
)
);
EXE_INC = \
-I.. \
-I../BCs/lnInclude \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
-I../BCs/lnInclude \
-I$(LIB_SRC)/sampling/lnInclude \
-I$(LIB_SRC)/turbulenceModels/compressible/turbulenceModel \
-I$(LIB_SRC)/dynamicMesh/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude
......@@ -13,5 +13,8 @@ EXE_LIBS = \
-lbasicThermophysicalModels \
-lspecie \
-lrhoCentralFoam \
-lcompressibleTurbulenceModel \
-lcompressibleRASModels \
-lcompressibleLESModels \
-ldynamicMesh \
-lmeshTools
......@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2009-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
......@@ -22,7 +22,7 @@ License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Application
rhoCentralFoam
rhoCentralDyMFoam
Description
Density-based compressible flow solver based on central-upwind schemes of
......@@ -32,6 +32,7 @@ Description
#include "fvCFD.H"
#include "basicPsiThermo.H"
#include "turbulenceModel.H"
#include "zeroGradientFvPatchFields.H"
#include "fixedRhoFvPatchScalarField.H"
#include "motionSolver.H"
......@@ -185,12 +186,11 @@ int main(int argc, char *argv[])
+ aSf*p_pos - aSf*p_neg
);
volTensorField tauMC("tauMC", mu*dev2(Foam::T(fvc::grad(U))));
volScalarField muEff(turbulence->muEff());
volTensorField tauMC("tauMC", muEff*dev2(Foam::T(fvc::grad(U))));
// --- Solve density
Info<< max(rho) << " " << min(rho) << endl;
solve(fvm::ddt(rho) + fvc::div(phi));
Info<< max(rho) << " " << min(rho) << endl;
// --- Solve momentum
solve(fvm::ddt(rhoU) + fvc::div(phiUp));
......@@ -206,7 +206,7 @@ int main(int argc, char *argv[])
solve
(
fvm::ddt(rho, U) - fvc::ddt(rho, U)
- fvm::laplacian(mu, U)
- fvm::laplacian(muEff, U)
- fvc::div(tauMC)
);
rhoU = rho*U;
......@@ -216,7 +216,7 @@ int main(int argc, char *argv[])
surfaceScalarField sigmaDotU
(
(
fvc::interpolate(mu)*mesh.magSf()*fvc::snGrad(U)
fvc::interpolate(muEff)*mesh.magSf()*fvc::snGrad(U)
+ (mesh.Sf() & fvc::interpolate(tauMC))
)
& (a_pos*U_pos + a_neg*U_neg)
......@@ -240,12 +240,12 @@ int main(int argc, char *argv[])
if (!inviscid)
{
volScalarField k("k", thermo.Cp()*mu/Pr);
volScalarField k("k", thermo.Cp()*muEff/Pr);
solve
(
fvm::ddt(rho, e) - fvc::ddt(rho, e)
- fvm::laplacian(thermo.alpha(), e)
+ fvc::laplacian(thermo.alpha(), e)
- fvm::laplacian(turbulence->alphaEff(), e)
+ fvc::laplacian(turbulence->alpha(), e)
- fvc::laplacian(k, T)
);
thermo.correct();
......@@ -258,6 +258,8 @@ int main(int argc, char *argv[])
p.correctBoundaryConditions();
rho.boundaryField() = psi.boundaryField()*p.boundaryField();
turbulence->correct();
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
......
......@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2009-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
......@@ -32,6 +32,7 @@ Description
#include "fvCFD.H"
#include "basicPsiThermo.H"
#include "turbulenceModel.H"
#include "zeroGradientFvPatchFields.H"
#include "fixedRhoFvPatchScalarField.H"
......@@ -51,7 +52,7 @@ int main(int argc, char *argv[])
#include "readFluxScheme.H"
dimensionedScalar v_zero("v_zero",dimVolume/dimTime, 0.0);
dimensionedScalar v_zero("v_zero", dimVolume/dimTime, 0.0);
Info<< "\nStarting time loop\n" << endl;
......@@ -128,14 +129,6 @@ int main(int argc, char *argv[])
surfaceScalarField amaxSf("amaxSf", max(mag(am), mag(ap)));
#include "compressibleCourantNo.H"
#include "readTimeControls.H"
#include "setDeltaT.H"
runTime++;
Info<< "Time = " << runTime.timeName() << nl << endl;
surfaceScalarField aSf(am*a_pos);
if (fluxScheme == "Tadmor")
......@@ -152,6 +145,18 @@ int main(int argc, char *argv[])
surfaceScalarField aphiv_pos(phiv_pos - aSf);
surfaceScalarField aphiv_neg(phiv_neg + aSf);
// Reuse amaxSf for the maximum positive and negative fluxes
// estimated by the central scheme
amaxSf = max(mag(aphiv_pos), mag(aphiv_neg));
#include "compressibleCourantNo.H"
#include "readTimeControls.H"
#include "setDeltaT.H"
runTime++;
Info<< "Time = " << runTime.timeName() << nl << endl;
surfaceScalarField phi("phi", aphiv_pos*rho_pos + aphiv_neg*rho_neg);
surfaceVectorField phiUp
......@@ -167,7 +172,8 @@ int main(int argc, char *argv[])
+ aSf*p_pos - aSf*p_neg
);
volTensorField tauMC("tauMC", mu*dev2(Foam::T(fvc::grad(U))));
volScalarField muEff(turbulence->muEff());
volTensorField tauMC("tauMC", muEff*dev2(Foam::T(fvc::grad(U))));
// --- Solve density
solve(fvm::ddt(rho) + fvc::div(phi));
......@@ -188,7 +194,7 @@ int main(int argc, char *argv[])
solve
(
fvm::ddt(rho, U) - fvc::ddt(rho, U)
- fvm::laplacian(mu, U)
- fvm::laplacian(muEff, U)
- fvc::div(tauMC)
);
rhoU = rho*U;
......@@ -198,7 +204,7 @@ int main(int argc, char *argv[])
surfaceScalarField sigmaDotU
(
(
fvc::interpolate(mu)*mesh.magSf()*fvc::snGrad(U)
fvc::interpolate(muEff)*mesh.magSf()*fvc::snGrad(U)
+ (mesh.Sf() & fvc::interpolate(tauMC))
)
& (a_pos*U_pos + a_neg*U_neg)
......@@ -222,12 +228,12 @@ int main(int argc, char *argv[])
if (!inviscid)
{
volScalarField k("k", thermo.Cp()*mu/Pr);
volScalarField k("k", thermo.Cp()*muEff/Pr);
solve
(
fvm::ddt(rho, e) - fvc::ddt(rho, e)
- fvm::laplacian(thermo.alpha(), e)
+ fvc::laplacian(thermo.alpha(), e)
- fvm::laplacian(turbulence->alphaEff(), e)
+ fvc::laplacian(turbulence->alpha(), e)
- fvc::laplacian(k, T)
);
thermo.correct();
......@@ -240,6 +246,8 @@ int main(int argc, char *argv[])
p.correctBoundaryConditions();
rho.boundaryField() = psi.boundaryField()*p.boundaryField();
turbulence->correct();
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
......
#!/bin/sh
cd ${0%/*} || exit 1 # run from this directory
set -x
wmake
wmake rhoPorousMRFPimpleFoam
# ----------------------------------------------------------------- end-of-file
......@@ -39,6 +39,16 @@
#include "compressibleCreatePhi.H"
dimensionedScalar rhoMax
(
mesh.solutionDict().subDict("PIMPLE").lookup("rhoMax")
);
dimensionedScalar rhoMin
(
mesh.solutionDict().subDict("PIMPLE").lookup("rhoMin")
);
Info<< "Creating turbulence model\n" << endl;
autoPtr<compressible::turbulenceModel> turbulence
(
......
rho = thermo.rho();
rho = max(rho, rhoMin);
rho = min(rho, rhoMax);
rho.relax();
U = rAU*UEqn().H();
......@@ -98,6 +101,9 @@ p.relax();
// Recalculate density from the relaxed pressure
rho = thermo.rho();
rho = max(rho, rhoMin);
rho = min(rho, rhoMax);
rho.relax();
Info<< "rho max/min : " << max(rho).value()
<< " " << min(rho).value() << endl;
......
......@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2004-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
......@@ -75,6 +75,7 @@ int main(int argc, char *argv[])
if (nOuterCorr != 1)
{
p.storePrevIter();
rho.storePrevIter();
}
#include "UEqn.H"
......
EXE_INC = \
-I../rhoPimpleFoam \
-I.. \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/turbulenceModels/compressible/turbulenceModel \
-I$(LIB_SRC)/finiteVolume/cfdTools \
......
MRFZones mrfZones(mesh);
mrfZones.correctBoundaryVelocity(U);
porousZones pZones(mesh);
Switch pressureImplicitPorosity(false);
rho = thermo.rho();
rho = max(rho, rhoMin);
rho = min(rho, rhoMax);
rho.relax();
volScalarField rAU(1.0/UEqn().A());
U = rAU*UEqn().H();
......@@ -96,30 +99,17 @@ else
#include "rhoEqn.H"
#include "compressibleContinuityErrs.H"
//if (oCorr != nOuterCorr-1)
{
// Explicitly relax pressure for momentum corrector
p.relax();
// Explicitly relax pressure for momentum corrector
p.relax();
rho = thermo.rho();
rho = max(rho, rhoMin);
rho = min(rho, rhoMax);
rho.relax();
Info<< "rho max/min : " << max(rho).value()
<< " " << min(rho).value() << endl;
}
rho = thermo.rho();
rho = max(rho, rhoMin);
rho = min(rho, rhoMax);
rho.relax();
Info<< "rho max/min : " << max(rho).value()
<< " " << min(rho).value() << endl;
U -= rAU*fvc::grad(p);
U.correctBoundaryConditions();
DpDt = fvc::DDt(surfaceScalarField("phiU", phi/fvc::interpolate(rho)), p);
// For closed-volume cases adjust the pressure and density levels
// to obey overall mass continuity
/*
if (closedVolume)
{
p += (initialMass - fvc::domainIntegrate(psi*p))
/fvc::domainIntegrate(psi);
}
*/
......@@ -48,6 +48,7 @@ int main(int argc, char *argv[])
#include "createTime.H"
#include "createMesh.H"
#include "createFields.H"
#include "createZones.H"
#include "initContinuityErrs.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
......
Info<< "Reading thermophysical properties\n" << endl;
autoPtr<basicPsiThermo> pThermo
(
basicPsiThermo::New(mesh)
);
basicPsiThermo& thermo = pThermo();
volScalarField& p = thermo.p();
volScalarField& h = thermo.h();
const volScalarField& psi = thermo.psi();
volScalarField rho
(
IOobject
(
"rho",
runTime.timeName(),
mesh,
IOobject::READ_IF_PRESENT,
IOobject::AUTO_WRITE
),
thermo.rho()
);
Info<< "Reading field U\n" << endl;
volVectorField U
(
IOobject
(
"U",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
#include "compressibleCreatePhi.H"
dimensionedScalar rhoMax
(
mesh.solutionDict().subDict("PIMPLE").lookup("rhoMax")
);
dimensionedScalar rhoMin
(
mesh.solutionDict().subDict("PIMPLE").lookup("rhoMin")
);
Info<< "Creating turbulence model\n" << endl;
autoPtr<compressible::turbulenceModel> turbulence
(
compressible::turbulenceModel::New
(
rho,
U,
phi,
thermo
)
);
//dimensionedScalar initialMass = fvc::domainIntegrate(rho);
Info<< "Creating field DpDt\n" << endl;
volScalarField DpDt
(
fvc::DDt(surfaceScalarField("phiU", phi/fvc::interpolate(rho)), p)
);
MRFZones mrfZones(mesh);
mrfZones.correctBoundaryVelocity(U);
porousZones pZones(mesh);
Switch pressureImplicitPorosity(false);
#!/bin/sh
cd ${0%/*} || exit 1 # run from this directory
set -x
wmake
wmake rhoPorousSimpleFoam
# ----------------------------------------------------------------- end-of-file
EXE_INC = \
-I.. \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/thermalPorousZone/lnInclude \
-I$(LIB_SRC)/turbulenceModels \
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment