69int main(
int argc,
char *argv[])
73 "Solver for two incompressible, isothermal, immiscible fluids with"
75 "using VOF (volume of fluid) phase-fraction based interface capturing,"
76 " with optional dynamic mesh motion (including overset)\n"
77 "and mesh topology changes including adaptive re-meshing."
88 #include "createDyMControls.H"
89 #include "initContinuityErrs.H"
90 #include "createFields.H"
99 IOobject::READ_IF_PRESENT,
107 #include "CourantNo.H"
108 #include "setInitialDeltaT.H"
117 Info<<
"\nStarting time loop\n" <<
endl;
121 #include "readControls.H"
128 #include "CourantNo.H"
129 #include "setDeltaT.H"
140 scalar timeBeforeMeshUpdate =
runTime.elapsedCpuTime();
146 Info<<
"Execution time for mesh.update() = "
147 <<
runTime.elapsedCpuTime() - timeBeforeMeshUpdate
158 localMin<scalar>(
mesh).interpolate(cellMask.oldTime());
171 Uf.boundaryFieldRef()[patchI] =
172 Uint.boundaryField()[patchI];
179 #include "correctPhi.H"
185 faceMask = localMin<scalar>(
mesh).interpolate(cellMask);
191 fvc::makeRelative(
phi,
U);
200 #include "alphaControls.H"
216 #include "alphaEqnSubCycle.H"
241 runTime.printExecutionTime(Info);
CMULES: Multidimensional universal limiter for explicit corrected implicit solution.
const uniformDimensionedVectorField & g
const surfaceScalarField & ghf
const volScalarField & gh
surfaceScalarField faceMask(localMin< scalar >(mesh).interpolate(cellMask))
Read the control parameters used by setDeltaT.
autoPtr< surfaceVectorField > Uf
Creates and initialises the velocity velocity field Uf.
compressible::turbulenceModel & turbulence
tmp< volScalarField > rAU
Calculates and outputs the mean and maximum Courant Numbers.
GeometricField< vector, fvsPatchField, surfaceMesh > surfaceVectorField
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
messageStream Info
Information stream (stdout output on master, null elsewhere)
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Ostream & endl(Ostream &os)
Add newline and flush stream.
constexpr char nl
The newline '\n' character (0x0a)
Adjust the balance of fluxes on the faces between interpolated and calculated to obey continuity.
Execute application functionObjects to post-process existing results.
Sets blocked cells mask field.
Sets blocked cells mask field.
interfaceProperties interface(alpha1, U, thermo->transportPropertiesDict())
Info<< "Creating temperaturePhaseChangeTwoPhaseMixture\n"<< endl;autoPtr< temperaturePhaseChangeTwoPhaseMixture > mixture
#define forAll(list, i)
Loop across all elements in list.