compressibleMultiphaseVoF.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) 2023 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 
27 #include "geometricZeroField.H"
28 #include "fvcDdt.H"
29 #include "fvcDiv.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 namespace solvers
37 {
40 }
41 }
42 
43 
44 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
45 
47 (
48  fvMesh& mesh
49 )
50 :
52  (
53  mesh,
55  (
57  )
58  ),
59 
60  mixture
61  (
63  (
65  )
66  ),
67 
68  phases(mixture.phases()),
69 
70  p(mixture.p()),
71 
72  pressureReference_
73  (
74  p,
75  p_rgh,
76  pimple.dict(),
77  false
78  ),
79 
80  K("K", 0.5*magSqr(U)),
81 
82  momentumTransport_
83  (
84  compressible::momentumTransportModel::New
85  (
86  rho,
87  U,
88  rhoPhi,
89  mixture
90  )
91  ),
92 
93  momentumTransport(momentumTransport_())
94 {
95  // Read the controls
96  readControls();
97 
98  if (correctPhi || mesh.topoChanging())
99  {
100  rAU = new volScalarField
101  (
102  IOobject
103  (
104  "rAU",
105  runTime.name(),
106  mesh,
109  ),
110  mesh,
112  );
113  }
114 }
115 
116 
117 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
118 
120 {}
121 
122 
123 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
124 
126 {
128 
129  contErr = fvc::ddt(rho)()() + fvc::div(rhoPhi)()();
130 
131  forAll(mixture.phases(), phasei)
132  {
133  const volScalarField& rho = phases[phasei].thermo().rho();
134  contErr.ref() -= fvModels().source(phases[phasei], rho)&rho;
135  }
136 
137  if (pimple.predictTransport())
138  {
139  momentumTransport.predict();
140  }
141 }
142 
143 
145 {
146  if (pimple.correctTransport())
147  {
148  momentumTransport.correct();
149  }
150 }
151 
152 
153 // ************************************************************************* //
#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
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
Compressible multiphase mixture for interface-capturing simulations.
const word & name() const
Return const reference to name.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:101
bool topoChanging() const
Does the mesh topology change?
Definition: fvMesh.C:647
Abstract base class for turbulence models (RAS, LES and laminar).
Multiphase VoF mixture with support for interface properties.
Abstract base class for run-time selectable region solvers.
Definition: solver.H:55
const Time & runTime
Time.
Definition: solver.H:97
const fvMesh & mesh
Region mesh.
Definition: solver.H:94
autoPtr< volScalarField > rAU
Inverse momentum equation diagonal.
Definition: VoFSolver.H:129
Solver module for the solution of multiple compressible, isothermal immiscible fluids using a VOF (vo...
virtual void prePredictor()
Called at the start of the PIMPLE loop.
virtual void postCorrector()
Correct the momentum and thermophysical transport modelling.
compressibleMultiphaseVoF(fvMesh &mesh)
Construct from region mesh.
bool correctPhi
Switch to correct the flux after mesh change.
Definition: fluidSolver.H:95
void readControls()
Read controls.
Definition: fluidSolver.C:45
Base solver module for the solution of multiple immiscible fluids using a VOF (volume of fluid) phase...
virtual void prePredictor()
Called at the start of the PIMPLE loop.
Foam::fvModels & fvModels(Foam::fvModels::New(mesh))
pimpleControl pimple(mesh)
Calculate the first temporal derivative.
Calculate the divergence of the given field.
K
Definition: pEqn.H:75
U
Definition: pEqn.H:72
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
tmp< VolField< Type > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
Definition: fvcDdt.C:45
tmp< VolField< Type > > div(const SurfaceField< Type > &ssf)
Definition: fvcDiv.C:47
addToRunTimeSelectionTable(solver, compressibleMultiphaseVoF, fvMesh)
defineTypeNameAndDebug(compressibleMultiphaseVoF, 0)
Namespace for OpenFOAM.
To & refCast(From &r)
Reference type cast template function.
Definition: typeInfo.H:111
const dimensionSet dimTime
const dimensionSet dimDensity
VolField< scalar > volScalarField
Definition: volFieldsFwd.H:61
dimensioned< scalar > magSqr(const dimensioned< Type > &)
dictionary dict
volScalarField & p