Merge branch 'mergeDyM'

This commit is contained in:
Henry Weller
2018-07-11 09:48:17 +01:00
308 changed files with 873 additions and 4164 deletions

View File

@ -1,9 +1,6 @@
EXE_INC = \
-I. \
-I$(FOAM_SOLVERS)/combustion/reactingFoam \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/sampling/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
@ -12,20 +9,26 @@ EXE_INC = \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/chemistryModel/lnInclude \
-I$(LIB_SRC)/ODE/lnInclude \
-I$(LIB_SRC)/combustionModels/lnInclude
-I$(LIB_SRC)/combustionModels/lnInclude \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/finiteVolume/cfdTools \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/sampling/lnInclude \
-I$(LIB_SRC)/dynamicFvMesh/lnInclude
EXE_LIBS = \
-lfiniteVolume \
-lfvOptions \
-lmeshTools \
-lsampling \
-lturbulenceModels \
-lcompressibleTurbulenceModels \
-lreactionThermophysicalModels \
-lspecie \
-lcompressibleTransportModels \
-lfluidThermophysicalModels \
-lspecie \
-lchemistryModel \
-lODE \
-lcombustionModels
-lcombustionModels \
-lreactionThermophysicalModels \
-lturbulenceModels \
-lcompressibleTurbulenceModels \
-lfiniteVolume \
-ldynamicFvMesh \
-ltopoChangerFvMesh \
-lmeshTools \
-lsampling \
-lfvOptions

View File

@ -31,12 +31,14 @@ Description
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "dynamicFvMesh.H"
#include "rhoReactionThermo.H"
#include "CombustionModel.H"
#include "turbulentFluidThermoModel.H"
#include "multivariateScheme.H"
#include "pimpleControl.H"
#include "pressureControl.H"
#include "CorrectPhi.H"
#include "fvOptions.H"
#include "localEulerDdtScheme.H"
#include "fvcSmooth.H"
@ -49,13 +51,12 @@ int main(int argc, char *argv[])
#include "setRootCaseLists.H"
#include "createTime.H"
#include "createMesh.H"
#include "createControl.H"
#include "createDynamicFvMesh.H"
#include "createDyMControls.H"
#include "initContinuityErrs.H"
#include "createFields.H"
#include "createFieldRefs.H"
#include "createRhoUfIfPresent.H"
#include "createTimeControls.H"
turbulence->validate();
@ -71,7 +72,20 @@ int main(int argc, char *argv[])
while (runTime.run())
{
#include "readTimeControls.H"
#include "readDyMControls.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 = new volScalarField
(
"divrhoU",
fvc::div(fvc::absolute(phi, rho, U))
);
}
if (LTS)
{
@ -87,11 +101,49 @@ int main(int argc, char *argv[])
Info<< "Time = " << runTime.timeName() << nl << endl;
#include "rhoEqn.H"
// --- Pressure-velocity PIMPLE corrector loop
while (pimple.loop())
{
if (pimple.firstIter() || moveMeshOuterCorrectors)
{
// Store momentum to set rhoUf for introduced faces.
autoPtr<volVectorField> rhoU;
if (rhoUf.valid())
{
rhoU = new volVectorField("rhoU", rho*U);
}
// Do any mesh changes
mesh.update();
if (mesh.changing())
{
MRF.update();
if (correctPhi)
{
// Calculate absolute flux
// from the mapped surface velocity
phi = mesh.Sf() & rhoUf();
#include "../../../compressible/rhoPimpleFoam/correctPhi.H"
// Make the fluxes relative to the mesh-motion
fvc::makeRelative(phi, rho, U);
}
if (checkMeshCourantNo)
{
#include "meshCourantNo.H"
}
}
}
if (pimple.firstIter() && !pimple.simpleRho())
{
#include "rhoEqn.H"
}
#include "UEqn.H"
#include "YEqn.H"
#include "EEqn.H"

View File

@ -3,6 +3,5 @@ cd ${0%/*} || exit 1 # Run from this directory
wclean libso BCs
wclean
wclean rhoCentralDyMFoam
#------------------------------------------------------------------------------

View File

@ -4,6 +4,6 @@ cd ${0%/*} || exit 1 # Run from this directory
# Parse arguments for library compilation
. $WM_PROJECT_DIR/wmake/scripts/AllwmakeParseArguments
(wmake $targetType BCs && wmake $targetType && wmake $targetType rhoCentralDyMFoam)
(wmake $targetType BCs && wmake $targetType)
#------------------------------------------------------------------------------

View File

@ -6,7 +6,7 @@ EXE_INC = \
-I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \
-I$(LIB_SRC)/dynamicMesh/lnInclude \
-I$(LIB_SRC)/dynamicFvMesh/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude
EXE_LIBS = \
@ -18,4 +18,6 @@ EXE_LIBS = \
-lrhoCentralFoam \
-lturbulenceModels \
-lcompressibleTurbulenceModels \
-ldynamicFvMesh \
-ltopoChangerFvMesh \
-lmeshTools

View File

@ -1,3 +0,0 @@
rhoCentralDyMFoam.C
EXE = $(FOAM_APPBIN)/rhoCentralDyMFoam

View File

@ -1,26 +0,0 @@
EXE_INC = \
-I.. \
-I../BCs/lnInclude \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/transportModels/compressible/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \
-I$(LIB_SRC)/dynamicMesh/lnInclude \
-I$(LIB_SRC)/dynamicFvMesh/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude
EXE_LIBS = \
-lfiniteVolume \
-lfvOptions \
-lcompressibleTransportModels \
-lfluidThermophysicalModels \
-lspecie \
-lrhoCentralFoam \
-lturbulenceModels \
-lcompressibleTurbulenceModels \
-ldynamicMesh \
-ldynamicFvMesh \
-ltopoChangerFvMesh \
-lmeshTools

View File

@ -1,269 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation
\\/ 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
rhoCentralDyMFoam
Description
Density-based compressible flow solver based on central-upwind schemes of
Kurganov and Tadmor with support for mesh-motion and topology changes.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "dynamicFvMesh.H"
#include "psiThermo.H"
#include "turbulentFluidThermoModel.H"
#include "fixedRhoFvPatchScalarField.H"
#include "directionInterpolate.H"
#include "localEulerDdtScheme.H"
#include "fvcSmooth.H"
#include "motionSolver.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
#define NO_CONTROL
#include "postProcess.H"
#include "setRootCaseLists.H"
#include "createTime.H"
#include "createDynamicFvMesh.H"
#include "createFields.H"
#include "createFieldRefs.H"
#include "createTimeControls.H"
turbulence->validate();
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#include "readFluxScheme.H"
dimensionedScalar v_zero("v_zero", dimVolume/dimTime, 0.0);
// Courant numbers used to adjust the time-step
scalar CoNum = 0.0;
scalar meanCoNum = 0.0;
Info<< "\nStarting time loop\n" << endl;
while (runTime.run())
{
#include "readTimeControls.H"
#include "setDeltaT.H"
runTime++;
Info<< "Time = " << runTime.timeName() << nl << endl;
// Do any mesh changes
mesh.update();
// --- Directed interpolation of primitive fields onto faces
surfaceScalarField rho_pos(interpolate(rho, pos));
surfaceScalarField rho_neg(interpolate(rho, neg));
surfaceVectorField rhoU_pos(interpolate(rhoU, pos, U.name()));
surfaceVectorField rhoU_neg(interpolate(rhoU, neg, U.name()));
volScalarField rPsi("rPsi", 1.0/psi);
surfaceScalarField rPsi_pos(interpolate(rPsi, pos, T.name()));
surfaceScalarField rPsi_neg(interpolate(rPsi, neg, T.name()));
surfaceScalarField e_pos(interpolate(e, pos, T.name()));
surfaceScalarField e_neg(interpolate(e, neg, T.name()));
surfaceVectorField U_pos("U_pos", rhoU_pos/rho_pos);
surfaceVectorField U_neg("U_neg", rhoU_neg/rho_neg);
surfaceScalarField p_pos("p_pos", rho_pos*rPsi_pos);
surfaceScalarField p_neg("p_neg", rho_neg*rPsi_neg);
surfaceScalarField phiv_pos("phiv_pos", U_pos & mesh.Sf());
surfaceScalarField phiv_neg("phiv_neg", U_neg & mesh.Sf());
// Make fluxes relative to mesh-motion
if (mesh.moving())
{
phiv_pos -= mesh.phi();
phiv_neg -= mesh.phi();
}
volScalarField c("c", sqrt(thermo.Cp()/thermo.Cv()*rPsi));
surfaceScalarField cSf_pos
(
"cSf_pos",
interpolate(c, pos, T.name())*mesh.magSf()
);
surfaceScalarField cSf_neg
(
"cSf_neg",
interpolate(c, neg, T.name())*mesh.magSf()
);
surfaceScalarField ap
(
"ap",
max(max(phiv_pos + cSf_pos, phiv_neg + cSf_neg), v_zero)
);
surfaceScalarField am
(
"am",
min(min(phiv_pos - cSf_pos, phiv_neg - cSf_neg), v_zero)
);
surfaceScalarField a_pos("a_pos", ap/(ap - am));
surfaceScalarField amaxSf("amaxSf", max(mag(am), mag(ap)));
surfaceScalarField aSf("aSf", am*a_pos);
if (fluxScheme == "Tadmor")
{
aSf = -0.5*amaxSf;
a_pos = 0.5;
}
surfaceScalarField a_neg("a_neg", 1.0 - a_pos);
phiv_pos *= a_pos;
phiv_neg *= a_neg;
surfaceScalarField aphiv_pos("aphiv_pos", phiv_pos - aSf);
surfaceScalarField aphiv_neg("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 "centralCourantNo.H"
phi = aphiv_pos*rho_pos + aphiv_neg*rho_neg;
surfaceVectorField phiUp
(
(aphiv_pos*rhoU_pos + aphiv_neg*rhoU_neg)
+ (a_pos*p_pos + a_neg*p_neg)*mesh.Sf()
);
surfaceScalarField phiEp
(
"phiEp",
aphiv_pos*(rho_pos*(e_pos + 0.5*magSqr(U_pos)) + p_pos)
+ aphiv_neg*(rho_neg*(e_neg + 0.5*magSqr(U_neg)) + p_neg)
+ aSf*p_pos - aSf*p_neg
);
// Make flux for pressure-work absolute
if (mesh.moving())
{
phiEp += mesh.phi()*(a_pos*p_pos + a_neg*p_neg);
}
volScalarField muEff("muEff", turbulence->muEff());
volTensorField tauMC("tauMC", muEff*dev2(Foam::T(fvc::grad(U))));
// --- Solve density
solve(fvm::ddt(rho) + fvc::div(phi));
// --- Solve momentum
solve(fvm::ddt(rhoU) + fvc::div(phiUp));
U.ref() =
rhoU()
/rho();
U.correctBoundaryConditions();
rhoU.boundaryFieldRef() == rho.boundaryField()*U.boundaryField();
if (!inviscid)
{
solve
(
fvm::ddt(rho, U) - fvc::ddt(rho, U)
- fvm::laplacian(muEff, U)
- fvc::div(tauMC)
);
rhoU = rho*U;
}
// --- Solve energy
surfaceScalarField sigmaDotU
(
"sigmaDotU",
(
fvc::interpolate(muEff)*mesh.magSf()*fvc::snGrad(U)
+ fvc::dotInterpolate(mesh.Sf(), tauMC)
)
& (a_pos*U_pos + a_neg*U_neg)
);
solve
(
fvm::ddt(rhoE)
+ fvc::div(phiEp)
- fvc::div(sigmaDotU)
);
e = rhoE/rho - 0.5*magSqr(U);
e.correctBoundaryConditions();
thermo.correct();
rhoE.boundaryFieldRef() ==
rho.boundaryField()*
(
e.boundaryField() + 0.5*magSqr(U.boundaryField())
);
if (!inviscid)
{
solve
(
fvm::ddt(rho, e) - fvc::ddt(rho, e)
- fvm::laplacian(turbulence->alphaEff(), e)
);
thermo.correct();
rhoE = rho*(e + 0.5*magSqr(U));
}
p.ref() =
rho()
/psi();
p.correctBoundaryConditions();
rho.boundaryFieldRef() == psi.boundaryField()*p.boundaryField();
turbulence->correct();
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
}
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

View File

@ -26,11 +26,12 @@ Application
Description
Density-based compressible flow solver based on central-upwind schemes of
Kurganov and Tadmor.
Kurganov and Tadmor with support for mesh-motion and topology changes.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "dynamicFvMesh.H"
#include "psiThermo.H"
#include "turbulentFluidThermoModel.H"
#include "fixedRhoFvPatchScalarField.H"
@ -47,7 +48,7 @@ int main(int argc, char *argv[])
#include "setRootCaseLists.H"
#include "createTime.H"
#include "createMesh.H"
#include "createDynamicFvMesh.H"
#include "createFields.H"
#include "createFieldRefs.H"
#include "createTimeControls.H"
@ -68,6 +69,17 @@ int main(int argc, char *argv[])
while (runTime.run())
{
#include "readTimeControls.H"
if (!LTS)
{
#include "setDeltaT.H"
runTime++;
// Do any mesh changes
mesh.update();
}
// --- Directed interpolation of primitive fields onto faces
surfaceScalarField rho_pos(interpolate(rho, pos));
@ -92,6 +104,13 @@ int main(int argc, char *argv[])
surfaceScalarField phiv_pos("phiv_pos", U_pos & mesh.Sf());
surfaceScalarField phiv_neg("phiv_neg", U_neg & mesh.Sf());
// Make fluxes relative to mesh-motion
if (mesh.moving())
{
phiv_pos -= mesh.phi();
phiv_neg -= mesh.phi();
}
volScalarField c("c", sqrt(thermo.Cp()/thermo.Cv()*rPsi));
surfaceScalarField cSf_pos
(
@ -140,18 +159,12 @@ int main(int argc, char *argv[])
amaxSf = max(mag(aphiv_pos), mag(aphiv_neg));
#include "centralCourantNo.H"
#include "readTimeControls.H"
if (LTS)
{
#include "setRDeltaT.H"
runTime++;
}
else
{
#include "setDeltaT.H"
}
runTime++;
Info<< "Time = " << runTime.timeName() << nl << endl;
@ -171,6 +184,12 @@ int main(int argc, char *argv[])
+ aSf*p_pos - aSf*p_neg
);
// Make flux for pressure-work absolute
if (mesh.moving())
{
phiEp += mesh.phi()*(a_pos*p_pos + a_neg*p_neg);
}
volScalarField muEff("muEff", turbulence->muEff());
volTensorField tauMC("tauMC", muEff*dev2(Foam::T(fvc::grad(U))));

View File

@ -7,8 +7,7 @@ EXE_INC = \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/sampling/lnInclude \
-I$(LIB_SRC)/dynamicFvMesh/lnInclude \
-I$(LIB_SRC)/dynamicMesh/lnInclude
-I$(LIB_SRC)/dynamicFvMesh/lnInclude
EXE_LIBS = \
-lcompressibleTransportModels \
@ -17,9 +16,8 @@ EXE_LIBS = \
-lturbulenceModels \
-lcompressibleTurbulenceModels \
-lfiniteVolume \
-lmeshTools \
-lsampling \
-lfvOptions \
-ldynamicFvMesh \
-ltopoChangerFvMesh \
-ldynamicMesh
-lmeshTools \
-lsampling \
-lfvOptions

View File

@ -37,7 +37,6 @@ Description
#include "dynamicFvMesh.H"
#include "fluidThermo.H"
#include "turbulentFluidThermoModel.H"
#include "bound.H"
#include "pimpleControl.H"
#include "pressureControl.H"
#include "CorrectPhi.H"

View File

@ -1,21 +0,0 @@
{
fvScalarMatrix EEqn
(
fvm::ddt(rho, e) + fvm::div(phi, e)
+ fvc::ddt(rho, K) + fvc::div(phi, K)
+ fvc::div(fvc::absolute(phi/fvc::interpolate(rho), U), p, "div(phiv,p)")
- fvm::laplacian(turbulence->alphaEff(), e)
==
fvOptions(rho, e)
);
EEqn.relax();
fvOptions.constrain(EEqn);
EEqn.solve();
fvOptions.correct(e);
thermo.correct();
}

View File

@ -1,3 +0,0 @@
sonicFoam.C
EXE = $(FOAM_APPBIN)/sonicFoam

View File

@ -1,19 +0,0 @@
EXE_INC = \
-I$(LIB_SRC)/transportModels/compressible/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/sampling/lnInclude \
-I$(LIB_SRC)/finiteVolume/lnInclude
EXE_LIBS = \
-lcompressibleTransportModels \
-lfluidThermophysicalModels \
-lspecie \
-lturbulenceModels \
-lcompressibleTurbulenceModels \
-lmeshTools \
-lsampling \
-lfvOptions \
-lfiniteVolume

View File

@ -1,24 +0,0 @@
// Solve the Momentum equation
MRF.correctBoundaryVelocity(U);
fvVectorMatrix UEqn
(
fvm::ddt(rho, U) + fvm::div(phi, U)
+ MRF.DDt(rho, U)
+ turbulence->divDevRhoReff(U)
==
fvOptions(rho, U)
);
UEqn.relax();
fvOptions.constrain(UEqn);
if (pimple.momentumPredictor())
{
solve(UEqn == -fvc::grad(p));
fvOptions.correct(U);
K = 0.5*magSqr(U);
}

View File

@ -1,2 +0,0 @@
volScalarField& e = thermo.he();
const volScalarField& psi = thermo.psi();

View File

@ -1,57 +0,0 @@
Info<< "Reading thermophysical properties\n" << endl;
autoPtr<psiThermo> pThermo
(
psiThermo::New(mesh)
);
psiThermo& thermo = pThermo();
thermo.validate(args.executable(), "e");
volScalarField& p = thermo.p();
volScalarField rho
(
IOobject
(
"rho",
runTime.timeName(),
mesh
),
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"
mesh.setFluxRequired(p.name());
Info<< "Creating turbulence model\n" << endl;
autoPtr<compressible::turbulenceModel> turbulence
(
compressible::turbulenceModel::New
(
rho,
U,
phi,
thermo
)
);
Info<< "Creating field kinetic energy K\n" << endl;
volScalarField K("K", 0.5*magSqr(U));
#include "createMRF.H"
#include "createFvOptions.H"

View File

@ -1,44 +0,0 @@
rho = thermo.rho();
volScalarField rAU(1.0/UEqn.A());
surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rho*rAU));
volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p));
surfaceScalarField phid
(
"phid",
fvc::interpolate(psi)
*(
fvc::flux(HbyA)
+ MRF.zeroFilter(rhorAUf*fvc::ddtCorr(rho, U, phi)/fvc::interpolate(rho))
)
);
MRF.makeRelative(fvc::interpolate(psi), phid);
// Non-orthogonal pressure corrector loop
while (pimple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
(
fvm::ddt(psi, p)
+ fvm::div(phid, p)
- fvm::laplacian(rhorAUf, p)
==
fvOptions(psi, p, rho.name())
);
pEqn.solve();
if (pimple.finalNonOrthogonalIter())
{
phi = pEqn.flux();
}
}
#include "rhoEqn.H"
#include "compressibleContinuityErrs.H"
U = HbyA - rAU*fvc::grad(p);
U.correctBoundaryConditions();
fvOptions.correct(U);
K = 0.5*magSqr(U);

View File

@ -1,3 +0,0 @@
sonicDyMFoam.C
EXE = $(FOAM_APPBIN)/sonicDyMFoam

View File

@ -1,28 +0,0 @@
EXE_INC = \
-I.. \
-I../../rhoPimpleFoam \
-I$(LIB_SRC)/transportModels/compressible/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \
-I$(LIB_SRC)/finiteVolume/cfdTools \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/sampling/lnInclude \
-I$(LIB_SRC)/dynamicFvMesh/lnInclude \
-I$(LIB_SRC)/dynamicMesh/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
EXE_LIBS = \
-lcompressibleTransportModels \
-lfluidThermophysicalModels \
-lspecie \
-lturbulenceModels \
-lcompressibleTurbulenceModels \
-lfiniteVolume \
-lmeshTools \
-lsampling \
-lfvOptions \
-ldynamicFvMesh \
-ltopoChangerFvMesh \
-ldynamicMesh \
-lmeshTools

View File

@ -1,51 +0,0 @@
rho = thermo.rho();
volScalarField rAU(1.0/UEqn.A());
surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rho*rAU));
volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p));
surfaceScalarField phid
(
"phid",
fvc::interpolate(psi)
*(
fvc::flux(HbyA)
+ rhorAUf*fvc::ddtCorr(rho, U, rhoUf)/fvc::interpolate(rho)
)
);
fvc::makeRelative(phid, psi, U);
MRF.makeRelative(fvc::interpolate(psi), phid);
// Non-orthogonal pressure corrector loop
while (pimple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
(
fvm::ddt(psi, p)
+ fvm::div(phid, p)
- fvm::laplacian(rhorAUf, p)
==
fvOptions(psi, p, rho.name())
);
pEqn.solve();
if (pimple.finalNonOrthogonalIter())
{
phi = pEqn.flux();
}
}
#include "rhoEqn.H"
#include "compressibleContinuityErrs.H"
U = HbyA - rAU*fvc::grad(p);
U.correctBoundaryConditions();
fvOptions.correct(U);
K = 0.5*magSqr(U);
{
rhoUf = fvc::interpolate(rho*U);
surfaceVectorField n(mesh.Sf()/mesh.magSf());
rhoUf += n*(fvc::absolute(phi, rho, U)/mesh.magSf() - (n & rhoUf));
}

View File

@ -1,151 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation
\\/ 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
sonicDyMFoam
Group
grpCompressibleSolvers grpMovingMeshSolvers
Description
Transient solver for trans-sonic/supersonic, turbulent flow of a
compressible gas, with optional mesh motion and mesh topology changes.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "dynamicFvMesh.H"
#include "psiThermo.H"
#include "turbulentFluidThermoModel.H"
#include "pimpleControl.H"
#include "CorrectPhi.H"
#include "fvOptions.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
#include "postProcess.H"
#include "setRootCaseLists.H"
#include "createTime.H"
#include "createDynamicFvMesh.H"
#include "createDyMControls.H"
#include "createFields.H"
#include "createFieldRefs.H"
#include "createRhoUf.H"
#include "compressibleCourantNo.H"
#include "setInitialDeltaT.H"
#include "initContinuityErrs.H"
turbulence->validate();
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
while (runTime.run())
{
#include "readDyMControls.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
volScalarField divrhoU
(
"divrhoU",
fvc::div(fvc::absolute(phi, rho, U))
);
#include "compressibleCourantNo.H"
#include "setDeltaT.H"
runTime++;
Info<< "Time = " << runTime.timeName() << nl << endl;
// Store momentum to set rhoUf for introduced faces.
volVectorField rhoU("rhoU", rho*U);
// Do any mesh changes
mesh.update();
if (mesh.changing())
{
MRF.update();
if (correctPhi)
{
// Calculate absolute flux from the mapped surface velocity
phi = mesh.Sf() & rhoUf;
#include "correctPhi.H"
// Make the fluxes relative to the mesh-motion
fvc::makeRelative(phi, rho, U);
}
}
if (checkMeshCourantNo)
{
#include "meshCourantNo.H"
}
}
#include "rhoEqn.H"
Info<< "rhoEqn max/min : " << max(rho).value()
<< " " << min(rho).value() << endl;
// --- Pressure-velocity PIMPLE corrector loop
while (pimple.loop())
{
#include "UEqn.H"
#include "EEqn.H"
// --- Pressure corrector loop
while (pimple.correct())
{
#include "pEqn.H"
}
if (pimple.turbCorr())
{
turbulence->correct();
}
}
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
}
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

View File

@ -1,103 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation
\\/ 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
sonicFoam
Group
grpCompressibleSolvers
Description
Transient solver for trans-sonic/supersonic, turbulent flow of a
compressible gas.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "psiThermo.H"
#include "turbulentFluidThermoModel.H"
#include "pimpleControl.H"
#include "fvOptions.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
#include "postProcess.H"
#include "setRootCaseLists.H"
#include "createTime.H"
#include "createMesh.H"
#include "createControl.H"
#include "createFields.H"
#include "createFieldRefs.H"
#include "initContinuityErrs.H"
turbulence->validate();
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
while (runTime.loop())
{
Info<< "Time = " << runTime.timeName() << nl << endl;
#include "compressibleCourantNo.H"
#include "rhoEqn.H"
// --- Pressure-velocity PIMPLE corrector loop
while (pimple.loop())
{
#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();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
}
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

View File

@ -1,3 +0,0 @@
sonicLiquidFoam.C
EXE = $(FOAM_APPBIN)/sonicLiquidFoam

View File

@ -1,7 +0,0 @@
EXE_INC = \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude
EXE_LIBS = \
-lfiniteVolume \
-lmeshTools

View File

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

View File

@ -1,49 +0,0 @@
#include "readThermodynamicProperties.H"
#include "readTransportProperties.H"
Info<< "Reading field p\n" << endl;
volScalarField p
(
IOobject
(
"p",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
Info<< "Reading field U\n" << endl;
volVectorField U
(
IOobject
(
"U",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField rho
(
IOobject
(
"rho",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
rhoO + psi*p
);
#include "compressibleCreatePhi.H"
mesh.setFluxRequired(p.name());

View File

@ -1,37 +0,0 @@
Info<< "Reading thermodynamicProperties\n" << endl;
IOdictionary thermodynamicProperties
(
IOobject
(
"thermodynamicProperties",
runTime.constant(),
mesh,
IOobject::MUST_READ_IF_MODIFIED,
IOobject::NO_WRITE
)
);
dimensionedScalar rho0
(
"rho0",
dimDensity,
thermodynamicProperties
);
dimensionedScalar p0
(
"p0",
dimPressure,
thermodynamicProperties
);
dimensionedScalar psi
(
"psi",
dimCompressibility,
thermodynamicProperties
);
// Density offset, i.e. the constant part of the density
dimensionedScalar rhoO("rhoO", rho0 - psi*p0);

View File

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

View File

@ -1,132 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation
\\/ 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
sonicLiquidFoam
Description
Transient solver for trans-sonic/supersonic, laminar flow of a
compressible liquid.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "pimpleControl.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
#include "postProcess.H"
#include "setRootCaseLists.H"
#include "createTime.H"
#include "createMesh.H"
#include "createControl.H"
#include "createFields.H"
#include "initContinuityErrs.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
while (runTime.loop())
{
Info<< "Time = " << runTime.timeName() << nl << endl;
#include "compressibleCourantNo.H"
solve(fvm::ddt(rho) + fvc::div(phi));
// --- Pressure-velocity PIMPLE corrector loop
while (pimple.loop())
{
fvVectorMatrix UEqn
(
fvm::ddt(rho, U)
+ fvm::div(phi, U)
- fvm::laplacian(mu, U)
);
solve(UEqn == -fvc::grad(p));
// --- Pressure corrector loop
while (pimple.correct())
{
volScalarField rAU("rAU", 1.0/UEqn.A());
surfaceScalarField rhorAUf
(
"rhorAUf",
fvc::interpolate(rho*rAU)
);
U = rAU*UEqn.H();
surfaceScalarField phid
(
"phid",
psi
*(
fvc::flux(U)
+ rhorAUf*fvc::ddtCorr(rho, U, phi)/fvc::interpolate(rho)
)
);
phi = (rhoO/psi)*phid;
fvScalarMatrix pEqn
(
fvm::ddt(psi, p)
+ fvc::div(phi)
+ fvm::div(phid, p)
- fvm::laplacian(rhorAUf, p)
);
pEqn.solve();
phi += pEqn.flux();
solve(fvm::ddt(rho) + fvc::div(phi));
#include "compressibleContinuityErrs.H"
U -= rAU*fvc::grad(p);
U.correctBoundaryConditions();
}
}
rho = rhoO + psi*p;
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
}
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

View File

@ -15,7 +15,7 @@ volScalarField h
);
Info<< "Reading field h0 if present\n" << endl;
volScalarField h0
const volScalarField h0
(
IOobject
(

View File

@ -5,7 +5,5 @@ wclean libso DPMTurbulenceModels
wclean
wclean MPPICFoam
wclean DPMDyMFoam
wclean DPMDyMFoam/MPPICDyMFoam
#------------------------------------------------------------------------------

View File

@ -8,7 +8,5 @@ wmake $targetType DPMTurbulenceModels
wmake $targetType
wmake $targetType MPPICFoam
wmake $targetType DPMDyMFoam
wmake $targetType DPMDyMFoam/MPPICDyMFoam
#------------------------------------------------------------------------------

View File

@ -1,166 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2017-2018 OpenFOAM Foundation
\\/ 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
DPMDyMFoam
Description
Transient solver for the coupled transport of a single kinematic particle
cloud including the effect of the volume fraction of particles on the
continuous phase, with optional mesh motion and mesh topology changes.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "dynamicFvMesh.H"
#include "singlePhaseTransportModel.H"
#include "PhaseIncompressibleTurbulenceModel.H"
#include "pimpleControl.H"
#include "CorrectPhi.H"
#ifdef MPPIC
#include "basicKinematicMPPICCloud.H"
#define basicKinematicTypeCloud basicKinematicMPPICCloud
#else
#include "basicKinematicCollidingCloud.H"
#define basicKinematicTypeCloud basicKinematicCollidingCloud
#endif
int main(int argc, char *argv[])
{
argList::addOption
(
"cloudName",
"name",
"specify alternative cloud name. default is 'kinematicCloud'"
);
#include "postProcess.H"
#include "setRootCaseLists.H"
#include "createTime.H"
#include "createDynamicFvMesh.H"
#include "createDyMControls.H"
#include "createFields.H"
#include "createUcf.H"
#include "initContinuityErrs.H"
Info<< "\nStarting time loop\n" << endl;
while (runTime.run())
{
#include "readDyMControls.H"
#include "CourantNo.H"
#include "setDeltaT.H"
runTime++;
Info<< "Time = " << runTime.timeName() << nl << endl;
// Store the particle positions
kinematicCloud.storeGlobalPositions();
mesh.update();
// Calculate absolute flux from the mapped surface velocity
phic = mesh.Sf() & Ucf;
if (mesh.changing() && correctPhi)
{
#include "correctPhic.H"
}
// Make the flux relative to the mesh motion
fvc::makeRelative(phic, Uc);
if (mesh.changing() && checkMeshCourantNo)
{
#include "meshCourantNo.H"
}
continuousPhaseTransport.correct();
muc = rhoc*continuousPhaseTransport.nu();
Info<< "Evolving " << kinematicCloud.name() << endl;
kinematicCloud.evolve();
// Update continuous phase volume fraction field
alphac = max(1.0 - kinematicCloud.theta(), alphacMin);
alphac.correctBoundaryConditions();
alphacf = fvc::interpolate(alphac);
alphaPhic = alphacf*phic;
fvVectorMatrix cloudSU(kinematicCloud.SU(Uc));
volVectorField cloudVolSUSu
(
IOobject
(
"cloudVolSUSu",
runTime.timeName(),
mesh
),
mesh,
dimensionedVector
(
"0",
cloudSU.dimensions()/dimVolume,
Zero
),
zeroGradientFvPatchVectorField::typeName
);
cloudVolSUSu.primitiveFieldRef() = -cloudSU.source()/mesh.V();
cloudVolSUSu.correctBoundaryConditions();
cloudSU.source() = Zero;
// --- Pressure-velocity PIMPLE corrector loop
while (pimple.loop())
{
#include "UcEqn.H"
// --- PISO loop
while (pimple.correct())
{
#include "pEqn.H"
}
if (pimple.turbCorr())
{
continuousPhaseTurbulence->correct();
}
}
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
}
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

View File

@ -1,3 +0,0 @@
MPPICDyMFoam.C
EXE = $(FOAM_APPBIN)/MPPICDyMFoam

View File

@ -1,42 +0,0 @@
EXE_INC = \
-I.. \
-I../.. \
-I../../DPMTurbulenceModels/lnInclude \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/lagrangian/basic/lnInclude \
-I$(LIB_SRC)/lagrangian/intermediate/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
-I$(LIB_SRC)/transportModels/compressible/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/radiation/lnInclude \
-I$(LIB_SRC)/transportModels \
-I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \
-I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/phaseIncompressible/lnInclude \
-I$(LIB_SRC)/regionModels/regionModel/lnInclude \
-I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \
-I$(LIB_SRC)/dynamicFvMesh/lnInclude \
-I$(LIB_SRC)/dynamicMesh/lnInclude
EXE_LIBS = \
-lfiniteVolume \
-lfvOptions \
-lmeshTools \
-llagrangian \
-llagrangianIntermediate \
-llagrangianTurbulence \
-lspecie \
-lradiationModels \
-lincompressibleTransportModels \
-lturbulenceModels \
-lincompressibleTurbulenceModels \
-lDPMTurbulenceModels \
-lregionModels \
-lsurfaceFilmModels \
-lsampling \
-ldynamicFvMesh \
-ltopoChangerFvMesh \
-ldynamicMesh

View File

@ -1,3 +0,0 @@
DPMDyMFoam.C
EXE = $(FOAM_APPBIN)/DPMDyMFoam

View File

@ -1,41 +0,0 @@
EXE_INC = \
-I.. \
-I../DPMTurbulenceModels/lnInclude \
-I$(LIB_SRC)/lagrangian/basic/lnInclude \
-I$(LIB_SRC)/lagrangian/intermediate/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
-I$(LIB_SRC)/transportModels/compressible/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/radiation/lnInclude \
-I$(LIB_SRC)/transportModels \
-I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \
-I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/phaseIncompressible/lnInclude \
-I$(LIB_SRC)/regionModels/regionModel/lnInclude \
-I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/dynamicFvMesh/lnInclude \
-I$(LIB_SRC)/dynamicMesh/lnInclude
EXE_LIBS = \
-llagrangian \
-llagrangianIntermediate \
-llagrangianTurbulence \
-lspecie \
-lradiationModels \
-lincompressibleTransportModels \
-lturbulenceModels \
-lincompressibleTurbulenceModels \
-lDPMTurbulenceModels \
-lregionModels \
-lsurfaceFilmModels \
-lsampling \
-lfiniteVolume \
-lfvOptions \
-lmeshTools \
-ldynamicFvMesh \
-ltopoChangerFvMesh \
-ldynamicMesh

View File

@ -1,62 +0,0 @@
{
volVectorField HbyA(constrainHbyA(rAUc*UcEqn.H(), Uc, p));
surfaceScalarField phiHbyA
(
"phiHbyA",
(
fvc::flux(HbyA)
+ alphacf*rAUcf*fvc::ddtCorr(Uc, Ucf)
)
);
if (p.needReference())
{
fvc::makeRelative(phiHbyA, Uc);
adjustPhi(phiHbyA, Uc, p);
fvc::makeAbsolute(phiHbyA, Uc);
}
phiHbyA += phicForces;
// Update the pressure BCs to ensure flux consistency
constrainPressure(p, Uc, phiHbyA, rAUcf);
// Non-orthogonal pressure corrector loop
while (pimple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
(
fvm::laplacian(alphacf*rAUcf, p)
==
fvc::ddt(alphac) + fvc::div(alphacf*phiHbyA)
);
pEqn.setReference(pRefCell, pRefValue);
pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter())));
if (pimple.finalNonOrthogonalIter())
{
phic = phiHbyA - pEqn.flux()/alphacf;
p.relax();
Uc = HbyA
+ rAUc
*fvc::reconstruct((phicForces - pEqn.flux()/alphacf)/rAUcf);
Uc.correctBoundaryConditions();
{
Ucf = fvc::interpolate(Uc);
surfaceVectorField n(mesh.Sf()/mesh.magSf());
Ucf += n*(phic/mesh.magSf() - (n & Ucf));
}
// Make the fluxes relative to the mesh motion
fvc::makeRelative(phic, Uc);
}
}
}
#include "continuityErrs.H"

View File

@ -27,14 +27,16 @@ Application
Description
Transient solver for the coupled transport of a single kinematic particle
cloud including the effect of the volume fraction of particles on the
continuous phase.
continuous phase, with optional mesh motion and mesh topology changes.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "dynamicFvMesh.H"
#include "singlePhaseTransportModel.H"
#include "PhaseIncompressibleTurbulenceModel.H"
#include "pimpleControl.H"
#include "CorrectPhi.H"
#ifdef MPPIC
#include "basicKinematicMPPICCloud.H"
@ -57,17 +59,17 @@ int main(int argc, char *argv[])
#include "setRootCaseLists.H"
#include "createTime.H"
#include "createMesh.H"
#include "createControl.H"
#include "createTimeControls.H"
#include "createDynamicFvMesh.H"
#include "createDyMControls.H"
#include "createFields.H"
#include "createUcfIfPresent.H"
#include "initContinuityErrs.H"
Info<< "\nStarting time loop\n" << endl;
while (runTime.run())
{
#include "readTimeControls.H"
#include "readDyMControls.H"
#include "CourantNo.H"
#include "setDeltaT.H"
@ -75,6 +77,30 @@ int main(int argc, char *argv[])
Info<< "Time = " << runTime.timeName() << nl << endl;
// Store the particle positions
kinematicCloud.storeGlobalPositions();
mesh.update();
if (mesh.changing())
{
if (correctPhi)
{
// Calculate absolute flux from the mapped surface velocity
phic = mesh.Sf() & Ucf();
#include "correctPhic.H"
// Make the flux relative to the mesh motion
fvc::makeRelative(phic, Uc);
}
if (checkMeshCourantNo)
{
#include "meshCourantNo.H"
}
}
continuousPhaseTransport.correct();
muc = rhoc*continuousPhaseTransport.nu();

View File

@ -1,8 +1,6 @@
EXE_INC = \
-I.. \
-I../DPMTurbulenceModels/lnInclude \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/lagrangian/basic/lnInclude \
-I$(LIB_SRC)/lagrangian/intermediate/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
@ -17,20 +15,20 @@ EXE_INC = \
-I$(LIB_SRC)/TurbulenceModels/phaseIncompressible/lnInclude \
-I$(LIB_SRC)/regionModels/regionModel/lnInclude \
-I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/dynamicFvMesh/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude
EXE_LIBS = \
-lfiniteVolume \
-lfvOptions \
-lmeshTools \
-llagrangian \
-llagrangianIntermediate \
-llagrangianTurbulence \
-lspecie \
-lradiationModels \
-lincompressibleTransportModels \
-lturbulenceModels \
-lincompressibleTurbulenceModels \
-lDPMTurbulenceModels \
-lregionModels \
-lsurfaceFilmModels \
-lsampling
-lfiniteVolume \
-ldynamicFvMesh \
-ltopoChangerFvMesh \
-lfvOptions \
-lmeshTools

View File

@ -15,21 +15,19 @@ EXE_INC = \
-I$(LIB_SRC)/regionModels/regionModel/lnInclude \
-I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/dynamicFvMesh/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude
EXE_LIBS = \
-llagrangian \
-llagrangianIntermediate \
-llagrangianTurbulence \
-lspecie \
-lradiationModels \
-lincompressibleTransportModels \
-lturbulenceModels \
-lincompressibleTurbulenceModels \
-lDPMTurbulenceModels \
-lregionModels \
-lsurfaceFilmModels \
-lsampling \
-lfiniteVolume \
-ldynamicFvMesh \
-ltopoChangerFvMesh \
-lfvOptions \
-lmeshTools

View File

@ -21,21 +21,35 @@ License
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Application
MPPICDyMFoam
Global
createUcfIfPresent
Description
Transient solver for the coupled transport of a single kinematic particle
cloud including the effect of the volume fraction of particles on the
continuous phase. Multi-Phase Particle In Cell (MPPIC) modeling is used to
represent collisions without resolving particle-particle interactions,
with optional mesh motion and mesh topology changes.
Creates and initialises the continuous phase face velocity field Ufc
if required.
\*---------------------------------------------------------------------------*/
#define MPPIC
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#include "DPMDyMFoam.C"
autoPtr<surfaceVectorField> Ucf;
if (mesh.dynamic())
{
Info<< "Constructing continuous phase face velocity Ucf\n" << endl;
Ucf = new surfaceVectorField
(
IOobject
(
"Ucf",
runTime.timeName(),
mesh,
IOobject::READ_IF_PRESENT,
IOobject::AUTO_WRITE
),
fvc::interpolate(Uc)
);
}
// ************************************************************************* //

View File

@ -6,13 +6,15 @@
"phiHbyA",
(
fvc::flux(HbyA)
+ alphacf*rAUcf*fvc::ddtCorr(Uc, phic)
+ alphacf*rAUcf*fvc::ddtCorr(Uc, phic, Ucf)
)
);
if (p.needReference())
{
fvc::makeRelative(phiHbyA, Uc);
adjustPhi(phiHbyA, Uc, p);
fvc::makeAbsolute(phiHbyA, Uc);
}
phiHbyA += phicForces;
@ -41,8 +43,15 @@
p.relax();
Uc = HbyA
+ rAUc*fvc::reconstruct((phicForces - pEqn.flux()/alphacf)/rAUcf);
+ rAUc
*fvc::reconstruct((phicForces - pEqn.flux()/alphacf)/rAUcf);
Uc.correctBoundaryConditions();
// Correct Ucf if the mesh is moving
fvc::correctUf(Ucf, Uc, phic);
// Make the fluxes relative to the mesh motion
fvc::makeRelative(phic, Uc);
}
}
}

View File

@ -1,33 +1,26 @@
EXE_INC = \
-I$(LIB_SRC)/lagrangian/basic/lnInclude \
-I$(LIB_SRC)/lagrangian/intermediate/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
-I$(LIB_SRC)/transportModels/compressible/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/radiation/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \
-I$(LIB_SRC)/transportModels \
-I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/regionModels/regionModel/lnInclude \
-I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude
-I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \
-I$(LIB_SRC)/dynamicFvMesh/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude
EXE_LIBS = \
-llagrangian \
-llagrangianIntermediate \
-llagrangianTurbulence \
-lcompressibleTransportModels \
-lfluidThermophysicalModels \
-lspecie \
-lradiationModels \
-lturbulenceModels \
-lincompressibleTurbulenceModels \
-lincompressibleTransportModels \
-lfiniteVolume \
-lfvOptions \
-lmeshTools \
-lregionModels \
-lsurfaceFilmModels
-lsurfaceFilmModels \
-ldynamicFvMesh \
-ltopoChangerFvMesh \
-lmeshTools

View File

@ -1,3 +0,0 @@
icoUncoupledKinematicParcelDyMFoam.C
EXE = $(FOAM_APPBIN)/icoUncoupledKinematicParcelDyMFoam

View File

@ -1,40 +0,0 @@
EXE_INC = \
-I.. \
-I$(LIB_SRC)/lagrangian/basic/lnInclude \
-I$(LIB_SRC)/lagrangian/intermediate/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
-I$(LIB_SRC)/transportModels/compressible/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/radiation/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \
-I$(LIB_SRC)/transportModels \
-I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/regionModels/regionModel/lnInclude \
-I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \
-I$(LIB_SRC)/dynamicMesh/lnInclude \
-I$(LIB_SRC)/dynamicFvMesh/lnInclude
EXE_LIBS = \
-llagrangian \
-llagrangianIntermediate \
-llagrangianTurbulence \
-lcompressibleTransportModels \
-lfluidThermophysicalModels \
-lspecie \
-lradiationModels \
-lturbulenceModels \
-lincompressibleTurbulenceModels \
-lincompressibleTransportModels \
-lfiniteVolume \
-lfvOptions \
-lmeshTools \
-lregionModels \
-lsurfaceFilmModels \
-ldynamicMesh \
-ldynamicFvMesh \
-ltopoChangerFvMesh

View File

@ -1,95 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation
\\/ 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
uncoupledKinematicParcelDyMFoam
Description
Transient solver for the passive transport of a single kinematic
particle cloud, with optional mesh motion and mesh topology changes.
Uses a pre-calculated velocity field to evolve the cloud.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "dynamicFvMesh.H"
#include "singlePhaseTransportModel.H"
#include "turbulentTransportModel.H"
#include "basicKinematicCollidingCloud.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
argList::addOption
(
"cloudName",
"name",
"specify alternative cloud name. default is 'kinematicCloud'"
);
#include "postProcess.H"
#include "setRootCaseLists.H"
#include "createTime.H"
#include "createDynamicFvMesh.H"
#include "createControl.H"
#include "createFields.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
while (runTime.loop())
{
Info<< "Time = " << runTime.timeName() << nl << endl;
kinematicCloud.storeGlobalPositions();
mesh.update();
U.correctBoundaryConditions();
Info<< "Evolving " << kinematicCloud.name() << endl;
laminarTransport.correct();
mu = laminarTransport.nu()*rhoInfValue;
kinematicCloud.evolve();
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
}
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

View File

@ -22,17 +22,18 @@ License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Application
icoUncoupledKinematicParcelFoam
uncoupledKinematicParcelFoam
Description
Transient solver for the passive transport of a single kinematic
particle cloud.
particle cloud, with optional mesh motion and mesh topology changes.
Uses a pre-calculated velocity field to evolve the cloud.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "dynamicFvMesh.H"
#include "singlePhaseTransportModel.H"
#include "turbulentTransportModel.H"
#include "basicKinematicCollidingCloud.H"
@ -48,13 +49,14 @@ int main(int argc, char *argv[])
"specify alternative cloud name. default is 'kinematicCloud'"
);
#define NO_CONTROL
#include "postProcess.H"
#include "setRootCaseLists.H"
#include "createTime.H"
#include "createMesh.H"
#include "createControl.H"
#include "createDynamicFvMesh.H"
#include "createFields.H"
#include "CourantNo.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -64,10 +66,18 @@ int main(int argc, char *argv[])
{
Info<< "Time = " << runTime.timeName() << nl << endl;
Info<< "Evolving " << kinematicCloud.name() << endl;
kinematicCloud.storeGlobalPositions();
mesh.update();
if (mesh.changing())
{
U.correctBoundaryConditions();
}
laminarTransport.correct();
Info<< "Evolving " << kinematicCloud.name() << endl;
mu = laminarTransport.nu()*rhoInfValue;
kinematicCloud.evolve();

View File

@ -1,6 +1,7 @@
EXE_INC = \
-I. \
-I../reactingParcelFoam \
-I../../compressible/rhoPimpleFoam \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/sampling/lnInclude \
@ -21,7 +22,9 @@ EXE_INC = \
-I$(LIB_SRC)/ODE/lnInclude \
-I$(LIB_SRC)/regionModels/regionModel/lnInclude \
-I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \
-I$(LIB_SRC)/combustionModels/lnInclude
-I$(LIB_SRC)/combustionModels/lnInclude \
-I$(LIB_SRC)/dynamicFvMesh/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude
EXE_LIBS = \
-lturbulenceModels \
@ -43,6 +46,8 @@ EXE_LIBS = \
-lsurfaceFilmModels \
-lcombustionModels \
-lfiniteVolume \
-ldynamicFvMesh \
-ltopoChangerFvMesh \
-lfvOptions \
-lmeshTools \
-lsampling

View File

@ -1,6 +1,5 @@
EXE_INC = \
-I. \
-I../sprayDyMFoam \
-I.. \
-I../../reactingParcelFoam \
-I$(LIB_SRC)/finiteVolume/lnInclude \

View File

@ -57,10 +57,10 @@ int main(int argc, char *argv[])
#include "readEngineTimeControls.H"
#include "createFields.H"
#include "createFieldRefs.H"
#include "createRhoUf.H"
#include "compressibleCourantNo.H"
#include "setInitialDeltaT.H"
#include "initContinuityErrs.H"
#include "createRhoUfIfPresent.H"
#include "startSummary.H"
turbulence->validate();

View File

@ -22,11 +22,12 @@ if (pimple.transonic())
fvc::flux(HbyA)
+ MRF.zeroFilter
(
rhorAUf*fvc::ddtCorr(rho, U, phi)/fvc::interpolate(rho)
rhorAUf*fvc::ddtCorr(rho, U, phi, rhoUf)/fvc::interpolate(rho)
)
)
);
fvc::makeRelative(phid, psi, U);
MRF.makeRelative(fvc::interpolate(psi), phid);
while (pimple.correctNonOrthogonal())
@ -54,12 +55,11 @@ else
surfaceScalarField phiHbyA
(
"phiHbyA",
(
fvc::flux(rho*HbyA)
+ rhorAUf*fvc::ddtCorr(rho, U, phi)
)
fvc::flux(rho*HbyA)
+ MRF.zeroFilter(rhorAUf*fvc::ddtCorr(rho, U, phi, rhoUf))
);
fvc::makeRelative(phiHbyA, rho, U);
MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
// Update the pressure BCs to ensure flux consistency
@ -105,7 +105,15 @@ U.correctBoundaryConditions();
fvOptions.correct(U);
K = 0.5*magSqr(U);
// Correct rhoUf if the mesh is moving
fvc::correctRhoUf(rhoUf, rho, U, phi);
if (thermo.dpdt())
{
dpdt = fvc::ddt(p);
if (mesh.moving())
{
dpdt -= fvc::div(fvc::meshPhi(rho, U), p);
}
}

View File

@ -1,3 +0,0 @@
sprayDyMFoam.C
EXE = $(FOAM_APPBIN)/sprayDyMFoam

View File

@ -1,54 +0,0 @@
EXE_INC = \
-I.. \
-I../../reactingParcelFoam \
-I../../../compressible/rhoPimpleFoam \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/sampling/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \
-I$(LIB_SRC)/lagrangian/basic/lnInclude \
-I$(LIB_SRC)/lagrangian/intermediate/lnInclude \
-I$(LIB_SRC)/lagrangian/spray/lnInclude \
-I$(LIB_SRC)/lagrangian/distributionModels/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
-I$(LIB_SRC)/transportModels/compressible/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/thermophysicalProperties/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/SLGThermo/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/chemistryModel/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/radiation/lnInclude \
-I$(LIB_SRC)/ODE/lnInclude \
-I$(LIB_SRC)/regionModels/regionModel/lnInclude \
-I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \
-I$(LIB_SRC)/combustionModels/lnInclude \
-I$(LIB_SRC)/dynamicFvMesh/lnInclude \
-I$(LIB_SRC)/dynamicMesh/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude
EXE_LIBS = \
-lturbulenceModels \
-lcompressibleTurbulenceModels \
-llagrangian \
-llagrangianIntermediate \
-llagrangianTurbulence \
-llagrangianSpray \
-lspecie \
-lcompressibleTransportModels \
-lfluidThermophysicalModels \
-lthermophysicalProperties \
-lreactionThermophysicalModels \
-lSLGThermo \
-lchemistryModel \
-lradiationModels \
-lODE \
-lregionModels \
-lsurfaceFilmModels \
-lcombustionModels \
-lfiniteVolume \
-lfvOptions \
-ldynamicFvMesh \
-ltopoChangerFvMesh \
-ldynamicMesh \
-lmeshTools

View File

@ -1,122 +0,0 @@
rho = thermo.rho();
rho = max(rho, rhoMin);
rho = min(rho, rhoMax);
rho.relax();
volScalarField rAU(1.0/UEqn.A());
surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rho*rAU));
volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p));
if (pimple.nCorrPISO() <= 1)
{
tUEqn.clear();
}
if (pimple.transonic())
{
surfaceScalarField phid
(
"phid",
fvc::interpolate(psi)
*(
fvc::flux(HbyA)
+ MRF.zeroFilter
(
rhorAUf*fvc::ddtCorr(rho, U, rhoUf)/fvc::interpolate(rho)
)
)
);
fvc::makeRelative(phid, psi, U);
MRF.makeRelative(fvc::interpolate(psi), phid);
while (pimple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
(
fvm::ddt(psi, p)
+ fvm::div(phid, p)
- fvm::laplacian(rhorAUf, p)
==
parcels.Srho()
+ fvOptions(psi, p, rho.name())
);
pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter())));
if (pimple.finalNonOrthogonalIter())
{
phi == pEqn.flux();
}
}
}
else
{
surfaceScalarField phiHbyA
(
"phiHbyA",
fvc::flux(rho*HbyA)
+ MRF.zeroFilter(rhorAUf*fvc::ddtCorr(rho, U, rhoUf))
);
fvc::makeRelative(phiHbyA, rho, U);
MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
// Update the pressure BCs to ensure flux consistency
constrainPressure(p, rho, U, phiHbyA, rhorAUf, MRF);
while (pimple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
(
fvm::ddt(psi, p)
+ fvc::div(phiHbyA)
- fvm::laplacian(rhorAUf, p)
==
parcels.Srho()
+ fvOptions(psi, p, rho.name())
);
pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter())));
if (pimple.finalNonOrthogonalIter())
{
phi = phiHbyA + pEqn.flux();
}
}
}
#include "rhoEqn.H"
#include "compressibleContinuityErrs.H"
// Explicitly relax pressure for momentum corrector
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;
U = HbyA - rAU*fvc::grad(p);
U.correctBoundaryConditions();
fvOptions.correct(U);
K = 0.5*magSqr(U);
{
rhoUf = fvc::interpolate(rho*U);
surfaceVectorField n(mesh.Sf()/mesh.magSf());
rhoUf += n*(fvc::absolute(phi, rho, U)/mesh.magSf() - (n & rhoUf));
}
if (thermo.dpdt())
{
dpdt = fvc::ddt(p);
if (mesh.moving())
{
dpdt -= fvc::div(fvc::meshPhi(rho, U), p);
}
}

View File

@ -1,159 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2015-2018 OpenFOAM Foundation
\\/ 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
sprayDyMFoam
Description
Transient solver for compressible, turbulent flow with a spray particle
cloud, with optional mesh motion and mesh topology changes.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "dynamicFvMesh.H"
#include "turbulenceModel.H"
#include "basicSprayCloud.H"
#include "psiReactionThermo.H"
#include "CombustionModel.H"
#include "radiationModel.H"
#include "SLGThermo.H"
#include "pimpleControl.H"
#include "CorrectPhi.H"
#include "fvOptions.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
#include "postProcess.H"
#include "setRootCaseLists.H"
#include "createTime.H"
#include "createDynamicFvMesh.H"
#include "createDyMControls.H"
#include "createFields.H"
#include "createFieldRefs.H"
#include "createRhoUf.H"
#include "compressibleCourantNo.H"
#include "setInitialDeltaT.H"
#include "initContinuityErrs.H"
turbulence->validate();
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
while (runTime.run())
{
#include "readDyMControls.H"
{
// Store divrhoU from the previous time-step/mesh for the correctPhi
volScalarField divrhoU
(
"divrhoU",
fvc::div(fvc::absolute(phi, rho, U))
);
#include "compressibleCourantNo.H"
#include "setDeltaT.H"
runTime++;
Info<< "Time = " << runTime.timeName() << nl << endl;
// Store momentum to set rhoUf for introduced faces.
volVectorField rhoU("rhoU", rho*U);
// Store the particle positions
parcels.storeGlobalPositions();
// Do any mesh changes
mesh.update();
if (mesh.changing())
{
MRF.update();
if (correctPhi)
{
// Calculate absolute flux from the mapped surface velocity
phi = mesh.Sf() & rhoUf;
#include "correctPhi.H"
// Make the fluxes relative to the mesh-motion
fvc::makeRelative(phi, rho, U);
}
if (checkMeshCourantNo)
{
#include "meshCourantNo.H"
}
}
}
parcels.evolve();
#include "rhoEqn.H"
// --- Pressure-velocity PIMPLE corrector loop
while (pimple.loop())
{
#include "UEqn.H"
#include "YEqn.H"
#include "EEqn.H"
// --- Pressure corrector loop
while (pimple.correct())
{
#include "pEqn.H"
}
if (pimple.turbCorr())
{
turbulence->correct();
}
}
rho = thermo.rho();
if (runTime.write())
{
combustion->Qdot()().write();
}
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
}
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

View File

@ -26,18 +26,20 @@ Application
Description
Transient solver for compressible, turbulent flow with a spray particle
cloud.
cloud, with optional mesh motion and mesh topology changes.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "turbulentFluidThermoModel.H"
#include "dynamicFvMesh.H"
#include "turbulenceModel.H"
#include "basicSprayCloud.H"
#include "psiReactionThermo.H"
#include "CombustionModel.H"
#include "radiationModel.H"
#include "SLGThermo.H"
#include "pimpleControl.H"
#include "CorrectPhi.H"
#include "fvOptions.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -48,14 +50,14 @@ int main(int argc, char *argv[])
#include "setRootCaseLists.H"
#include "createTime.H"
#include "createMesh.H"
#include "createControl.H"
#include "createTimeControls.H"
#include "createDynamicFvMesh.H"
#include "createDyMControls.H"
#include "createFields.H"
#include "createFieldRefs.H"
#include "compressibleCourantNo.H"
#include "setInitialDeltaT.H"
#include "initContinuityErrs.H"
#include "createRhoUfIfPresent.H"
turbulence->validate();
@ -65,7 +67,21 @@ int main(int argc, char *argv[])
while (runTime.run())
{
#include "readTimeControls.H"
#include "readDyMControls.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 = new volScalarField
(
"divrhoU",
fvc::div(fvc::absolute(phi, rho, U))
);
}
#include "compressibleCourantNo.H"
#include "setDeltaT.H"
@ -73,44 +89,68 @@ int main(int argc, char *argv[])
Info<< "Time = " << runTime.timeName() << nl << endl;
parcels.evolve();
if (!pimple.frozenFlow())
// Store momentum to set rhoUf for introduced faces.
autoPtr<volVectorField> rhoU;
if (rhoUf.valid())
{
#include "rhoEqn.H"
rhoU = new volVectorField("rhoU", rho*U);
}
// --- Pressure-velocity PIMPLE corrector loop
while (pimple.loop())
// Store the particle positions
parcels.storeGlobalPositions();
// Do any mesh changes
mesh.update();
if (mesh.changing())
{
MRF.update();
if (correctPhi)
{
#include "UEqn.H"
#include "YEqn.H"
#include "EEqn.H"
// Calculate absolute flux from the mapped surface velocity
phi = mesh.Sf() & rhoUf();
// --- Pressure corrector loop
while (pimple.correct())
{
#include "pEqn.H"
}
#include "correctPhi.H"
if (pimple.turbCorr())
{
turbulence->correct();
}
// Make the fluxes relative to the mesh-motion
fvc::makeRelative(phi, rho, U);
}
rho = thermo.rho();
if (runTime.write())
if (checkMeshCourantNo)
{
combustion->Qdot()().write();
#include "meshCourantNo.H"
}
}
else
parcels.evolve();
#include "rhoEqn.H"
// --- Pressure-velocity PIMPLE corrector loop
while (pimple.loop())
{
if (runTime.writeTime())
#include "UEqn.H"
#include "YEqn.H"
#include "EEqn.H"
// --- Pressure corrector loop
while (pimple.correct())
{
parcels.write();
#include "pEqn.H"
}
if (pimple.turbCorr())
{
turbulence->correct();
}
}
rho = thermo.rho();
if (runTime.write())
{
combustion->Qdot()().write();
}
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"

View File

@ -11,7 +11,8 @@ EXE_INC = \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/regionModels/regionModel/lnInclude \
-I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude
-I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \
-I$(LIB_SRC)/dynamicFvMesh/lnInclude
EXE_LIBS = \
-llagrangian \
@ -27,4 +28,5 @@ EXE_LIBS = \
-lfvOptions \
-lmeshTools \
-lregionModels \
-lsurfaceFilmModels
-lsurfaceFilmModels \
-ldynamicFvMesh

View File

@ -1,3 +0,0 @@
uncoupledKinematicParcelDyMFoam.C
EXE = $(FOAM_APPBIN)/uncoupledKinematicParcelDyMFoam

View File

@ -1,36 +0,0 @@
EXE_INC = \
-I.. \
-I$(LIB_SRC)/lagrangian/basic/lnInclude \
-I$(LIB_SRC)/lagrangian/intermediate/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
-I$(LIB_SRC)/transportModels/compressible/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/radiation/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/regionModels/regionModel/lnInclude \
-I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \
-I$(LIB_SRC)/dynamicMesh/lnInclude \
-I$(LIB_SRC)/dynamicFvMesh/lnInclude
EXE_LIBS = \
-llagrangian \
-llagrangianIntermediate \
-llagrangianTurbulence \
-lcompressibleTransportModels \
-lfluidThermophysicalModels \
-lspecie \
-lradiationModels \
-lturbulenceModels \
-lcompressibleTurbulenceModels \
-lfiniteVolume \
-lfvOptions \
-lmeshTools \
-lregionModels \
-lsurfaceFilmModels \
-ldynamicMesh \
-ldynamicFvMesh \
-ltopoChangerFvMesh

View File

@ -1,90 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2017-2018 OpenFOAM Foundation
\\/ 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
uncoupledKinematicParcelDyMFoam
Description
Transient solver for the passive transport of a particle cloud.
Uses a pre- calculated velocity field to evolve the cloud.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "dynamicFvMesh.H"
#include "psiThermo.H"
#include "turbulentFluidThermoModel.H"
#include "basicKinematicCloud.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
argList::addOption
(
"cloudName",
"name",
"specify alternative cloud name. default is 'kinematicCloud'"
);
#define NO_CONTROL
#include "postProcess.H"
#include "setRootCaseLists.H"
#include "createTime.H"
#include "createDynamicFvMesh.H"
#include "createFields.H"
#include "compressibleCourantNo.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
while (runTime.loop())
{
Info<< "Time = " << runTime.timeName() << nl << endl;
kinematicCloud.storeGlobalPositions();
mesh.update();
U.correctBoundaryConditions();
Info<< "Evolving " << kinematicCloud.name() << endl;
kinematicCloud.evolve();
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
}
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

View File

@ -32,6 +32,7 @@ Description
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "dynamicFvMesh.H"
#include "psiThermo.H"
#include "turbulentFluidThermoModel.H"
#include "basicKinematicCloud.H"
@ -52,7 +53,7 @@ int main(int argc, char *argv[])
#include "setRootCaseLists.H"
#include "createTime.H"
#include "createMesh.H"
#include "createDynamicFvMesh.H"
#include "createFields.H"
#include "compressibleCourantNo.H"
@ -64,6 +65,15 @@ int main(int argc, char *argv[])
{
Info<< "Time = " << runTime.timeName() << nl << endl;
kinematicCloud.storeGlobalPositions();
mesh.update();
if (mesh.changing())
{
U.correctBoundaryConditions();
}
Info<< "Evolving " << kinematicCloud.name() << endl;
kinematicCloud.evolve();

View File

@ -7,6 +7,7 @@ EXE_INC = \
-I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/barotropicCompressibilityModel/lnInclude \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/dynamicFvMesh/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude
EXE_LIBS = \
@ -14,5 +15,6 @@ EXE_LIBS = \
-lincompressibleTurbulenceModels \
-lbarotropicCompressibilityModel \
-lfiniteVolume \
-ldynamicFvMesh \
-lfvOptions \
-lmeshTools

View File

@ -1,3 +0,0 @@
cavitatingDyMFoam.C
EXE = $(FOAM_APPBIN)/cavitatingDyMFoam

View File

@ -1,23 +0,0 @@
EXE_INC = \
-I.. \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/transportModels/twoPhaseMixture/lnInclude \
-I$(LIB_SRC)/transportModels \
-I$(LIB_SRC)/transportModels/incompressible/lnInclude \
-I$(LIB_SRC)/transportModels/interfaceProperties/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/barotropicCompressibilityModel/lnInclude \
-I$(LIB_SRC)/dynamicMesh/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/dynamicFvMesh/lnInclude
EXE_LIBS = \
-lturbulenceModels \
-lincompressibleTurbulenceModels \
-lbarotropicCompressibilityModel \
-lfiniteVolume \
-lfvOptions \
-ldynamicMesh \
-ldynamicFvMesh \
-lmeshTools

View File

@ -1,125 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2012-2018 OpenFOAM Foundation
\\/ 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
cavitatingFoam
Description
Transient cavitation code based on the homogeneous equilibrium model
from which the compressibility of the liquid/vapour "mixture" is obtained,
with optional mesh motion and mesh topology changes.
Turbulence modelling is generic, i.e. laminar, RAS or LES may be selected.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "dynamicFvMesh.H"
#include "barotropicCompressibilityModel.H"
#include "incompressibleTwoPhaseMixture.H"
#include "turbulentTransportModel.H"
#include "CorrectPhi.H"
#include "pimpleControl.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
#include "postProcess.H"
#include "setRootCaseLists.H"
#include "createTime.H"
#include "createDynamicFvMesh.H"
#include "createControls.H"
#include "createFields.H"
#include "createUf.H"
#include "createPcorrTypes.H"
#include "CourantNo.H"
#include "setInitialDeltaT.H"
turbulence->validate();
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
while (runTime.run())
{
#include "readControls.H"
{
#include "CourantNo.H"
#include "setDeltaT.H"
runTime++;
Info<< "Time = " << runTime.timeName() << nl << endl;
// Do any mesh changes
mesh.update();
if (mesh.changing() && correctPhi)
{
// Calculate absolute flux from the mapped surface velocity
phi = mesh.Sf() & Uf;
#include "correctPhi.H"
// Make the flux relative to the mesh motion
fvc::makeRelative(phi, U);
}
}
// --- Pressure-velocity PIMPLE corrector loop
while (pimple.loop())
{
#include "rhoEqn.H"
#include "alphavPsi.H"
#include "UEqn.H"
// --- Pressure corrector loop
while (pimple.correct())
{
#include "pEqn.H"
}
if (pimple.turbCorr())
{
turbulence->correct();
}
}
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
}
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

View File

@ -1,6 +0,0 @@
#include "createDyMControls.H"
scalar maxAcousticCo
(
readScalar(runTime.controlDict().lookup("maxAcousticCo"))
);

View File

@ -1,88 +0,0 @@
{
if (pimple.nCorrPimple() == 1)
{
p =
(
rho
- alphal*rhol0
- ((alphav*psiv + alphal*psil) - psi)*pSat
)/psi;
}
surfaceScalarField rhof("rhof", fvc::interpolate(rho));
volScalarField rAU(1.0/UEqn.A());
surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rho*rAU));
volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p));
phi = fvc::flux(HbyA)
+ rhorAUf*fvc::ddtCorr(U, Uf);
fvc::makeRelative(phi, U);
surfaceScalarField phiGradp(rhorAUf*mesh.magSf()*fvc::snGrad(p));
phi -= phiGradp/rhof;
volScalarField rho0(rho - psi*p);
while (pimple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
(
fvc::ddt(rho)
+ psi*correction(fvm::ddt(p))
+ fvc::div(phi, rho)
+ fvc::div(phiGradp)
- fvm::laplacian(rhorAUf, p)
);
pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter())));
if (pimple.finalNonOrthogonalIter())
{
phi += (phiGradp + pEqn.flux())/rhof;
}
}
Info<< "Predicted p max-min : " << max(p).value()
<< " " << min(p).value() << endl;
rho == max(rho0 + psi*p, rhoMin);
#include "alphavPsi.H"
p =
(
rho
- alphal*rhol0
- ((alphav*psiv + alphal*psil) - psi)*pSat
)/psi;
p.correctBoundaryConditions();
Info<< "Phase-change corrected p max-min : " << max(p).value()
<< " " << min(p).value() << endl;
// Correct velocity
U = HbyA - rAU*fvc::grad(p);
// Remove the swirl component of velocity for "wedge" cases
if (pimple.dict().found("removeSwirl"))
{
label swirlCmpt(readLabel(pimple.dict().lookup("removeSwirl")));
Info<< "Removing swirl component-" << swirlCmpt << " of U" << endl;
U.field().replace(swirlCmpt, 0.0);
}
U.correctBoundaryConditions();
Info<< "max(U) " << max(mag(U)).value() << endl;
{
Uf = fvc::interpolate(U);
surfaceVectorField n(mesh.Sf()/mesh.magSf());
Uf += n*(phi/mesh.magSf() - (n & Uf));
}
}

View File

@ -1,3 +0,0 @@
#include "readDyMControls.H"
maxAcousticCo = readScalar(runTime.controlDict().lookup("maxAcousticCo"));

View File

@ -26,16 +26,19 @@ Application
Description
Transient cavitation code based on the homogeneous equilibrium model
from which the compressibility of the liquid/vapour "mixture" is obtained.
from which the compressibility of the liquid/vapour "mixture" is obtained,
with optional mesh motion and mesh topology changes.
Turbulence modelling is generic, i.e. laminar, RAS or LES may be selected.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "dynamicFvMesh.H"
#include "barotropicCompressibilityModel.H"
#include "incompressibleTwoPhaseMixture.H"
#include "turbulentTransportModel.H"
#include "CorrectPhi.H"
#include "pimpleControl.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -46,10 +49,11 @@ int main(int argc, char *argv[])
#include "setRootCaseLists.H"
#include "createTime.H"
#include "createMesh.H"
#include "createControl.H"
#include "createDynamicFvMesh.H"
#include "createControls.H"
#include "createFields.H"
#include "createUfIfPresent.H"
#include "createPcorrTypes.H"
#include "CourantNo.H"
#include "setInitialDeltaT.H"
@ -62,11 +66,29 @@ int main(int argc, char *argv[])
while (runTime.run())
{
#include "readControls.H"
#include "CourantNo.H"
#include "setDeltaT.H"
runTime++;
Info<< "Time = " << runTime.timeName() << nl << endl;
{
#include "CourantNo.H"
#include "setDeltaT.H"
runTime++;
Info<< "Time = " << runTime.timeName() << nl << endl;
// Do any mesh changes
mesh.update();
if (mesh.changing() && correctPhi)
{
// Calculate absolute flux from the mapped surface velocity
phi = mesh.Sf() & Uf();
#include "correctPhi.H"
// Make the flux relative to the mesh motion
fvc::makeRelative(phi, U);
}
}
// --- Pressure-velocity PIMPLE corrector loop
while (pimple.loop())

View File

@ -1,4 +1,4 @@
#include "createTimeControls.H"
#include "createDyMControls.H"
scalar maxAcousticCo
(

View File

@ -15,19 +15,21 @@
surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rho*rAU));
volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p));
phi = fvc::flux(HbyA)
+ rhorAUf*fvc::ddtCorr(U, phi);
phi = fvc::flux(HbyA) + rhorAUf*fvc::ddtCorr(U, phi, Uf);
fvc::makeRelative(phi, U);
surfaceScalarField phiGradp(rhorAUf*mesh.magSf()*fvc::snGrad(p));
phi -= phiGradp/rhof;
volScalarField rho0(rho - psi*p);
while (pimple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
(
fvm::ddt(psi, p)
- (rhol0 + (psil - psiv)*pSat)*fvc::ddt(alphav) - pSat*fvc::ddt(psi)
fvc::ddt(rho)
+ psi*correction(fvm::ddt(p))
+ fvc::div(phi, rho)
+ fvc::div(phiGradp)
- fvm::laplacian(rhorAUf, p)
@ -44,13 +46,7 @@
Info<< "Predicted p max-min : " << max(p).value()
<< " " << min(p).value() << endl;
rho == max
(
psi*p
+ alphal*rhol0
+ ((alphav*psiv + alphal*psil) - psi)*pSat,
rhoMin
);
rho == max(rho0 + psi*p, rhoMin);
#include "alphavPsi.H"
@ -82,4 +78,7 @@
U.correctBoundaryConditions();
Info<< "max(U) " << max(mag(U)).value() << endl;
// Correct Uf if the mesh is moving
fvc::correctUf(Uf, U, phi);
}

View File

@ -1,3 +1,3 @@
#include "readTimeControls.H"
#include "readDyMControls.H"
maxAcousticCo = readScalar(runTime.controlDict().lookup("maxAcousticCo"));

View File

@ -6,7 +6,6 @@ wclean libso surfaceTensionModels
wclean libso VoFphaseCompressibleTurbulenceModels
wclean
wclean compressibleInterDyMFoam
wclean compressibleInterFilmFoam
#------------------------------------------------------------------------------

View File

@ -9,7 +9,6 @@ wmake $targetType surfaceTensionModels
wmake $targetType VoFphaseCompressibleTurbulenceModels
wmake $targetType
wmake $targetType compressibleInterDyMFoam
wmake $targetType compressibleInterFilmFoam
#------------------------------------------------------------------------------

View File

@ -2,6 +2,7 @@ EXE_INC = \
-I. \
-I../VoF \
-ItwoPhaseMixtureThermo \
-IVoFphaseCompressibleTurbulenceModels/lnInclude \
-I$(LIB_SRC)/transportModels/compressible/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/transportModels/twoPhaseMixture/lnInclude \
@ -9,9 +10,10 @@ EXE_INC = \
-I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/phaseCompressible/lnInclude \
-IVoFphaseCompressibleTurbulenceModels/lnInclude \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/dynamicMesh/lnInclude \
-I$(LIB_SRC)/dynamicFvMesh/lnInclude \
-I$(LIB_SRC)/finiteVolume/lnInclude
EXE_LIBS = \
-ltwoPhaseMixtureThermo \
@ -27,4 +29,6 @@ EXE_LIBS = \
-lVoFphaseCompressibleTurbulenceModels \
-lfiniteVolume \
-lfvOptions \
-lmeshTools
-lmeshTools \
-ldynamicMesh \
-ldynamicFvMesh

View File

@ -4,7 +4,7 @@
fvm::ddt(rho, T) + fvm::div(rhoPhi, T) - fvm::Sp(contErr, T)
- fvm::laplacian(turbulence.alphaEff(), T)
+ (
fvc::div(fvc::absolute(phi, U), p)()() - contErr/rho*p
fvc::div(fvc::absolute(phi, U), p)()() // - contErr/rho*p
+ (fvc::ddt(rho, K) + fvc::div(rhoPhi, K))()() - contErr*K
)
*(

View File

@ -1,3 +0,0 @@
compressibleInterDyMFoam.C
EXE = $(FOAM_APPBIN)/compressibleInterDyMFoam

View File

@ -1,35 +0,0 @@
EXE_INC = \
-I. \
-I.. \
-I../../VoF \
-I../twoPhaseMixtureThermo \
-I../VoFphaseCompressibleTurbulenceModels/lnInclude \
-I$(LIB_SRC)/transportModels/compressible/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/transportModels/twoPhaseMixture/lnInclude \
-I$(LIB_SRC)/transportModels/interfaceProperties/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/phaseCompressible/lnInclude \
-I$(LIB_SRC)/dynamicMesh/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/dynamicFvMesh/lnInclude \
-I$(LIB_SRC)/finiteVolume/lnInclude
EXE_LIBS = \
-ltwoPhaseMixtureThermo \
-ltwoPhaseSurfaceTension \
-lcompressibleTransportModels \
-lfluidThermophysicalModels \
-lspecie \
-ltwoPhaseMixture \
-ltwoPhaseProperties \
-linterfaceProperties \
-lturbulenceModels \
-lcompressibleTurbulenceModels \
-lVoFphaseCompressibleTurbulenceModels \
-ldynamicMesh \
-lmeshTools \
-ldynamicFvMesh \
-lfiniteVolume \
-lfvOptions

View File

@ -1,176 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation
\\/ 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
compressibleInterDyMFoam
Description
Solver for 2 compressible, non-isothermal immiscible fluids using a VOF
(volume of fluid) phase-fraction based interface capturing approach,
with optional mesh motion and mesh topology changes including adaptive
re-meshing.
The momentum and other fluid properties are of the "mixture" and a single
momentum equation is solved.
Turbulence modelling is generic, i.e. laminar, RAS or LES may be selected.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "dynamicFvMesh.H"
#include "CMULES.H"
#include "EulerDdtScheme.H"
#include "localEulerDdtScheme.H"
#include "CrankNicolsonDdtScheme.H"
#include "subCycle.H"
#include "compressibleInterPhaseTransportModel.H"
#include "pimpleControl.H"
#include "fvOptions.H"
#include "CorrectPhi.H"
#include "fvcSmooth.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
#include "postProcess.H"
#include "setRootCaseLists.H"
#include "createTime.H"
#include "createDynamicFvMesh.H"
#include "initContinuityErrs.H"
#include "createDyMControls.H"
#include "createFields.H"
#include "createUf.H"
#include "CourantNo.H"
#include "setInitialDeltaT.H"
volScalarField& p = mixture.p();
volScalarField& T = mixture.T();
const volScalarField& psi1 = mixture.thermo1().psi();
const volScalarField& psi2 = mixture.thermo2().psi();
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
while (runTime.run())
{
#include "readDyMControls.H"
// Store divU from the previous mesh so that it can be mapped
// and used in correctPhi to ensure the corrected phi has the
// same divergence
volScalarField divU("divU0", fvc::div(fvc::absolute(phi, U)));
if (LTS)
{
#include "setRDeltaT.H"
}
else
{
#include "CourantNo.H"
#include "alphaCourantNo.H"
#include "setDeltaT.H"
}
runTime++;
Info<< "Time = " << runTime.timeName() << nl << endl;
// --- Pressure-velocity PIMPLE corrector loop
while (pimple.loop())
{
if (pimple.firstIter() || moveMeshOuterCorrectors)
{
scalar timeBeforeMeshUpdate = runTime.elapsedCpuTime();
mesh.update();
if (mesh.changing())
{
MRF.update();
Info<< "Execution time for mesh.update() = "
<< runTime.elapsedCpuTime() - timeBeforeMeshUpdate
<< " s" << endl;
gh = (g & mesh.C()) - ghRef;
ghf = (g & mesh.Cf()) - ghRef;
if (correctPhi)
{
// Calculate absolute flux
// from the mapped surface velocity
phi = mesh.Sf() & Uf;
#include "correctPhi.H"
// Make the fluxes relative to the mesh motion
fvc::makeRelative(phi, U);
mixture.correct();
}
if (checkMeshCourantNo)
{
#include "meshCourantNo.H"
}
}
}
#include "alphaControls.H"
#include "compressibleAlphaEqnSubCycle.H"
turbulence.correctPhasePhi();
#include "UEqn.H"
#include "TEqn.H"
// --- Pressure corrector loop
while (pimple.correct())
{
#include "pEqn.H"
}
if (pimple.turbCorr())
{
turbulence.correct();
}
}
runTime.write();
Info<< "ExecutionTime = "
<< runTime.elapsedCpuTime()
<< " s\n\n" << endl;
}
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

View File

@ -1,156 +0,0 @@
{
volScalarField rAU("rAU", 1.0/UEqn.A());
surfaceScalarField rAUf("rAUf", fvc::interpolate(rAU));
volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p_rgh));
surfaceScalarField phiHbyA
(
"phiHbyA",
fvc::flux(HbyA)
+ MRF.zeroFilter(fvc::interpolate(rho*rAU)*fvc::ddtCorr(U, Uf))
);
MRF.makeRelative(phiHbyA);
surfaceScalarField phig
(
(
mixture.surfaceTensionForce()
- ghf*fvc::snGrad(rho)
)*rAUf*mesh.magSf()
);
phiHbyA += phig;
// Update the pressure BCs to ensure flux consistency
constrainPressure(p_rgh, U, phiHbyA, rAUf, MRF);
// Make the fluxes relative to the mesh motion
fvc::makeRelative(phiHbyA, U);
tmp<fvScalarMatrix> p_rghEqnComp1;
tmp<fvScalarMatrix> p_rghEqnComp2;
if (pimple.transonic())
{
#include "rhofs.H"
surfaceScalarField phid1("phid1", fvc::interpolate(psi1)*phi);
surfaceScalarField phid2("phid2", fvc::interpolate(psi2)*phi);
p_rghEqnComp1 =
pos(alpha1)
*(
(
fvc::ddt(alpha1, rho1) + fvc::div(alphaPhi1*rho1f)
- (fvOptions(alpha1, mixture.thermo1().rho())&rho1)
)/rho1
- fvc::ddt(alpha1) - fvc::div(alphaPhi1)
+ (alpha1/rho1)
*correction
(
psi1*fvm::ddt(p_rgh)
+ fvm::div(phid1, p_rgh) - fvm::Sp(fvc::div(phid1), p_rgh)
)
);
p_rghEqnComp1.ref().relax();
p_rghEqnComp2 =
pos(alpha2)
*(
(
fvc::ddt(alpha2, rho2) + fvc::div(alphaPhi2*rho2f)
- (fvOptions(alpha2, mixture.thermo2().rho())&rho2)
)/rho2
- fvc::ddt(alpha2) - fvc::div(alphaPhi2)
+ (alpha2/rho2)
*correction
(
psi2*fvm::ddt(p_rgh)
+ fvm::div(phid2, p_rgh) - fvm::Sp(fvc::div(phid2), p_rgh)
)
);
p_rghEqnComp2.ref().relax();
}
else
{
#include "rhofs.H"
p_rghEqnComp1 =
pos(alpha1)
*(
(
fvc::ddt(alpha1, rho1) + fvc::div(alphaPhi1*rho1f)
- (fvOptions(alpha1, mixture.thermo1().rho())&rho1)
)/rho1
- fvc::ddt(alpha1) - fvc::div(alphaPhi1)
+ fvc::div(mesh.phi())*alpha1
+ (alpha1*psi1/rho1)*correction(fvm::ddt(p_rgh))
);
p_rghEqnComp2 =
pos(alpha2)
*(
(
fvc::ddt(alpha2, rho2) + fvc::div(alphaPhi2*rho2f)
- (fvOptions(alpha2, mixture.thermo2().rho())&rho2)
)/rho2
- fvc::ddt(alpha2) - fvc::div(alphaPhi2)
+ fvc::div(mesh.phi())*alpha2
+ (alpha2*psi2/rho2)*correction(fvm::ddt(p_rgh))
);
}
// Cache p_rgh prior to solve for density update
volScalarField p_rgh_0(p_rgh);
while (pimple.correctNonOrthogonal())
{
fvScalarMatrix p_rghEqnIncomp
(
fvc::div(phiHbyA)
- fvm::laplacian(rAUf, p_rgh)
);
solve
(
p_rghEqnComp1() + p_rghEqnComp2() + p_rghEqnIncomp,
mesh.solver(p_rgh.select(pimple.finalInnerIter()))
);
if (pimple.finalNonOrthogonalIter())
{
p = max(p_rgh + (alpha1*rho1 + alpha2*rho2)*gh, pMin);
p_rgh = p - (alpha1*rho1 + alpha2*rho2)*gh;
dgdt =
(
alpha1*(p_rghEqnComp2 & p_rgh)
- alpha2*(p_rghEqnComp1 & p_rgh)
);
phi = phiHbyA + p_rghEqnIncomp.flux();
U = HbyA
+ rAU*fvc::reconstruct((phig + p_rghEqnIncomp.flux())/rAUf);
U.correctBoundaryConditions();
fvOptions.correct(U);
}
}
{
Uf = fvc::interpolate(U);
surfaceVectorField n(mesh.Sf()/mesh.magSf());
Uf += n*(fvc::absolute(phi, U)/mesh.magSf() - (n & Uf));
}
// Update densities from change in p_rgh
mixture.thermo1().correctRho(psi1*(p_rgh - p_rgh_0));
mixture.thermo2().correctRho(psi2*(p_rgh - p_rgh_0));
rho = alpha1*rho1 + alpha2*rho2;
// Correct p_rgh for consistency with p and the updated densities
p_rgh = p - rho*gh;
p_rgh.correctBoundaryConditions();
K = 0.5*magSqr(U);
}

View File

@ -26,7 +26,9 @@ Application
Description
Solver for 2 compressible, non-isothermal immiscible fluids using a VOF
(volume of fluid) phase-fraction based interface capturing approach.
(volume of fluid) phase-fraction based interface capturing approach,
with optional mesh motion and mesh topology changes including adaptive
re-meshing.
The momentum and other fluid properties are of the "mixture" and a single
momentum equation is solved.
@ -39,6 +41,7 @@ Description
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "dynamicFvMesh.H"
#include "CMULES.H"
#include "EulerDdtScheme.H"
#include "localEulerDdtScheme.H"
@ -47,6 +50,7 @@ Description
#include "compressibleInterPhaseTransportModel.H"
#include "pimpleControl.H"
#include "fvOptions.H"
#include "CorrectPhi.H"
#include "fvcSmooth.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -57,30 +61,30 @@ int main(int argc, char *argv[])
#include "setRootCaseLists.H"
#include "createTime.H"
#include "createMesh.H"
#include "createControl.H"
#include "createTimeControls.H"
#include "createDynamicFvMesh.H"
#include "initContinuityErrs.H"
#include "createDyMControls.H"
#include "createFields.H"
#include "CourantNo.H"
#include "setInitialDeltaT.H"
#include "createUfIfPresent.H"
volScalarField& p = mixture.p();
volScalarField& T = mixture.T();
const volScalarField& psi1 = mixture.thermo1().psi();
const volScalarField& psi2 = mixture.thermo2().psi();
if (!LTS)
{
#include "readTimeControls.H"
#include "CourantNo.H"
#include "setInitialDeltaT.H"
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
while (runTime.run())
{
#include "readTimeControls.H"
#include "readDyMControls.H"
// Store divU from the previous mesh so that it can be mapped
// and used in correctPhi to ensure the corrected phi has the
// same divergence
volScalarField divU("divU0", fvc::div(fvc::absolute(phi, U)));
if (LTS)
{
@ -100,6 +104,40 @@ int main(int argc, char *argv[])
// --- Pressure-velocity PIMPLE corrector loop
while (pimple.loop())
{
if (pimple.firstIter() || moveMeshOuterCorrectors)
{
mesh.update();
if (mesh.changing())
{
gh = (g & mesh.C()) - ghRef;
ghf = (g & mesh.Cf()) - ghRef;
MRF.update();
if (correctPhi)
{
// Calculate absolute flux
// from the mapped surface velocity
phi = mesh.Sf() & Uf();
#include "correctPhi.H"
// Make the fluxes relative to the mesh motion
fvc::makeRelative(phi, U);
mixture.correct();
}
if (checkMeshCourantNo)
{
#include "meshCourantNo.H"
}
}
}
divU = fvc::div(fvc::absolute(phi, U));
#include "alphaControls.H"
#include "compressibleAlphaEqnSubCycle.H"

View File

@ -8,6 +8,4 @@ CorrectPhi
pimple
);
//***HGW phi.oldTime() = phi;
#include "continuityErrs.H"

View File

@ -6,7 +6,7 @@
(
"phiHbyA",
fvc::flux(HbyA)
+ MRF.zeroFilter(fvc::interpolate(rho*rAU)*fvc::ddtCorr(U, phi))
+ MRF.zeroFilter(fvc::interpolate(rho*rAU)*fvc::ddtCorr(U, phi, Uf))
);
MRF.makeRelative(phiHbyA);
@ -23,6 +23,9 @@
// Update the pressure BCs to ensure flux consistency
constrainPressure(p_rgh, U, phiHbyA, rAUf, MRF);
// Make the fluxes relative to the mesh motion
fvc::makeRelative(phiHbyA, U);
tmp<fvScalarMatrix> p_rghEqnComp1;
tmp<fvScalarMatrix> p_rghEqnComp2;
@ -34,8 +37,7 @@
surfaceScalarField phid2("phid2", fvc::interpolate(psi2)*phi);
p_rghEqnComp1 =
pos(alpha1)
*(
(
(
fvc::ddt(alpha1, rho1) + fvc::div(alphaPhi1*rho1f)
- (fvOptions(alpha1, mixture.thermo1().rho())&rho1)
@ -48,11 +50,9 @@
+ fvm::div(phid1, p_rgh) - fvm::Sp(fvc::div(phid1), p_rgh)
)
);
p_rghEqnComp1.ref().relax();
p_rghEqnComp2 =
pos(alpha2)
*(
(
(
fvc::ddt(alpha2, rho2) + fvc::div(alphaPhi2*rho2f)
- (fvOptions(alpha2, mixture.thermo2().rho())&rho2)
@ -65,7 +65,6 @@
+ fvm::div(phid2, p_rgh) - fvm::Sp(fvc::div(phid2), p_rgh)
)
);
p_rghEqnComp2.ref().relax();
}
else
{
@ -94,6 +93,21 @@
);
}
if (mesh.moving())
{
p_rghEqnComp1.ref() += fvc::div(mesh.phi())*alpha1;
p_rghEqnComp2.ref() += fvc::div(mesh.phi())*alpha2;
}
p_rghEqnComp1.ref() *= pos(alpha1);
p_rghEqnComp2.ref() *= pos(alpha2);
if (pimple.transonic())
{
p_rghEqnComp1.ref().relax();
p_rghEqnComp2.ref().relax();
}
// Cache p_rgh prior to solve for density update
volScalarField p_rgh_0(p_rgh);
@ -131,6 +145,9 @@
}
}
// Correct Uf if the mesh is moving
fvc::correctUf(Uf, U, fvc::absolute(phi, U));
// Update densities from change in p_rgh
mixture.thermo1().correctRho(psi1*(p_rgh - p_rgh_0));
mixture.thermo2().correctRho(psi2*(p_rgh - p_rgh_0));

View File

@ -18,6 +18,6 @@ EXE_LIBS = \
-lincompressibleTurbulenceModels \
-lfiniteVolume \
-ldynamicFvMesh \
-ltopoChangerFvMesh \
-lfvOptions \
-lmeshTools \
-lsampling
-lmeshTools

View File

@ -3,6 +3,5 @@ cd ${0%/*} || exit 1 # Run from this directory
wclean libso phaseChangeTwoPhaseMixtures
wclean
wclean interPhaseChangeDyMFoam
#------------------------------------------------------------------------------

View File

@ -6,6 +6,5 @@ cd ${0%/*} || exit 1 # Run from this directory
wmake $targetType phaseChangeTwoPhaseMixtures
wmake $targetType
wmake $targetType interPhaseChangeDyMFoam
#------------------------------------------------------------------------------

View File

@ -1,14 +1,17 @@
EXE_INC = \
-I$(LIB_SRC)/transportModels/twoPhaseMixture/lnInclude \
-I. \
-I../VoF \
-I$(LIB_SRC)/transportModels \
-I$(LIB_SRC)/transportModels/incompressible/lnInclude \
-I$(LIB_SRC)/transportModels/twoPhaseMixture/lnInclude \
-I$(LIB_SRC)/transportModels/interfaceProperties/lnInclude \
-I$(LIB_SRC)/transportModels/immiscibleIncompressibleTwoPhaseMixture/lnInclude \
-IphaseChangeTwoPhaseMixtures/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \
-IphaseChangeTwoPhaseMixtures/lnInclude \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/sampling/lnInclude
-I$(LIB_SRC)/dynamicFvMesh/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude
EXE_LIBS = \
-lphaseChangeTwoPhaseMixtures \
@ -19,6 +22,7 @@ EXE_LIBS = \
-lturbulenceModels \
-lincompressibleTurbulenceModels \
-lfiniteVolume \
-lmeshTools \
-ldynamicFvMesh \
-ltopoChangerFvMesh \
-lfvOptions \
-lsampling
-lmeshTools

View File

@ -1,13 +1,16 @@
fvVectorMatrix UEqn
(
fvm::ddt(rho, U)
+ fvm::div(rhoPhi, U)
fvm::ddt(rho, U) + fvm::div(rhoPhi, U)
- fvm::Sp(fvc::ddt(rho) + fvc::div(rhoPhi), U)
+ turbulence->divDevRhoReff(rho, U)
==
fvOptions(rho, U)
);
UEqn.relax();
fvOptions.constrain(UEqn);
if (pimple.momentumPredictor())
{
solve
@ -17,10 +20,12 @@
fvc::reconstruct
(
(
interface.surfaceTensionForce()
mixture.surfaceTensionForce()
- ghf*fvc::snGrad(rho)
- fvc::snGrad(p_rgh)
) * mesh.magSf()
)
);
fvOptions.correct(U);
}

View File

@ -2,10 +2,9 @@
word alphaScheme("div(phi,alpha)");
word alpharScheme("div(phirb,alpha)");
surfaceScalarField phir("phir", phic*interface.nHatf());
surfaceScalarField phir("phir", phic*mixture.nHatf());
Pair<tmp<volScalarField>> vDotAlphal =
mixture->vDotAlphal();
Pair<tmp<volScalarField>> vDotAlphal = mixture.vDotAlphal();
const volScalarField& vDotcAlphal = vDotAlphal[0]();
const volScalarField& vDotvAlphal = vDotAlphal[1]();
const volScalarField vDotvmcAlphal(vDotvAlphal - vDotcAlphal);

View File

@ -1,16 +1,14 @@
{
// Standard face-flux compression coefficient
surfaceScalarField phic(interface.cAlpha()*mag(phi/mesh.magSf()));
surfaceScalarField phic(mixture.cAlpha()*mag(phi/mesh.magSf()));
// Add the optional isotropic compression contribution
if (icAlpha > 0)
{
phic *= (1.0 - icAlpha);
phic += (interface.cAlpha()*icAlpha)*fvc::interpolate(mag(U));
phic += (mixture.cAlpha()*icAlpha)*fvc::interpolate(mag(U));
}
volScalarField divU(fvc::div(phi));
if (nAlphaSubCycles > 1)
{
dimensionedScalar totalDeltaT = runTime.deltaT();

View File

@ -3,7 +3,7 @@ CorrectPhi
U,
phi,
p_rgh,
surfaceScalarField("rAUf", fvc::interpolate(rAU)),
surfaceScalarField("rAUf", fvc::interpolate(rAU())),
divU,
pimple
);

View File

@ -30,14 +30,18 @@ volVectorField U
Info<< "Creating phaseChangeTwoPhaseMixture\n" << endl;
autoPtr<phaseChangeTwoPhaseMixture> mixture =
phaseChangeTwoPhaseMixture::New(U, phi);
autoPtr<phaseChangeTwoPhaseMixture> mixturePtr
(
phaseChangeTwoPhaseMixture::New(U, phi)
);
volScalarField& alpha1(mixture->alpha1());
volScalarField& alpha2(mixture->alpha2());
phaseChangeTwoPhaseMixture& mixture = mixturePtr();
const dimensionedScalar& rho1 = mixture->rho1();
const dimensionedScalar& rho2 = mixture->rho2();
volScalarField& alpha1(mixture.alpha1());
volScalarField& alpha2(mixture.alpha2());
const dimensionedScalar& rho1 = mixture.rho1();
const dimensionedScalar& rho2 = mixture.rho2();
// Need to store rho for ddt(rho, U)
@ -55,13 +59,10 @@ volScalarField rho
rho.oldTime();
// Construct interface from alpha1 distribution
interfaceProperties interface(alpha1, U, mixture());
// Construct incompressible turbulence model
autoPtr<incompressible::turbulenceModel> turbulence
(
incompressible::turbulenceModel::New(U, phi, mixture())
incompressible::turbulenceModel::New(U, phi, mixture)
);

View File

@ -0,0 +1,36 @@
tmp<volScalarField> rAU;
if (correctPhi)
{
rAU = new volScalarField
(
IOobject
(
"rAU",
runTime.timeName(),
mesh,
IOobject::READ_IF_PRESENT,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("rAU", dimTime/dimDensity, 1)
);
volScalarField divU("divU0", fvc::div(fvc::absolute(phi, U)));
#include "correctPhi.H"
}
else
{
CorrectPhi
(
U,
phi,
p_rgh,
dimensionedScalar("rAUf", dimTime/rho.dimensions(), 1),
geometricZeroField(),
pimple
);
#include "continuityErrs.H"
}

View File

@ -1,3 +0,0 @@
interPhaseChangeDyMFoam.C
EXE = $(FOAM_APPBIN)/interPhaseChangeDyMFoam

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