All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
multiphaseEulerFoam.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2021 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 Application
25  multiphaseEulerFoam
26 
27 Description
28  Solver for a system of any number of compressible fluid phases with a
29  common pressure, but otherwise separate properties. The type of phase model
30  is run time selectable and can optionally represent multiple species and
31  in-phase reactions. The phase system is also run time selectable and can
32  optionally represent different types of momentum, heat and mass transfer.
33 
34 \*---------------------------------------------------------------------------*/
35 
36 #include "fvCFD.H"
37 #include "dynamicFvMesh.H"
38 #include "phaseSystem.H"
40 #include "pimpleControl.H"
41 #include "pressureReference.H"
42 #include "localEulerDdtScheme.H"
43 #include "fvcSmooth.H"
44 
45 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
46 
47 int main(int argc, char *argv[])
48 {
49  #include "postProcess.H"
50 
51  #include "setRootCaseLists.H"
52  #include "createTime.H"
53  #include "createDynamicFvMesh.H"
54  #include "createDyMControls.H"
55  #include "createFields.H"
56  #include "createFieldRefs.H"
57 
58  if (!LTS)
59  {
60  #include "CourantNo.H"
61  #include "setInitialDeltaT.H"
62  }
63 
64  Switch faceMomentum
65  (
66  pimple.dict().lookupOrDefault<Switch>("faceMomentum", false)
67  );
68  Switch partialElimination
69  (
70  pimple.dict().lookupOrDefault<Switch>("partialElimination", false)
71  );
72 
73  #include "createRDeltaTf.H"
74 
75  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
76 
77  Info<< "\nStarting time loop\n" << endl;
78 
79  while (pimple.run(runTime))
80  {
81  #include "readDyMControls.H"
82 
83  int nEnergyCorrectors
84  (
85  pimple.dict().lookupOrDefault<int>("nEnergyCorrectors", 1)
86  );
87 
88  if (LTS)
89  {
90  #include "setRDeltaT.H"
91  if (faceMomentum)
92  {
93  #include "setRDeltaTf.H"
94  }
95  }
96  else
97  {
98  #include "CourantNo.H"
99  #include "setDeltaT.H"
100  }
101 
102  runTime++;
103  Info<< "Time = " << runTime.timeName() << nl << endl;
104 
105  // --- Pressure-velocity PIMPLE corrector loop
106  while (pimple.loop())
107  {
108  if (!pimple.flow())
109  {
110  if (pimple.models())
111  {
112  fvModels.correct();
113  }
114 
115  if (pimple.thermophysics())
116  {
117  fluid.solve(rAUs, rAUfs);
118  fluid.correct();
119  fluid.correctContinuityError();
120 
121  #include "YEqns.H"
122  #include "EEqns.H"
123  #include "pEqnComps.H"
124 
126  {
127  phases[phasei].divU(-pEqnComps[phasei] & p_rgh);
128  }
129  }
130  }
131  else
132  {
133  if (pimple.firstPimpleIter() || moveMeshOuterCorrectors)
134  {
135  // Store divU from the previous mesh so that it can be
136  // mapped and used in correctPhi to ensure the corrected phi
137  // has the same divergence
138  tmp<volScalarField> divU;
139 
140  if
141  (
142  correctPhi
143  )
144  {
145  // Construct and register divU for mapping
146  divU = new volScalarField
147  (
148  "divU0",
149  fvc::div
150  (
151  fvc::absolute(phi, fluid.movingPhases()[0].U())
152  )
153  );
154  }
155 
157 
158  mesh.update();
159 
160  if (mesh.changing())
161  {
162  gh = (g & mesh.C()) - ghRef;
163  ghf = (g & mesh.Cf()) - ghRef;
164 
165  fluid.meshUpdate();
166 
167  if (correctPhi)
168  {
169  fluid.correctPhi
170  (
171  p_rgh,
172  divU,
174  pimple
175  );
176  }
177 
178  if (checkMeshCourantNo)
179  {
180  #include "meshCourantNo.H"
181  }
182  }
183  }
184 
185  if (pimple.models())
186  {
187  fvModels.correct();
188  }
189 
190  fluid.solve(rAUs, rAUfs);
191  fluid.correct();
192  fluid.correctContinuityError();
193 
194  if (pimple.thermophysics())
195  {
196  #include "YEqns.H"
197  }
198 
199  if (faceMomentum)
200  {
201  #include "pUf/UEqns.H"
202 
203  if (pimple.thermophysics())
204  {
205  #include "EEqns.H"
206  }
207 
208  #include "pUf/pEqn.H"
209  }
210  else
211  {
212  #include "pU/UEqns.H"
213 
214  if (pimple.thermophysics())
215  {
216  #include "EEqns.H"
217  }
218 
219  #include "pU/pEqn.H"
220  }
221 
222  fluid.correctKinematics();
223 
224  if (pimple.turbCorr())
225  {
226  fluid.correctTurbulence();
227  }
228  }
229  }
230 
231  runTime.write();
232 
233  Info<< "ExecutionTime = "
234  << runTime.elapsedCpuTime()
235  << " s\n\n" << endl;
236  }
237 
238  Info<< "End\n" << endl;
239 
240  return 0;
241 }
242 
243 
244 // ************************************************************************* //
pressureReference & pressureReference
volScalarField & p_rgh
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
pimpleNoLoopControl & pimple
label phasei
Definition: pEqn.H:27
Info<< "Reading field p_rgh\"<< endl;volScalarField p_rgh(IOobject("p_rgh", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);if(fluid.found("pMin")){ IOWarningInFunction(fluid)<< "Pressure limits, pMin and pMax, are now read from "<< pimple.dict().name()<< endl;}pressureReference pressureReference(p, p_rgh, pimple.dict(), fluid.incompressible());if(fluid.incompressible()){ p=p_rgh+fluid.rho() *gh;}if(p_rgh.needReference() &&fluid.incompressible()){ p+=dimensionedScalar("p", p.dimensions(), pressureReference.refValue() - getRefCellValue(p, pressureReference.refCell()));}p_rgh=p - fluid.rho() *gh;mesh.setFluxRequired(p_rgh.name());PtrList< volScalarField > rAUs
Definition: createFields.H:67
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
virtual void correct()
Correct the fvModels.
Definition: fvModels.C:316
engineTime & runTime
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
virtual void preUpdateMesh()
Prepare for mesh update.
Definition: fvModels.C:240
correctPhi
checkMeshCourantNo
PtrList< fvScalarMatrix > pEqnComps(phases.size())
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
dynamicFvMesh & mesh
phaseSystem & fluid
Definition: createFields.H:11
phi
Definition: correctPhi.H:3
const surfaceScalarField & ghf
static const char nl
Definition: Ostream.H:260
bool LTS
Definition: createRDeltaT.H:1
PtrList< surfaceScalarField > rAUfs
Definition: createFields.H:68
moveMeshOuterCorrectors
tmp< volScalarField > divU
Definition: alphaSuSp.H:9
Foam::fvModels & fvModels
Calculates and outputs the mean and maximum Courant Numbers.
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
Definition: fvcMeshPhi.C:188
phaseSystem::phaseModelList & phases
Definition: createFields.H:12
const volScalarField & gh
messageStream Info
const dimensionedVector & g
Execute application functionObjects to post-process existing results.
Provides functions smooth spread and sweep which use the FaceCellWave algorithm to smooth and redistr...