incompressibleMultiphaseVoF.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 "fvCorrectPhi.H"
28 #include "geometricZeroField.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace solvers
36 {
39 }
40 }
41 
42 
43 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
44 
46 (
47  fvMesh& mesh
48 )
49 :
51  (
52  mesh,
54  (
56  )
57  ),
58 
59  mixture
60  (
62  (
64  )
65  ),
66 
67  phases(mixture.phases()),
68 
69  p
70  (
71  IOobject
72  (
73  "p",
74  runTime.name(),
75  mesh,
76  IOobject::NO_READ,
77  IOobject::AUTO_WRITE
78  ),
79  p_rgh + rho*buoyancy.gh
80  ),
81 
82  pressureReference_
83  (
84  p,
85  p_rgh,
86  pimple.dict(),
87  false
88  ),
89 
90  momentumTransport_
91  (
92  incompressible::momentumTransportModel::New
93  (
94  U,
95  phi,
96  mixture
97  )
98  ),
99 
100  momentumTransport(momentumTransport_())
101 {
102  if (correctPhi || mesh.topoChanging())
103  {
104  rAU = new volScalarField
105  (
106  IOobject
107  (
108  "rAU",
109  runTime.name(),
110  mesh,
113  ),
114  mesh,
116  );
117  }
118 
119  if (!runTime.restart() || !divergent())
120  {
121  correctUphiBCs(U_, phi_, true);
122 
124  (
125  phi_,
126  U,
127  p_rgh,
128  rAU,
131  pimple
132  );
133  }
134 }
135 
136 
137 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
138 
140 {}
141 
142 
143 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
144 
146 {
148 
149  if (pimple.predictTransport())
150  {
151  momentumTransport.predict();
152  }
153 }
154 
155 
157 {}
158 
159 
161 {
162  if (pimple.correctTransport())
163  {
164  momentumTransport.correct();
165  }
166 }
167 
168 
169 // ************************************************************************* //
Macros for easy insertion into run-time selection tables.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
bool restart() const
Return true if the run is a restart, i.e. startTime != beginTime.
Definition: Time.H:432
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
const word & name() const
Return const reference to name.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
bool topoChanging() const
Does the mesh topology change?
Definition: fvMesh.C:598
Incompressible multiphase mixture for interface-capturing simulations.
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:56
pimpleNoLoopControl pimple
PIMPLE inner-loop controls.
Definition: solver.H:107
const Time & runTime
Time.
Definition: solver.H:104
const fvMesh & mesh
Region mesh.
Definition: solver.H:101
volScalarField & p_rgh
Reference to the buoyant pressure for buoyant cases.
Definition: VoFSolver.H:107
volVectorField U_
Velocity field.
Definition: VoFSolver.H:94
autoPtr< volScalarField > rAU
Inverse momentum equation diagonal.
Definition: VoFSolver.H:129
const volVectorField & U
Reference to the velocity field.
Definition: VoFSolver.H:209
surfaceScalarField phi_
Volumetric flux field.
Definition: VoFSolver.H:97
Buoyancy related data for the Foam::solvers::isothermalFluid solver module when solving buoyant cases...
Definition: buoyancy.H:70
bool correctPhi
Switch to correct the flux after mesh change.
Definition: fluidSolver.H:101
Solver module for the solution of multiple incompressible, isothermal immiscible fluids using a VOF (...
virtual void thermophysicalPredictor()
Construct and solve the energy equation,.
virtual const Foam::pressureReference & pressureReference() const
Return the pressure reference.
virtual void prePredictor()
Called at the start of the PIMPLE loop.
virtual bool divergent() const
The flow is not incompressible and hence not divergent.
virtual void postCorrector()
Correct the momentum and thermophysical transport modelling.
incompressibleMultiphaseVoF(fvMesh &mesh)
Construct from region mesh.
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.
pimpleControl pimple(mesh)
Flux correction functions to ensure continuity.
U
Definition: pEqn.H:72
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
void correctPhi(surfaceScalarField &phi, const volVectorField &U, const volScalarField &p, const autoPtr< volScalarField > &rAU, const autoPtr< volScalarField > &divU, const pressureReference &pressureReference, nonOrthogonalSolutionControl &pcorrControl)
Definition: fvCorrectPhi.C:32
addToRunTimeSelectionTable(solver, compressibleMultiphaseVoF, fvMesh)
defineTypeNameAndDebug(compressibleMultiphaseVoF, 0)
Namespace for OpenFOAM.
To & refCast(From &r)
Reference type cast template function.
Definition: typeInfo.H:129
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:64
void correctUphiBCs(volVectorField &U, surfaceScalarField &phi, const bool evaluateUBCs)
If the mesh is moving correct the velocity BCs on the moving walls to.
dictionary dict
volScalarField & p