overInterDyMFoam.C
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2014 OpenFOAM Foundation
9 Copyright (C) 2016-2017 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27Application
28 overInterDyMFoam
29
30Group
31 grpMultiphaseSolvers grpMovingMeshSolvers
32
33Description
34 Solver for two incompressible, isothermal immiscible fluids using a VOF
35 (volume of fluid) phase-fraction based interface capturing approach,
36 with optional mesh motion and mesh topology changes including adaptive
37 re-meshing.
38
39\*---------------------------------------------------------------------------*/
40
41#include "fvCFD.H"
42#include "dynamicFvMesh.H"
43#include "CMULES.H"
44#include "EulerDdtScheme.H"
45#include "localEulerDdtScheme.H"
47#include "subCycle.H"
50#include "pimpleControl.H"
51#include "fvOptions.H"
52#include "CorrectPhi.H"
53#include "fvcSmooth.H"
55#include "localMin.H"
57#include "transform.H"
58#include "fvMeshSubset.H"
59#include "oversetAdjustPhi.H"
60
61// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
62
63int main(int argc, char *argv[])
64{
65 argList::addNote
66 (
67 "Solver for two incompressible, isothermal immiscible fluids using"
68 " VOF phase-fraction based interface capturing\n"
69 "With optional mesh motion and mesh topology changes including"
70 " adaptive re-meshing."
71 );
72
73 #include "postProcess.H"
74
75 #include "setRootCaseLists.H"
76 #include "createTime.H"
77 #include "createDynamicFvMesh.H"
78 #include "initContinuityErrs.H"
79 pimpleControl pimple(mesh);
80 #include "createTimeControls.H"
81 #include "createDyMControls.H"
82 #include "createFields.H"
83 #include "createAlphaFluxes.H"
84 #include "createFvOptions.H"
85
87 (
88 IOobject
89 (
90 "rAU",
91 runTime.timeName(),
92 mesh,
93 IOobject::READ_IF_PRESENT,
94 IOobject::AUTO_WRITE
95 ),
96 mesh,
97 dimensionedScalar("rAUf", dimTime/rho.dimensions(), 1.0)
98 );
99
100 if (correctPhi)
101 {
102 #include "correctPhi.H"
103 }
104 #include "createUf.H"
105
106 #include "setCellMask.H"
107 #include "setInterpolatedCells.H"
108
109 turbulence->validate();
110
111 if (!LTS)
112 {
113 #include "CourantNo.H"
114 #include "setInitialDeltaT.H"
115 }
116
117 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
118 Info<< "\nStarting time loop\n" << endl;
119
120 while (runTime.run())
121 {
122 #include "readControls.H"
123
124 if (LTS)
125 {
126 #include "setRDeltaT.H"
127 }
128 else
129 {
130 #include "CourantNo.H"
131 #include "alphaCourantNo.H"
132 #include "setDeltaT.H"
133 }
134
135 ++runTime;
136
137 Info<< "Time = " << runTime.timeName() << nl << endl;
138
139 // --- Pressure-velocity PIMPLE corrector loop
140 while (pimple.loop())
141 {
142 if (pimple.firstIter() || moveMeshOuterCorrectors)
143 {
144 scalar timeBeforeMeshUpdate = runTime.elapsedCpuTime();
145
146 mesh.update();
147
148 if (mesh.changing())
149 {
150 Info<< "Execution time for mesh.update() = "
151 << runTime.elapsedCpuTime() - timeBeforeMeshUpdate
152 << " s" << endl;
153
154 // Do not apply previous time-step mesh compression flux
155 // if the mesh topology changed
156 if (mesh.topoChanging())
157 {
158 talphaPhi1Corr0.clear();
159 }
160
161 gh = (g & mesh.C()) - ghRef;
162 ghf = (g & mesh.Cf()) - ghRef;
163
164 // Update cellMask field for blocking out hole cells
165 #include "setCellMask.H"
166 #include "setInterpolatedCells.H"
167
168 const surfaceScalarField faceMaskOld
169 (
170 localMin<scalar>(mesh).interpolate(cellMask.oldTime())
171 );
172
173 // Zero Uf on old faceMask (H-I)
174 Uf *= faceMaskOld;
175
176 const surfaceVectorField Uint(fvc::interpolate(U));
177 // Update Uf and phi on new C-I faces
178 Uf += (1-faceMaskOld)*Uint;
179
180 // Update Uf boundary
181 forAll(Uf.boundaryField(), patchI)
182 {
183 Uf.boundaryFieldRef()[patchI] =
184 Uint.boundaryField()[patchI];
185 }
186
187 phi = mesh.Sf() & Uf;
188
189 // Correct phi on individual regions
190 if (correctPhi)
191 {
192 #include "correctPhi.H"
193 }
194
195 mixture.correct();
196
197 // Zero phi on current H-I
199 (
200 localMin<scalar>(mesh).interpolate(cellMask)
201 );
202 phi *= faceMask;
203 U *= cellMask;
204
205 // Make the flux relative to the mesh motion
206 fvc::makeRelative(phi, U);
207
208 }
209
210 if (mesh.changing() && checkMeshCourantNo)
211 {
212 #include "meshCourantNo.H"
213 }
214 }
215
216
217 #include "alphaControls.H"
218 #include "alphaEqnSubCycle.H"
219
221 (
222 localMin<scalar>(mesh).interpolate(cellMask)
223 );
224 rhoPhi *= faceMask;
225
226 mixture.correct();
227
228 #include "UEqn.H"
229
230 // --- Pressure corrector loop
231 while (pimple.correct())
232 {
233 #include "pEqn.H"
234 }
235
236 if (pimple.turbCorr())
237 {
238 turbulence->correct();
239 }
240 }
241
242 runTime.write();
243
244 runTime.printExecutionTime(Info);
245 }
246
247 Info<< "End\n" << endl;
248
249 return 0;
250}
251
252
253// ************************************************************************* //
CMULES: Multidimensional universal limiter for explicit corrected implicit solution.
rhoPhi
Definition: rhoEqn.H:10
const uniformDimensionedVectorField & g
surfaceScalarField & phi
const surfaceScalarField & ghf
const volScalarField & gh
pimpleControl & pimple
U
Definition: pEqn.H:72
surfaceScalarField faceMask(localMin< scalar >(mesh).interpolate(cellMask))
mesh interpolate(rAU)
tmp< surfaceScalarField > talphaPhi1Corr0
dynamicFvMesh & mesh
engineTime & runTime
bool LTS
Definition: createRDeltaT.H:1
Read the control parameters used by setDeltaT.
autoPtr< surfaceVectorField > Uf
Creates and initialises the velocity velocity field Uf.
compressible::turbulenceModel & turbulence
Provides functions smooth spread and sweep which use the FaceCellWave algorithm to smooth and redistr...
tmp< volScalarField > rAU
Definition: initCorrectPhi.H:1
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
Definition: volFieldsFwd.H:82
messageStream Info
Information stream (stdout output on master, null elsewhere)
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
Adjust the balance of fluxes on the faces between interpolated and calculated to obey continuity.
Execute application functionObjects to post-process existing results.
checkMeshCourantNo
correctPhi
moveMeshOuterCorrectors
Sets blocked cells mask field.
Sets blocked cells mask field.
Info<< "Creating temperaturePhaseChangeTwoPhaseMixture\n"<< endl;autoPtr< temperaturePhaseChangeTwoPhaseMixture > mixture
Definition: createFields.H:39
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
3D tensor transformation operations.