All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
multiphaseEuler.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) 2022-2024 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 \*---------------------------------------------------------------------------*/
25 
26 #include "multiphaseEuler.H"
27 #include "localEulerDdtScheme.H"
28 #include "surfaceFields.H"
29 #include "fvcDiv.H"
30 #include "fvcSurfaceIntegrate.H"
31 #include "fvcMeshPhi.H"
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38 namespace solvers
39 {
42 }
43 }
44 
45 
46 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
47 
49 {
51 
53  pimple.dict().lookupOrDefault<bool>("momentumPredictor", false);
54 
55  faceMomentum =
56  pimple.dict().lookupOrDefault<Switch>("faceMomentum", false);
57 
59  pimple.dict().lookupOrDefault<Switch>("dragCorrection", false);
60 
62  pimple.dict().lookupOrDefault<int>("nEnergyCorrectors", 1);
63 
64  return true;
65 }
66 
67 
68 void Foam::solvers::multiphaseEuler::correctCoNum()
69 {
70  scalarField sumPhi
71  (
72  fvc::surfaceSum(mag(phi))().primitiveField()
73  );
74 
75  forAll(movingPhases, movingPhasei)
76  {
77  sumPhi = max
78  (
79  sumPhi,
80  fvc::surfaceSum(mag(movingPhases[movingPhasei].phi()))()
81  .primitiveField()
82  );
83  }
84 
85  CoNum_ = 0.5*gMax(sumPhi/mesh.V().primitiveField())*runTime.deltaTValue();
86 
87  const scalar meanCoNum =
88  0.5
89  *(gSum(sumPhi)/gSum(mesh.V().primitiveField()))
90  *runTime.deltaTValue();
91 
92  Info<< "Courant Number mean: " << meanCoNum
93  << " max: " << CoNum << endl;
94 }
95 
96 
97 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
98 
100 :
101  fluidSolver(mesh),
102 
103  predictMomentum
104  (
105  pimple.dict().lookupOrDefault<Switch>("momentumPredictor", false)
106  ),
107 
108  faceMomentum
109  (
110  pimple.dict().lookupOrDefault<Switch>("faceMomentum", false)
111  ),
112 
113  dragCorrection
114  (
115  pimple.dict().lookupOrDefault<Switch>("dragCorrection", false)
116  ),
117 
118  nEnergyCorrectors
119  (
120  pimple.dict().lookupOrDefault<int>("nEnergyCorrectors", 1)
121  ),
122 
123  trDeltaT
124  (
125  LTS
126  ? new volScalarField
127  (
128  IOobject
129  (
130  fv::localEulerDdt::rDeltaTName,
131  runTime.name(),
132  mesh,
133  IOobject::READ_IF_PRESENT,
134  IOobject::AUTO_WRITE
135  ),
136  mesh,
138  extrapolatedCalculatedFvPatchScalarField::typeName
139  )
140  : nullptr
141  ),
142 
143  trDeltaTf
144  (
145  LTS && faceMomentum
146  ? new surfaceScalarField
147  (
148  IOobject
149  (
150  fv::localEulerDdt::rDeltaTfName,
151  runTime.name(),
152  mesh,
153  IOobject::READ_IF_PRESENT,
154  IOobject::AUTO_WRITE
155  ),
156  mesh,
158  )
159  : nullptr
160  ),
161 
162  buoyancy(mesh),
163 
164  fluidPtr_(phaseSystem::New(mesh)),
165 
166  fluid_(fluidPtr_()),
167 
168  phases_(fluid_.phases()),
169 
170  movingPhases_(fluid_.movingPhases()),
171 
172  phi_(fluid_.phi()),
173 
174  p_(movingPhases_[0].fluidThermo().p()),
175 
176  p_rgh(buoyancy.p_rgh),
177 
179  (
180  p_,
181  p_rgh,
182  pimple.dict(),
183  fluid_.incompressible()
184  ),
185 
186  MRF(fluid_.MRF()),
187 
188  fluid(fluid_),
189  phases(phases_),
190  movingPhases(movingPhases_),
191  p(p_),
192  phi(phi_)
193 {
194  // Read the controls
195  read();
196 
198 
199  if (transient())
200  {
201  correctCoNum();
202  }
203 }
204 
205 
206 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
207 
209 {}
210 
211 
212 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
213 
215 {
216  if (transient())
217  {
218  correctCoNum();
219  }
220  else if (LTS)
221  {
222  setRDeltaT();
223  }
224 
225  // Store divU from the previous mesh so that it can be
226  // mapped and used in correctPhi to ensure the corrected phi
227  // has the same divergence
228  if (correctPhi || mesh.topoChanging())
229  {
230  // Construct and register divU for mapping
231  divU = new volScalarField
232  (
233  "divU0",
234  fvc::div(fvc::absolute(phi, movingPhases[0].U()))
235  );
236  }
237 
239 
240  // Update the mesh for topology change, mesh to mesh mapping
241  mesh_.update();
242 }
243 
244 
246 {
247  if (pimple.thermophysics() || pimple.flow())
248  {
249  fluid_.solve(rAs);
250  fluid_.correct();
251  fluid_.correctContinuityError();
252  }
253 
254  if (pimple.flow() && pimple.predictTransport())
255  {
256  fluid_.predictMomentumTransport();
257  }
258 }
259 
260 
262 {
263  if (pimple.flow() && pimple.correctTransport())
264  {
265  fluid_.correctMomentumTransport();
266  fluid_.correctThermophysicalTransport();
267  }
268 }
269 
270 
272 {
273  divU.clear();
274 }
275 
276 
277 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Macros for easy insertion into run-time selection tables.
Generic GeometricField class.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
const word & name() const
Return name.
Definition: IOobject.H:310
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
Definition: Switch.H:61
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T, if not found return the given default.
Base-class for fluid thermodynamic properties.
Definition: fluidThermo.H:57
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
const fvSchemes & schemes() const
Return the fvSchemes.
Definition: fvMesh.C:1751
virtual void preUpdateMesh()
Prepare for mesh update.
Definition: fvModels.C:260
void setFluxRequired(const word &name) const
Definition: fvSchemes.C:503
Class to represent a system of phases and model interfacial transfers between them.
Definition: phaseSystem.H:73
Provides controls for the pressure reference in closed-volume simulations.
virtual const dictionary & dict() const
Return the solution dictionary.
Abstract base class for run-time selectable region solvers.
Definition: solver.H:56
pimpleNoLoopControl pimple
PIMPLE inner-loop controls.
Definition: solver.H:107
const fvMesh & mesh
Region mesh.
Definition: solver.H:101
Buoyancy related data for the Foam::solvers::isothermalFluid solver module when solving buoyant cases...
Definition: buoyancy.H:70
Base solver module for fluid solvers.
Definition: fluidSolver.H:62
virtual bool read()
Read controls.
Definition: fluidSolver.C:51
Solver module for steady or transient turbulent flow of compressible fluids with heat-transfer for HV...
Definition: fluid.H:70
Solver module for a system of any number of compressible fluid phases with a common pressure,...
volScalarField & p_rgh
Reference to the buoyant pressure for buoyant cases.
virtual void prePredictor()
Called at the start of the PIMPLE loop.
virtual void postSolve()
Called after the PIMPLE loop at the end of the time-step.
int nEnergyCorrectors
Number of energy correctors.
virtual void postCorrector()
Correct the momentum and thermophysical transport modelling.
Switch predictMomentum
Momentum equation predictor switch.
Switch faceMomentum
Cell/face momentum equation switch.
Switch dragCorrection
Cell/face drag correction for cell momentum corrector.
virtual void preSolve()
Called at the start of the time-step, before the PIMPLE loop.
virtual ~multiphaseEuler()
Destructor.
multiphaseEuler(fvMesh &mesh)
Construct from region mesh.
virtual bool read()
Read controls.
scalar CoNum
scalar meanCoNum
Foam::fvModels & fvModels(Foam::fvModels::New(mesh))
IOMRFZoneList MRF(mesh)
pimpleControl pimple(mesh)
Calculate the divergence of the given field.
Calculate the mesh motion flux and convert fluxes from absolute to relative and back.
Surface integrate surfaceField creating a volField. Surface sum a surfaceField creating a volField.
U
Definition: pEqn.H:72
tmp< VolField< Type > > surfaceSum(const SurfaceField< Type > &ssf)
tmp< VolField< Type > > div(const SurfaceField< Type > &ssf)
Definition: fvcDiv.C:47
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
Definition: fvcMeshPhi.C:202
addToRunTimeSelectionTable(solver, compressibleMultiphaseVoF, fvMesh)
defineTypeNameAndDebug(compressibleMultiphaseVoF, 0)
Namespace for OpenFOAM.
Type gSum(const FieldField< Field, Type > &f)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
const dimensionSet dimless
messageStream Info
const dimensionSet dimTime
VolField< scalar > volScalarField
Definition: volFieldsFwd.H:64
dimensioned< scalar > mag(const dimensioned< Type > &)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
Type gMax(const FieldField< Field, Type > &f)
correctPhi
labelList fv(nPoints)
dictionary dict
volScalarField & p
Foam::surfaceFields.