mirror of
https://develop.openfoam.com/Development/openfoam.git
synced 2025-11-28 03:28:01 +00:00
ENH: Adding overBuoyantPimpleFoam and tutorial
This commit is contained in:
sergio
@ -0,0 +1,3 @@
|
||||
overBuoyantPimpleDyMFoam.C
|
||||
|
||||
EXE = $(FOAM_APPBIN)/overBuoyantPimpleDyMFoam
|
||||
@ -0,0 +1,30 @@
|
||||
EXE_INC = \
|
||||
-I$(FOAM_SOLVERS)/heatTransfer/buoyantPimpleFoam \
|
||||
-I$(FOAM_SOLVERS)/compressible/rhoPimpleFoam/overRhoPimpleDyMFoam \
|
||||
-I$(LIB_SRC)/finiteVolume/lnInclude \
|
||||
-I$(LIB_SRC)/sampling/lnInclude \
|
||||
-I$(LIB_SRC)/dynamicFvMesh/lnInclude \
|
||||
-I$(LIB_SRC)/dynamicMesh/lnInclude \
|
||||
-I$(LIB_SRC)/overset/lnInclude \
|
||||
-I$(LIB_SRC)/meshTools/lnInclude \
|
||||
-I$(LIB_SRC)/transportModels/compressible/lnInclude \
|
||||
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
|
||||
-I$(LIB_SRC)/thermophysicalModels/radiation/lnInclude \
|
||||
-I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
|
||||
-I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude
|
||||
|
||||
EXE_LIBS = \
|
||||
-lfiniteVolume \
|
||||
-lsampling \
|
||||
-lmeshTools \
|
||||
-lfvOptions \
|
||||
-lcompressibleTransportModels \
|
||||
-lfluidThermophysicalModels \
|
||||
-lradiationModels \
|
||||
-lspecie \
|
||||
-lturbulenceModels \
|
||||
-lcompressibleTurbulenceModels \
|
||||
-latmosphericModels \
|
||||
-loverset \
|
||||
-ldynamicFvMesh \
|
||||
-ltopoChangerFvMesh
|
||||
@ -0,0 +1,91 @@
|
||||
if (mesh.changing())
|
||||
{
|
||||
volVectorField::Boundary& bfld = U.boundaryFieldRef();
|
||||
forAll(bfld, patchi)
|
||||
{
|
||||
if (bfld[patchi].fixesValue())
|
||||
{
|
||||
bfld[patchi].initEvaluate();
|
||||
}
|
||||
}
|
||||
|
||||
surfaceScalarField::Boundary& phiBfld = phi.boundaryFieldRef();
|
||||
forAll(bfld, patchi)
|
||||
{
|
||||
if (bfld[patchi].fixesValue())
|
||||
{
|
||||
bfld[patchi].evaluate();
|
||||
|
||||
phiBfld[patchi] =
|
||||
rho.boundaryField()[patchi]
|
||||
* (
|
||||
bfld[patchi]
|
||||
& mesh.Sf().boundaryField()[patchi]
|
||||
);
|
||||
}
|
||||
}
|
||||
}
|
||||
// Initialize BCs list for pcorr to zero-gradient
|
||||
wordList pcorrTypes
|
||||
(
|
||||
p.boundaryField().size(),
|
||||
zeroGradientFvPatchScalarField::typeName
|
||||
);
|
||||
|
||||
// Set BCs of pcorr to fixed-value for patches at which p is fixed
|
||||
forAll(p.boundaryField(), patchi)
|
||||
{
|
||||
if (p.boundaryField()[patchi].fixesValue())
|
||||
{
|
||||
pcorrTypes[patchi] = fixedValueFvPatchScalarField::typeName;
|
||||
}
|
||||
}
|
||||
|
||||
volScalarField pcorr
|
||||
(
|
||||
IOobject
|
||||
(
|
||||
"pcorr",
|
||||
runTime.timeName(),
|
||||
mesh,
|
||||
IOobject::NO_READ,
|
||||
IOobject::NO_WRITE
|
||||
),
|
||||
mesh,
|
||||
dimensionedScalar(p.dimensions(), Zero),
|
||||
pcorrTypes
|
||||
);
|
||||
|
||||
mesh.setFluxRequired(pcorr.name());
|
||||
|
||||
{
|
||||
dimensionedScalar rAUf("rAUf", dimTime, 1.0);
|
||||
|
||||
while (pimple.correctNonOrthogonal())
|
||||
{
|
||||
fvScalarMatrix pcorrEqn
|
||||
(
|
||||
fvm::ddt(psi, pcorr)
|
||||
+ fvc::div(phi)
|
||||
- fvm::laplacian(rAUf, pcorr)
|
||||
==
|
||||
divrhoU()
|
||||
);
|
||||
|
||||
//pcorrEqn.solve(mesh.solver(pcorr.select(pimple.finalInnerIter())));
|
||||
//Bypass virtual layer
|
||||
const dictionary& d = mesh.solver
|
||||
(
|
||||
pcorr.select
|
||||
(
|
||||
pimple.finalInnerIter()
|
||||
)
|
||||
);
|
||||
mesh.fvMesh::solve(pcorrEqn, d);
|
||||
|
||||
if (pimple.finalNonOrthogonalIter())
|
||||
{
|
||||
phi += pcorrEqn.flux();
|
||||
}
|
||||
}
|
||||
}
|
||||
@ -0,0 +1 @@
|
||||
const volScalarField& psi = thermo.psi();
|
||||
@ -0,0 +1,126 @@
|
||||
Info<< "Reading thermophysical properties\n" << endl;
|
||||
|
||||
autoPtr<rhoThermo> pThermo(rhoThermo::New(mesh));
|
||||
rhoThermo& thermo = pThermo();
|
||||
thermo.validate(args.executable(), "h", "e");
|
||||
|
||||
volScalarField rho
|
||||
(
|
||||
IOobject
|
||||
(
|
||||
"rho",
|
||||
runTime.timeName(),
|
||||
mesh,
|
||||
IOobject::NO_READ,
|
||||
IOobject::NO_WRITE
|
||||
),
|
||||
thermo.rho()
|
||||
);
|
||||
|
||||
volScalarField& p = thermo.p();
|
||||
|
||||
Info<< "Reading field U\n" << endl;
|
||||
volVectorField U
|
||||
(
|
||||
IOobject
|
||||
(
|
||||
"U",
|
||||
runTime.timeName(),
|
||||
mesh,
|
||||
IOobject::MUST_READ,
|
||||
IOobject::AUTO_WRITE
|
||||
),
|
||||
mesh
|
||||
);
|
||||
|
||||
#include "compressibleCreatePhi.H"
|
||||
|
||||
|
||||
Info<< "Creating turbulence model\n" << endl;
|
||||
autoPtr<compressible::turbulenceModel> turbulence
|
||||
(
|
||||
compressible::turbulenceModel::New
|
||||
(
|
||||
rho,
|
||||
U,
|
||||
phi,
|
||||
thermo
|
||||
)
|
||||
);
|
||||
|
||||
|
||||
#include "readGravitationalAcceleration.H"
|
||||
#include "readhRef.H"
|
||||
#include "gh.H"
|
||||
|
||||
|
||||
Info<< "Reading field p_rgh\n" << endl;
|
||||
volScalarField p_rgh
|
||||
(
|
||||
IOobject
|
||||
(
|
||||
"p_rgh",
|
||||
runTime.timeName(),
|
||||
mesh,
|
||||
IOobject::MUST_READ,
|
||||
IOobject::AUTO_WRITE
|
||||
),
|
||||
mesh
|
||||
);
|
||||
|
||||
// Force p_rgh to be consistent with p
|
||||
p_rgh = p - rho*gh;
|
||||
|
||||
mesh.setFluxRequired(p_rgh.name());
|
||||
|
||||
label pRefCell = 0;
|
||||
scalar pRefValue = 0.0;
|
||||
|
||||
if (p_rgh.needReference())
|
||||
{
|
||||
setRefCell
|
||||
(
|
||||
p,
|
||||
p_rgh,
|
||||
pimple.dict(),
|
||||
pRefCell,
|
||||
pRefValue
|
||||
);
|
||||
|
||||
p += dimensionedScalar
|
||||
(
|
||||
"p",
|
||||
p.dimensions(),
|
||||
pRefValue - getRefCellValue(p, pRefCell)
|
||||
);
|
||||
}
|
||||
|
||||
dimensionedScalar initialMass("initialMass", fvc::domainIntegrate(rho));
|
||||
|
||||
#include "createDpdt.H"
|
||||
|
||||
#include "createK.H"
|
||||
|
||||
#include "createMRF.H"
|
||||
#include "createRadiationModel.H"
|
||||
#include "createFvOptions.H"
|
||||
|
||||
//- Overset specific
|
||||
|
||||
// Add solver-specific interpolations
|
||||
{
|
||||
wordHashSet& nonInt =
|
||||
const_cast<wordHashSet&>(Stencil::New(mesh).nonInterpolatedFields());
|
||||
|
||||
nonInt.insert("HbyA");
|
||||
nonInt.insert("grad(p_rgh)");
|
||||
nonInt.insert("surfaceIntegrate(phi)");
|
||||
nonInt.insert("surfaceIntegrate(phiHbyA)");
|
||||
nonInt.insert("cellMask");
|
||||
nonInt.insert("cellDisplacement");
|
||||
nonInt.insert("interpolatedCells");
|
||||
nonInt.insert("cellInterpolationWeight");
|
||||
}
|
||||
|
||||
// Mask field for zeroing out contributions on hole cells
|
||||
#include "createCellMask.H"
|
||||
@ -0,0 +1,212 @@
|
||||
/*---------------------------------------------------------------------------*\
|
||||
========= |
|
||||
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
|
||||
\\ / O peration |
|
||||
\\ / A nd | Copyright (C) 2019 OpenCFD Ltd.
|
||||
\\/ M anipulation |
|
||||
-------------------------------------------------------------------------------
|
||||
|
||||
License
|
||||
This file is part of OpenFOAM.
|
||||
|
||||
OpenFOAM is free software: you can redistribute it and/or modify it
|
||||
under the terms of the GNU General Public License as published by
|
||||
the Free Software Foundation, either version 3 of the License, or
|
||||
(at your option) any later version.
|
||||
|
||||
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
|
||||
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
|
||||
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
|
||||
for more details.
|
||||
|
||||
You should have received a copy of the GNU General Public License
|
||||
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
|
||||
|
||||
Application
|
||||
overBuoyantPimpleDymFoam
|
||||
|
||||
Group
|
||||
grpHeatTransferSolvers
|
||||
|
||||
Description
|
||||
Transient solver for buoyant, turbulent flow of compressible fluids for
|
||||
ventilation and heat-transfer with overset feature
|
||||
|
||||
Turbulence is modelled using a run-time selectable compressible RAS or
|
||||
LES model.
|
||||
|
||||
\*---------------------------------------------------------------------------*/
|
||||
|
||||
#include "fvCFD.H"
|
||||
#include "dynamicFvMesh.H"
|
||||
#include "rhoThermo.H"
|
||||
#include "turbulentFluidThermoModel.H"
|
||||
#include "radiationModel.H"
|
||||
#include "fvOptions.H"
|
||||
#include "pimpleControl.H"
|
||||
#include "pressureControl.H"
|
||||
|
||||
#include "CorrectPhi.H"
|
||||
#include "cellCellStencilObject.H"
|
||||
#include "localMin.H"
|
||||
|
||||
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
|
||||
|
||||
int main(int argc, char *argv[])
|
||||
{
|
||||
argList::addNote
|
||||
(
|
||||
"Transient solver for buoyant, turbulent fluid flow"
|
||||
" of compressible fluids, including radiation."
|
||||
);
|
||||
|
||||
#include "postProcess.H"
|
||||
|
||||
#include "addCheckCaseOptions.H"
|
||||
#include "setRootCaseLists.H"
|
||||
#include "createTime.H"
|
||||
#include "createDynamicFvMesh.H"
|
||||
#include "createDyMControls.H"
|
||||
#include "createFields.H"
|
||||
#include "createFieldRefs.H"
|
||||
#include "initContinuityErrs.H"
|
||||
|
||||
#include "createRhoUfIfPresent.H"
|
||||
#include "createControls.H"
|
||||
|
||||
#include "compressibleCourantNo.H"
|
||||
#include "setInitialDeltaT.H"
|
||||
|
||||
turbulence->validate();
|
||||
|
||||
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
|
||||
|
||||
Info<< "\nStarting time loop\n" << endl;
|
||||
|
||||
while (runTime.run())
|
||||
{
|
||||
#include "readTimeControls.H"
|
||||
|
||||
#include "readControls.H"
|
||||
#include "readDyMControls.H"
|
||||
|
||||
#include "compressibleCourantNo.H"
|
||||
#include "setDeltaT.H"
|
||||
|
||||
// Store divrhoU from the previous mesh so that it can be mapped
|
||||
// and used in correctPhi to ensure the corrected phi has the
|
||||
// same divergence
|
||||
autoPtr<volScalarField> divrhoU;
|
||||
if (correctPhi)
|
||||
{
|
||||
divrhoU.reset
|
||||
(
|
||||
new volScalarField
|
||||
(
|
||||
"divrhoU",
|
||||
fvc::div(fvc::absolute(phi, rho, U))
|
||||
)
|
||||
);
|
||||
}
|
||||
|
||||
|
||||
++runTime;
|
||||
|
||||
Info<< "Time = " << runTime.timeName() << nl << endl;
|
||||
|
||||
// --- Pressure-velocity PIMPLE corrector loop
|
||||
while (pimple.loop())
|
||||
{
|
||||
if (pimple.firstIter() || moveMeshOuterCorrectors)
|
||||
{
|
||||
// Do any mesh changes
|
||||
mesh.update();
|
||||
|
||||
if (mesh.changing())
|
||||
{
|
||||
MRF.update();
|
||||
|
||||
#include "setCellMask.H"
|
||||
|
||||
const surfaceScalarField faceMaskOld
|
||||
(
|
||||
localMin<scalar>(mesh).interpolate(cellMask.oldTime())
|
||||
);
|
||||
|
||||
// Zero Uf on old faceMask (H-I)
|
||||
rhoUf() *= faceMaskOld;
|
||||
|
||||
//fvc::correctRhoUf(rhoUfint, rho, U, phi);
|
||||
surfaceVectorField rhoUfint(fvc::interpolate(rho*U));
|
||||
|
||||
// Update Uf and phi on new C-I faces
|
||||
rhoUf() += (1-faceMaskOld)*rhoUfint;
|
||||
|
||||
// Update Uf boundary
|
||||
forAll(rhoUf().boundaryField(), patchI)
|
||||
{
|
||||
rhoUf().boundaryFieldRef()[patchI] =
|
||||
rhoUfint.boundaryField()[patchI];
|
||||
}
|
||||
|
||||
// Calculate absolute flux from the mapped surface velocity
|
||||
phi = mesh.Sf() & rhoUf();
|
||||
|
||||
if (correctPhi)
|
||||
{
|
||||
#include "correctPhi.H"
|
||||
}
|
||||
|
||||
// Zero phi on current H-I
|
||||
const surfaceScalarField faceMask
|
||||
(
|
||||
localMin<scalar>(mesh).interpolate(cellMask)
|
||||
);
|
||||
|
||||
phi *= faceMask;
|
||||
U *= cellMask;
|
||||
|
||||
// Make the fluxes relative to the mesh-motion
|
||||
fvc::makeRelative(phi, rho, U);
|
||||
}
|
||||
|
||||
if (checkMeshCourantNo)
|
||||
{
|
||||
#include "meshCourantNo.H"
|
||||
}
|
||||
}
|
||||
|
||||
if (pimple.firstIter())
|
||||
{
|
||||
#include "rhoEqn.H"
|
||||
}
|
||||
|
||||
#include "UEqn.H"
|
||||
#include "EEqn.H"
|
||||
|
||||
// --- Pressure corrector loop
|
||||
while (pimple.correct())
|
||||
{
|
||||
#include "pEqn.H"
|
||||
}
|
||||
|
||||
if (pimple.turbCorr())
|
||||
{
|
||||
turbulence->correct();
|
||||
}
|
||||
}
|
||||
|
||||
rho = thermo.rho();
|
||||
|
||||
runTime.write();
|
||||
|
||||
runTime.printExecutionTime(Info);
|
||||
}
|
||||
|
||||
Info<< "End\n" << endl;
|
||||
|
||||
return 0;
|
||||
}
|
||||
|
||||
|
||||
// ************************************************************************* //
|
||||
@ -0,0 +1,130 @@
|
||||
|
||||
dimensionedScalar compressibility = fvc::domainIntegrate(psi);
|
||||
bool compressible = (compressibility.value() > SMALL);
|
||||
|
||||
rho = thermo.rho();
|
||||
|
||||
// Thermodynamic density needs to be updated by psi*d(p) after the
|
||||
// pressure solution
|
||||
const volScalarField psip0(psi*p);
|
||||
|
||||
volScalarField rAU(1.0/UEqn.A());
|
||||
mesh.interpolate(rAU);
|
||||
|
||||
surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rho*rAU));
|
||||
volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p_rgh));
|
||||
|
||||
surfaceScalarField phig(-rhorAUf*ghf*fvc::snGrad(rho)*mesh.magSf());
|
||||
|
||||
surfaceScalarField phiHbyA
|
||||
(
|
||||
"phiHbyA",
|
||||
fvc::flux(rho*HbyA) + phig
|
||||
);
|
||||
|
||||
if (ddtCorr)
|
||||
{
|
||||
surfaceScalarField faceMaskOld
|
||||
(
|
||||
localMin<scalar>(mesh).interpolate(cellMask.oldTime())
|
||||
);
|
||||
|
||||
phiHbyA +=
|
||||
faceMaskOld*MRF.zeroFilter(rhorAUf*fvc::ddtCorr(rho, U, phi));
|
||||
}
|
||||
|
||||
MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
|
||||
|
||||
// Update the pressure BCs to ensure flux consistency
|
||||
constrainPressure(p_rgh, rho, U, phiHbyA, rhorAUf, MRF);
|
||||
|
||||
fvScalarMatrix p_rghDDtEqn
|
||||
(
|
||||
fvc::ddt(rho) + psi*correction(fvm::ddt(p_rgh))
|
||||
+ fvc::div(phiHbyA)
|
||||
==
|
||||
fvOptions(psi, p_rgh, rho.name())
|
||||
);
|
||||
|
||||
while (pimple.correctNonOrthogonal())
|
||||
{
|
||||
fvScalarMatrix p_rghEqn
|
||||
(
|
||||
p_rghDDtEqn
|
||||
- fvm::laplacian(rhorAUf, p_rgh)
|
||||
);
|
||||
|
||||
p_rghEqn.solve(mesh.solver(p_rgh.select(pimple.finalInnerIter())));
|
||||
|
||||
if (pimple.finalNonOrthogonalIter())
|
||||
{
|
||||
// Calculate the conservative fluxes
|
||||
phi = phiHbyA + p_rghEqn.flux();
|
||||
|
||||
// Explicitly relax pressure for momentum corrector
|
||||
p_rgh.relax();
|
||||
|
||||
// Correct the momentum source with the pressure gradient flux
|
||||
// calculated from the relaxed pressure
|
||||
U =
|
||||
cellMask*
|
||||
(
|
||||
HbyA + rAU*fvc::reconstruct((phig + p_rghEqn.flux())/rhorAUf)
|
||||
);
|
||||
U.correctBoundaryConditions();
|
||||
fvOptions.correct(U);
|
||||
K = 0.5*magSqr(U);
|
||||
}
|
||||
}
|
||||
|
||||
p = p_rgh + rho*gh;
|
||||
|
||||
#include "rhoEqn.H"
|
||||
#include "compressibleContinuityErrs.H"
|
||||
|
||||
if (p_rgh.needReference())
|
||||
{
|
||||
if (!compressible)
|
||||
{
|
||||
p += dimensionedScalar
|
||||
(
|
||||
"p",
|
||||
p.dimensions(),
|
||||
pRefValue - getRefCellValue(p, pRefCell)
|
||||
);
|
||||
}
|
||||
else
|
||||
{
|
||||
p += (initialMass - fvc::domainIntegrate(psi*p))
|
||||
/compressibility;
|
||||
thermo.correctRho(psi*p - psip0);
|
||||
rho = thermo.rho();
|
||||
p_rgh = p - rho*gh;
|
||||
}
|
||||
}
|
||||
else
|
||||
{
|
||||
thermo.correctRho(psi*p - psip0);
|
||||
}
|
||||
|
||||
rho = thermo.rho();
|
||||
|
||||
{
|
||||
fvc::correctRhoUf(rhoUf, rho, U, phi);
|
||||
}
|
||||
|
||||
if (thermo.dpdt())
|
||||
{
|
||||
dpdt = fvc::ddt(p);
|
||||
|
||||
if (mesh.moving())
|
||||
{
|
||||
dpdt -= fvc::div(fvc::meshPhi(rho, U), p);
|
||||
}
|
||||
}
|
||||
|
||||
surfaceScalarField faceMask
|
||||
(
|
||||
localMin<scalar>(mesh).interpolate(cellMask)
|
||||
);
|
||||
phi *= faceMask;
|
||||
@ -0,0 +1,9 @@
|
||||
#include "readTimeControls.H"
|
||||
|
||||
correctPhi = pimple.dict().lookupOrDefault("correctPhi", false);
|
||||
|
||||
checkMeshCourantNo =
|
||||
pimple.dict().lookupOrDefault("checkMeshCourantNo", false);
|
||||
|
||||
|
||||
ddtCorr = pimple.dict().lookupOrDefault("ddtCorr", true);
|
||||
Reference in New Issue
Block a user