compressibleVoF.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-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 
26 #include "compressibleVoF.H"
27 #include "localEulerDdtScheme.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 :
49  (
50  mesh,
52  ),
53 
54  mixture_
55  (
57  ),
58 
59  p(mixture_.p()),
60 
61  dgdt
62  (
63  IOobject
64  (
65  "dgdt",
66  runTime.name(),
67  mesh,
68  IOobject::READ_IF_PRESENT,
69  IOobject::AUTO_WRITE
70  ),
71  alpha1*fvc::div(phi)
72  ),
73 
74  pressureReference_
75  (
76  p,
77  p_rgh,
78  pimple.dict(),
79  false
80  ),
81 
82  alphaRhoPhi1
83  (
84  IOobject::groupName("alphaRhoPhi", alpha1.group()),
85  fvc::interpolate(mixture_.thermo1().rho())*alphaPhi1
86  ),
87 
88  alphaRhoPhi2
89  (
90  IOobject::groupName("alphaRhoPhi", alpha2.group()),
91  fvc::interpolate(mixture_.thermo2().rho())*(phi - alphaPhi1)
92  ),
93 
94  K("K", 0.5*magSqr(U)),
95 
96  momentumTransport
97  (
98  rho,
99  U,
100  phi,
101  rhoPhi,
102  alphaPhi1,
103  alphaRhoPhi1,
104  alphaRhoPhi2,
105  mixture_
106  ),
107 
108  thermophysicalTransport(momentumTransport),
109 
110  mixture(mixture_)
111 {
112  // Read the controls
113  readControls();
114 
115  if (correctPhi || mesh.topoChanging())
116  {
117  rAU = new volScalarField
118  (
119  IOobject
120  (
121  "rAU",
122  runTime.name(),
123  mesh,
126  ),
127  mesh,
129  );
130  }
131 }
132 
133 
134 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
135 
137 {}
138 
139 
140 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
141 
143 {
145 
146  const volScalarField& rho1 = mixture_.thermo1().rho();
147  const volScalarField& rho2 = mixture_.thermo2().rho();
148 
149  alphaRhoPhi1 = fvc::interpolate(rho1)*alphaPhi1;
150  alphaRhoPhi2 = fvc::interpolate(rho2)*(phi - alphaPhi1);
151 
152  rhoPhi = alphaRhoPhi1 + alphaRhoPhi2;
153 
154  contErr1 =
155  (
156  fvc::ddt(alpha1, rho1)()() + fvc::div(alphaRhoPhi1)()()
157  - (fvModels().source(alpha1, rho1)&rho1)()
158  );
159 
160  contErr2 =
161  (
162  fvc::ddt(alpha2, rho2)()() + fvc::div(alphaRhoPhi2)()()
163  - (fvModels().source(alpha2, rho2)&rho2)()
164  );
165 
166  if (pimple.predictTransport())
167  {
168  momentumTransport.predict();
169  }
170 
171  if (pimple.predictTransport())
172  {
173  thermophysicalTransport.predict();
174  }
175 }
176 
177 
179 {
180  if (pimple.correctTransport())
181  {
182  momentumTransport.correct();
183  thermophysicalTransport.correct();
184  }
185 }
186 
187 
188 // ************************************************************************* //
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
Class to represent a mixture of two rhoThermo-based phases.
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 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 for 2 compressible, non-isothermal immiscible fluids using a VOF (volume of fluid) ...
virtual void prePredictor()
Called at the start of the PIMPLE loop.
virtual void postCorrector()
Correct the momentum and thermophysical transport modelling.
compressibleVoF(fvMesh &mesh)
Construct from region mesh.
virtual ~compressibleVoF()
Destructor.
bool correctPhi
Switch to correct the flux after mesh change.
Definition: fluidSolver.H:95
void readControls()
Read controls.
Definition: fluidSolver.C:45
virtual void prePredictor()
Called at the start of the PIMPLE loop.
Solver module base-class for for 2 immiscible fluids using a VOF (volume of fluid) phase-fraction bas...
Class to represent a VoF mixture.
Foam::fvModels & fvModels(Foam::fvModels::New(mesh))
pimpleControl pimple(mesh)
Calculate the first temporal derivative.
Calculate the divergence of the given field.
Info<< "Creating thermophysical transport model\n"<< endl;turbulenceThermophysicalTransportModels::unityLewisEddyDiffusivity< RASThermophysicalTransportModel< ThermophysicalTransportModel< compressibleMomentumTransportModel, fluidThermo > >> thermophysicalTransport(turbulence(), thermo, true)
K
Definition: pEqn.H:75
U
Definition: pEqn.H:72
const char *const group
Group name for atomic constants.
static tmp< SurfaceField< Type > > interpolate(const VolField< Type > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
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
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
const dimensionSet dimTime
const dimensionSet dimDensity
VolField< scalar > volScalarField
Definition: volFieldsFwd.H:61
dimensioned< scalar > magSqr(const dimensioned< Type > &)
dictionary dict
volScalarField & p