phaseInterface.H
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 Class
25  Foam::phaseInterface
26 
27 Description
28  Class to represent an interface between phases. Derivations can further
29  specify the configuration of that interface; e.g., representing dispersal,
30  displacement or sidedness.
31 
32 SourceFiles
33  phaseInterface.C
34  phaseInterfaceI.H
35 
36 \*---------------------------------------------------------------------------*/
37 
38 #ifndef phaseInterface_H
39 #define phaseInterface_H
40 
41 #include "phaseModel.H"
42 #include "compressibleTwoPhases.H"
44 #include "runTimeSelectionTables.H"
45 
46 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
47 
48 namespace Foam
49 {
50 
51 class phaseInterfaceKey;
52 
53 /*---------------------------------------------------------------------------*\
54  Class phaseInterface Declaration
55 \*---------------------------------------------------------------------------*/
56 
57 class phaseInterface
58 :
60 {
61 private:
62 
63  // Private Static Data
64 
65  //- List of head separators. This is the first separator in an
66  // interface name, and it delimits the names of phases on either side
67  // of the interface. Other separators are modifiers (e.g.,
68  // displacedBy, inThe, ...).
69  static wordList headSeparators_;
70 
71  //- Map from old-format separators to separators
72  static HashTable<word> oldSeparatorToSeparator_;
73 
74 
75 public:
76 
77  // Public Static Functions
78 
79  //- Get a reference to phase1 after sorting the phases by index
80  static const phaseModel& getPhase1
81  (
82  const phaseModel& phase1,
83  const phaseModel& phase2
84  );
85 
86  //- Get a reference to phase2 after sorting the phases by index
87  static const phaseModel& getPhase2
88  (
89  const phaseModel& phase1,
90  const phaseModel& phase2
91  );
92 
93  //- Add a head separator to the list
94  static bool addHeadSeparator(const word& separator);
95 
96  //- Add a old separator to separator to the table
97  static bool addOldSeparatorToSeparator
98  (
99  const word& oldSeparator,
100  const word& separator
101  );
102 
103  //- Split an interface name and return all its parts
105  (
106  const phaseSystem& fluid,
107  const word& name
108  );
109 
110  //- Split an interface name and return its separators
112  (
113  const phaseSystem& fluid,
114  const word& name
115  );
116 
117  //- Convert a list of separators into a type name
119  (
120  const wordList& separators
121  );
122 
123  //- Convert an interface name into a type name. Essentially just
124  // replaces valid phase names with a "<phase>" placeholder.
125  static word nameToTypeName
126  (
127  const phaseSystem& fluid,
128  const word& name
129  );
130 
131  //- Convert interface name parts to an interface name
132  static word namePartsToName
133  (
134  const phaseSystem& fluid,
135  const wordList& nameParts
136  );
137 
138  //- Convert old-format interface name parts to an interface name. Used
139  // in phaseSystem to provide backwards compatible input.
140  static word oldNamePartsToName
141  (
142  const phaseSystem& fluid,
143  const wordList& oldNameParts
144  );
145 
146  //- Return references to the phases associated with a given name, and a
147  // list of valid separators
149  (
150  const phaseSystem& fluid,
151  const word& name,
152  const wordList& separators
153  );
154 
155 
156 private:
157 
158  // Private Data
159 
160  //- Phase 1
161  const phaseModel& phase1_;
162 
163  //- Phase 2
164  const phaseModel& phase2_;
165 
166  //- Gravitational acceleration
168 
169 
170 public:
171 
172  //- Runtime type information
173  TypeName("phaseInterface");
174 
175 
176  // Declare runtime construction
177 
179  (
180  autoPtr,
182  word,
183  (
184  const phaseSystem& fluid,
185  const word& name
186  ),
187  (fluid, name)
188  );
189 
190 
191  // Constructors
192 
193  //- Construct from phases
195  (
196  const phaseModel& phase1,
197  const phaseModel& phase2
198  );
199 
200  //- Construct from phases
202  (
204  );
205 
206  //- Construct from fluid and name
208  (
209  const phaseSystem& fluid,
210  const word& name
211  );
212 
213  //- Construct from fluid and key
215  (
216  const phaseSystem& fluid,
217  const phaseInterfaceKey& name
218  );
219 
220  //- Clone function
221  virtual autoPtr<phaseInterface> clone() const;
222 
223 
224  //- Destructor
225  virtual ~phaseInterface();
226 
227 
228  // Selectors
229 
230  //- Select given fluid and name
232  (
233  const phaseSystem& fluid,
234  const word& name
235  );
236 
237  //- Select by combining two interfaces
239  (
240  const phaseInterface& interface1,
241  const phaseInterface& interface2
242  );
243 
244  //- Class used for construction of PtrLists of phaseInterfaces
245  class iNew
246  {
247  const phaseSystem& fluid_;
248 
249  public:
250 
251  iNew(const phaseSystem& fluid)
252  :
253  fluid_(fluid)
254  {}
255 
257  {
258  return phaseInterface::New(fluid_, word(is));
259  }
260  };
261 
262 
263  // Static Member Functions
264 
265  //- Return the separator that delimits this interface's name
266  static word separator()
267  {
268  return word::null;
269  }
270 
271 
272  // Member Functions
273 
274  //- Name
275  virtual word name() const;
276 
277  //- Cast to derived type for use in a model
278  template<class ModelType, class Derived>
279  const Derived& modelCast() const
280  {
281  if (!isA<Derived>(*this))
282  {
284  << "Constructing " << ModelType::typeName
285  << " for interface " << name()
286  << " which is not of the required type "
287  << Derived::typeName << exit(FatalError);
288  }
289 
290  return refCast<const Derived>(*this);
291  }
292 
293 
294  // Access
295 
296  //- Return the volume fraction of phase 1
297  virtual inline const volScalarField& alpha1() const;
298 
299  //- Return the volume fraction of phase 2
300  virtual inline const volScalarField& alpha2() const;
301 
302  //- Return the thermo for phase 1
303  virtual inline const rhoFluidThermo& thermo1() const;
304 
305  //- Return the thermo for phase 2
306  virtual inline const rhoFluidThermo& thermo2() const;
307 
308  //- Return the density of phase 1
309  virtual inline const volScalarField& rho1() const;
310 
311  //- Return the density of phase 2
312  virtual inline const volScalarField& rho2() const;
313 
314  //- Return phase 1
315  inline const phaseModel& phase1() const;
316 
317  //- Return phase 2
318  inline const phaseModel& phase2() const;
319 
320  //- Return true if this phaseInterface contains the given phase
321  inline bool contains(const phaseModel& phase) const;
322 
323  //- Return the other phase relative to the given phase
324  // Generates a FatalError if this phaseInterface does not contain
325  // the given phase
326  inline const phaseModel& otherPhase(const phaseModel& phase) const;
327 
328  //- Return the index of the given phase. Generates a FatalError if
329  // this phaseInterface does not contain the given phase
330  inline label index(const phaseModel& phase) const;
331 
332  //- Return the phase system
333  inline const phaseSystem& fluid() const;
334 
335  //- Return the mesh
336  inline const fvMesh& mesh() const;
337 
338  //- Return gravitational acceleration
339  inline const uniformDimensionedVectorField& g() const;
340 
341 
342  // Properties
343 
344  //- Average density
345  tmp<volScalarField> rho() const;
346 
347  //- Relative velocity magnitude
348  tmp<volScalarField> magUr() const;
349 
350  //- Surface tension coefficient
351  tmp<volScalarField> sigma() const;
352 
353 
354  //- STL const_iterator
355  class const_iterator
356  {
357  // Private Data
358 
359  //- Reference to the pair for which this is an iterator
360  const phaseInterface& pair_;
361 
362  //- Current index
363  label index_;
364 
365  //- Construct an iterator with the given index
366  inline const_iterator(const phaseInterface&, const label index);
367 
368  public:
369 
370  friend class phaseInterface;
371 
372  // Constructors
373 
374  //- Construct from pair, moving to its 'begin' position
375  inline explicit const_iterator(const phaseInterface&);
376 
377 
378  // Access
379 
380  //- Return the current index
381  inline label index() const;
382 
383 
384  // Member Operators
385 
386  inline bool operator==(const const_iterator&) const;
387 
388  inline bool operator!=(const const_iterator&) const;
389 
390  inline const phaseModel& operator*() const;
391  inline const phaseModel& operator()() const;
392 
393  inline const phaseModel& otherPhase() const;
394 
395  inline const_iterator& operator++();
396  inline const_iterator operator++(int);
397  };
398 
399 
400  //- const_iterator set to the beginning of the pair
401  inline const_iterator cbegin() const;
402 
403  //- const_iterator set to beyond the end of the pair
404  inline const_iterator cend() const;
405 
406  //- const_iterator set to the beginning of the pair
407  inline const_iterator begin() const;
408 
409  //- const_iterator set to beyond the end of the pair
410  inline const_iterator end() const;
411 };
412 
413 
414 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
415 
416 } // End namespace Foam
417 
418 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
419 
420 #include "phaseInterfaceI.H"
421 
422 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
423 
424 #endif
425 
426 // ************************************************************************* //
Generic GeometricField class.
An STL-conforming hash table.
Definition: HashTable.H:127
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:60
A 2-tuple for storing two objects of different types.
Definition: Tuple2.H:66
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
Word-pair based class used for keying interface models in hash tables.
const phaseModel & operator*() const
bool operator!=(const const_iterator &) const
const phaseModel & otherPhase() const
const phaseModel & operator()() const
bool operator==(const const_iterator &) const
label index() const
Return the current index.
Class used for construction of PtrLists of phaseInterfaces.
iNew(const phaseSystem &fluid)
autoPtr< phaseInterface > operator()(Istream &is) const
Class to represent an interface between phases. Derivations can further specify the configuration of ...
static word oldNamePartsToName(const phaseSystem &fluid, const wordList &oldNameParts)
Convert old-format interface name parts to an interface name. Used.
static const phaseModel & getPhase1(const phaseModel &phase1, const phaseModel &phase2)
Get a reference to phase1 after sorting the phases by index.
static Tuple2< const phaseModel &, const phaseModel & > identifyPhases(const phaseSystem &fluid, const word &name, const wordList &separators)
Return references to the phases associated with a given name, and a.
const uniformDimensionedVectorField & g() const
Return gravitational acceleration.
const_iterator end() const
const_iterator set to beyond the end of the pair
const phaseSystem & fluid() const
Return the phase system.
TypeName("phaseInterface")
Runtime type information.
const phaseModel & otherPhase(const phaseModel &phase) const
Return the other phase relative to the given phase.
static word separator()
Return the separator that delimits this interface's name.
virtual autoPtr< phaseInterface > clone() const
Clone function.
virtual word name() const
Name.
static wordList nameToNameParts(const phaseSystem &fluid, const word &name)
Split an interface name and return all its parts.
static word nameToTypeName(const phaseSystem &fluid, const word &name)
Convert an interface name into a type name. Essentially just.
virtual const rhoFluidThermo & thermo2() const
Return the thermo for phase 2.
const_iterator begin() const
const_iterator set to the beginning of the pair
virtual const rhoFluidThermo & thermo1() const
Return the thermo for phase 1.
const_iterator cend() const
const_iterator set to beyond the end of the pair
static bool addHeadSeparator(const word &separator)
Add a head separator to the list.
virtual ~phaseInterface()
Destructor.
tmp< volScalarField > magUr() const
Relative velocity magnitude.
static word separatorsToTypeName(const wordList &separators)
Convert a list of separators into a type name.
static const phaseModel & getPhase2(const phaseModel &phase1, const phaseModel &phase2)
Get a reference to phase2 after sorting the phases by index.
virtual const volScalarField & rho2() const
Return the density of phase 2.
declareRunTimeSelectionTable(autoPtr, phaseInterface, word,(const phaseSystem &fluid, const word &name),(fluid, name))
tmp< volScalarField > sigma() const
Surface tension coefficient.
static wordList nameToSeparators(const phaseSystem &fluid, const word &name)
Split an interface name and return its separators.
phaseInterface(const phaseModel &phase1, const phaseModel &phase2)
Construct from phases.
label index(const phaseModel &phase) const
Return the index of the given phase. Generates a FatalError if.
static bool addOldSeparatorToSeparator(const word &oldSeparator, const word &separator)
Add a old separator to separator to the table.
const_iterator cbegin() const
const_iterator set to the beginning of the pair
const fvMesh & mesh() const
Return the mesh.
bool contains(const phaseModel &phase) const
Return true if this phaseInterface contains the given phase.
tmp< volScalarField > rho() const
Average density.
virtual const volScalarField & rho1() const
Return the density of phase 1.
const Derived & modelCast() const
Cast to derived type for use in a model.
static word namePartsToName(const phaseSystem &fluid, const wordList &nameParts)
Convert interface name parts to an interface name.
virtual const volScalarField & alpha2() const
Return the volume fraction of phase 2.
static autoPtr< phaseInterface > New(const phaseSystem &fluid, const word &name)
Select given fluid and name.
const phaseModel & phase1() const
Return phase 1.
const phaseModel & phase2() const
Return phase 2.
virtual const volScalarField & alpha1() const
Return the volume fraction of phase 1.
Class to represent a system of phases and model interfacial transfers between them.
Definition: phaseSystem.H:73
Base-class for fluid thermodynamic properties based on density.
A class for managing temporary objects.
Definition: tmp.H:55
A class for handling words, derived from string.
Definition: word.H:62
static const word null
An empty word.
Definition: word.H:77
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
Namespace for OpenFOAM.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
error FatalError
Macros to ease declaration of run-time selection tables.