VoFClouds.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) 2021-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 "VoFClouds.H"
28 #include "fvmSup.H"
31 
32 // * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36  namespace fv
37  {
39 
41  (
42  fvModel,
43  VoFClouds,
45  );
46  }
47 }
48 
49 
50 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
51 
53 (
54  const word& sourceName,
55  const word& modelType,
56  const fvMesh& mesh,
57  const dictionary& dict
58 )
59 :
60  fvModel(sourceName, modelType, mesh, dict),
61  phaseName_(dict.lookup("phase")),
62  carrierPhaseName_(dict.lookup("carrierPhase")),
63  thermo_
64  (
65  mesh.lookupObject<fluidThermo>
66  (
67  IOobject::groupName(physicalProperties::typeName, phaseName_)
68  )
69  ),
70  carrierThermo_
71  (
72  mesh.lookupObject<fluidThermo>
73  (
74  IOobject::groupName(physicalProperties::typeName, carrierPhaseName_)
75  )
76  ),
77  clouds_
78  (
79  carrierThermo_.rho(),
80  mesh.lookupObject<volVectorField>("U"),
81  mesh.lookupObject<uniformDimensionedVectorField>("g"),
82  carrierThermo_
83  ),
84  curTimeIndex_(-1)
85 {}
86 
87 
88 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
89 
91 {
92  return wordList({thermo_.rho()().name(), thermo_.he().name(), "U"});
93 }
94 
95 
97 {
98  if (curTimeIndex_ == mesh().time().timeIndex())
99  {
100  return;
101  }
102 
103  clouds_.evolve();
104 
105  curTimeIndex_ = mesh().time().timeIndex();
106 }
107 
108 
110 (
111  const volScalarField& alpha,
112  const volScalarField& rho,
113  fvMatrix<scalar>& eqn
114 ) const
115 {
116  if (debug)
117  {
118  Info<< type() << ": applying source to " << eqn.psi().name() << endl;
119  }
120 
121  if (&rho == &thermo_.rho()())
122  {
123  eqn += clouds_.Srho();
124  }
125  else
126  {
128  << "Support for field " << rho.name() << " is not implemented"
129  << exit(FatalError);
130  }
131 }
132 
133 
135 (
136  const volScalarField& alpha,
137  const volScalarField& rho,
138  const volScalarField& he,
139  fvMatrix<scalar>& eqn
140 ) const
141 {
142  if (debug)
143  {
144  Info<< type() << ": applying source to " << eqn.psi().name() << endl;
145  }
146 
147  if (&he == &thermo_.he())
148  {
149  eqn += clouds_.Sh(eqn.psi());
150  }
151  else
152  {
154  << "Support for field " << he.name() << " is not implemented"
155  << exit(FatalError);
156  }
157 }
158 
159 
161 (
162  const volScalarField& rho,
163  const volVectorField& U,
164  fvMatrix<vector>& eqn
165 ) const
166 {
167  if (debug)
168  {
169  Info<< type() << ": applying source to " << eqn.psi().name() << endl;
170  }
171 
172  if (U.name() == "U")
173  {
174  eqn += clouds_.SU(eqn.psi());
175  }
176  else
177  {
179  << "Support for field " << U.name() << " is not implemented"
180  << exit(FatalError);
181  }
182 }
183 
184 
186 {
187  // Store the particle positions
188  clouds_.storeGlobalPositions();
189 }
190 
191 
193 {
194  clouds_.topoChange(map);
195 }
196 
197 
199 {
200  clouds_.mapMesh(map);
201 }
202 
203 
205 {
206  clouds_.distribute(map);
207 }
208 
209 
211 {
212  return true;
213 }
214 
215 
216 // ************************************************************************* //
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
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Base-class for fluid thermodynamic properties.
Definition: fluidThermo.H:57
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:118
VolField< Type > & psi()
Definition: fvMatrix.H:289
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
Finite volume model abstract base class.
Definition: fvModel.H:59
Lagrangian clouds model for VoF simulations.
Definition: VoFClouds.H:65
virtual bool movePoints()
Update for mesh motion.
Definition: VoFClouds.C:210
virtual wordList addSupFields() const
Return the list of fields for which the option adds source term.
Definition: VoFClouds.C:90
virtual void correct()
Solve the Lagrangian clouds and update the sources.
Definition: VoFClouds.C:96
virtual void topoChange(const polyTopoChangeMap &)
Update topology using the given map.
Definition: VoFClouds.C:192
virtual void distribute(const polyDistributionMap &)
Redistribute or update using the given distribution map.
Definition: VoFClouds.C:204
VoFClouds(const word &sourceName, const word &modelType, const fvMesh &mesh, const dictionary &dict)
Construct from explicit source name and mesh.
Definition: VoFClouds.C:53
virtual void preUpdateMesh()
Prepare for mesh_.update.
Definition: VoFClouds.C:185
virtual void mapMesh(const polyMeshMap &)
Update from another mesh using the given map.
Definition: VoFClouds.C:198
virtual void addSup(const volScalarField &alpha, const volScalarField &rho, fvMatrix< scalar > &eqn) const
Add explicit contribution to phase continuity.
Definition: VoFClouds.C:110
A base class for physical properties.
Class containing mesh-to-mesh mapping information after a mesh distribution where we send parts of me...
Class containing mesh-to-mesh mapping information.
Definition: polyMeshMap.H:51
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
A class for handling words, derived from string.
Definition: word.H:62
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
Calculate the matrix for implicit and explicit sources.
U
Definition: pEqn.H:72
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
addToRunTimeSelectionTable(fvConstraint, bound, dictionary)
defineTypeNameAndDebug(bound, 0)
Namespace for OpenFOAM.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
List< word > wordList
A List of words.
Definition: fileName.H:54
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
messageStream Info
error FatalError
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
label timeIndex
Definition: getTimeIndex.H:4
thermo he()
labelList fv(nPoints)
dictionary dict