All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
phaseSystem.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) 2015-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::phaseSystem
26 
27 Description
28  Class to represent a system of phases and model interfacial transfers
29  between them.
30 
31 SourceFiles
32  phaseSystem.C
33 
34 \*---------------------------------------------------------------------------*/
35 
36 #ifndef phaseSystem_H
37 #define phaseSystem_H
38 
39 #include "IOdictionary.H"
40 
41 #include "phaseModel.H"
42 #include "phaseInterface.H"
43 #include "phaseInterfaceKey.H"
44 #include "HashPtrTable.H"
45 #include "PtrListDictionary.H"
46 #include "hashedWordList.H"
47 
48 #include "pimpleNoLoopControl.H"
49 
50 #include "IOMRFZoneList.H"
51 #include "fvModels.H"
52 #include "fvConstraints.H"
53 
54 #include "volFields.H"
55 #include "surfaceFields.H"
56 #include "fvMatricesFwd.H"
57 
58 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
59 
60 namespace Foam
61 {
62 
63 class interfaceSurfaceTensionModel;
64 class pressureReference;
65 class nonOrthogonalSolutionControl;
66 
67 /*---------------------------------------------------------------------------*\
68  Class phaseSystem Declaration
69 \*---------------------------------------------------------------------------*/
70 
71 class phaseSystem
72 :
73  public IOdictionary
74 {
75 public:
76 
77  // Public Typedefs
78 
80 
82 
84 
86 
88 
89  typedef
91  <
95  >
96  dmdtfTable;
97 
98  typedef
100  <
104  >
105  dmidtfTable;
106 
107 
108 
109 protected:
110 
111  // Protected typedefs
112 
113  typedef
114  HashTable
115  <
119  >
121 
123  cAlphaTable;
124 
125 
126  // Protected data
127 
128  //- Reference to the mesh
129  const fvMesh& mesh_;
130 
131  //- Reference to pimpleNoLoopControl
133 
134  //- Optional MRF zones
136 
137  //- Name of optional reference phase which is not solved for
138  // but obtained from the sum of the other phases
140 
141  //- Phase models
143 
144  //- Moving phase models
146 
147  //- Stationary phase models
149 
150  //- Anisothermal phase models
152 
153  //- Multi-component phase models
155 
156  //- Total volumetric flux
158 
159  //- Rate of change of pressure
161 
162  //- Interface compression coefficients
164 
165  //- Stabilisation for normalisation of the interface normal
167 
168 
169  // Sub Models
170 
171  //- Surface tension models
173 
174 
175  // Protected member functions
176 
177  //- Return the sum of the phase fractions of the moving phases
179 
180  //- Re-normalise the velocity of the phases
181  // around the specified mixture mean
182  void setMixtureU(const volVectorField& Um);
183 
184 
185  // Functions required for interface compression
186 
187  //- Normal to interface between two phases
188  // Used for interface compression
190  (
191  const volScalarField& alpha1,
192  const volScalarField& alpha2
193  ) const;
194 
195  //- Normal to interface between two phases dotted with face areas
196  // Used for interface compression
198  (
199  const volScalarField& alpha1,
200  const volScalarField& alpha2
201  ) const;
202 
203  //- Curvature of interface between two phases
204  // Used for interface compression
206  (
207  const phaseModel& alpha1,
208  const phaseModel& alpha2
209  ) const;
210 
211 
212  // Sub-model construction
213 
214  //- Return the dictionary containing interfacial model or value
215  // settings for the given name. Performs backwards compatibility
216  // conversions.
217  template<class Type>
218  dictionary interfacialDict(const word& name) const;
219 
220 
221 public:
222 
223  //- Runtime type information
224  TypeName("phaseSystem");
225 
226  //- Default name of the phase properties dictionary
227  static const word propertiesName;
228 
229 
230  // Declare runtime construction
231 
233  (
234  autoPtr,
235  phaseSystem,
236  dictionary,
237  (
238  const fvMesh& mesh
239  ),
240  (mesh)
241  );
242 
243 
244  // Constructors
245 
246  //- Construct from fvMesh
247  phaseSystem(const fvMesh& mesh);
248 
249 
250  // Selectors
251 
253  (
254  const fvMesh& mesh
255  );
256 
257 
258  //- Destructor
259  virtual ~phaseSystem();
260 
261 
262  // Member Functions
263 
264  // Access
265 
266  //- Return the mesh
267  inline const fvMesh& mesh() const;
268 
269  //- Return pimpleNoLoopControl
270  inline const pimpleNoLoopControl& pimple() const;
271 
272  //- Return the phase models
273  inline const phaseModelList& phases() const;
274 
275  //- Access the phase models
276  inline phaseModelList& phases();
277 
278  //- Return the models for phases that are moving
279  inline const phaseModelPartialList& movingPhases() const;
280 
281  //- Access the models for phases that are moving
283 
284  //- Return the models for phases that are stationary
285  inline const phaseModelPartialList& stationaryPhases() const;
286 
287  //- Access the models for phases that are stationary
289 
290  //- Return the models for phases that have variable temperature
291  inline const phaseModelPartialList& anisothermalPhases() const;
292 
293  //- Access the models for phases that have variable temperature
295 
296  //- Return the models for phases that have multiple species
297  inline const phaseModelPartialList& multicomponentPhases() const;
298 
299  //- Access the models for phases that have multiple species
301 
302  //- Return the phase not given as an argument in a two-phase system
303  // An error is generated if the system is not two-phase
304  inline const phaseModel& otherPhase(const phaseModel& phase) const;
305 
306  //- Return the mixture flux
307  inline const surfaceScalarField& phi() const;
308 
309  //- Access the mixture flux
310  inline surfaceScalarField& phi();
311 
312  //- Return the rate of change of the pressure
313  inline const volScalarField& dpdt() const;
314 
315  //- Access the rate of change of the pressure
316  inline volScalarField& dpdt();
317 
318  //- Return MRF zones
319  inline const IOMRFZoneList& MRF() const;
320 
321  //- Access the fvModels
323 
324  //- Access the fvModels
325  inline const Foam::fvModels& fvModels() const;
326 
327  //- Access the fvConstraints
329 
330  //- Access the fvConstraints
331  inline const Foam::fvConstraints& fvConstraints() const;
332 
333 
334  // Sub-model construction
335 
336  //- Return the model name. This is the same as the model's typename
337  // but without "Model" on the end.
338  template<class ModelType>
339  word modelName() const;
340 
341  //- Generate interfacial-model lists
342  template<class ModelType, class ... InterfaceTypes>
344  (
345  const dictionary& dict,
346  const phaseInterface& interface,
347  PtrList<phaseInterface>& interfaces,
348  PtrList<ModelType>& models
349  ) const;
350 
351  //- Generate interfacial-model tables
352  template<class ModelType>
354  (
355  const dictionary& dict,
356  HashTable
357  <
361  >& models
362  ) const;
363 
364  //- Generate interfacial-model tables
365  template<class ModelType>
367  (
368  HashTable
369  <
373  >& models
374  ) const;
375 
376  //- Generate interfacial-model tables
377  template<class ValueType>
379  (
380  const dictionary& dict,
381  HashTable
382  <
383  ValueType,
386  >& values
387  ) const;
388 
389  //- Generate interfacial-model tables
390  template<class ValueType>
392  (
393  const word& valueName,
394  HashTable
395  <
396  ValueType,
399  >& values
400  ) const;
401 
402  //- Return the dictionary from which to construct a low-level
403  // sub-model. Checks that there is just one sub-dictionary then
404  // returns it.
405  template<class ModelType>
406  static const dictionary& modelSubDict(const dictionary& dict);
407 
408  //- Check that mass transfer is supported across the given interface
409  template<class ModelType>
410  void validateMassTransfer(const phaseInterface& interface) const;
411 
412 
413  // Sub-model lookup
414 
415  //- Check availability of a sub model for a given interface
416  template<class ModelType>
417  bool foundInterfacialModel(const phaseInterface& interface) const;
418 
419  //- Return a sub model for an interface
420  template<class ModelType>
421  const ModelType& lookupInterfacialModel
422  (
423  const phaseInterface& interface
424  ) const;
425 
426 
427  // Properties
428 
429  //- Return the mixture density
430  tmp<volScalarField> rho() const;
431 
432  //- Return the mixture velocity
433  tmp<volVectorField> U() const;
434 
435  //- Return the surface tension coefficient for an interface
436  tmp<volScalarField> sigma(const phaseInterfaceKey& key) const;
437 
438  //- Return the surface tension coefficient for an interface on a
439  // patch
441  (
442  const phaseInterfaceKey& key,
443  const label patchi
444  ) const;
445 
446  //- Indicator of the proximity of the interface
447  // Field values are 1 near and 0 away for the interface.
449 
450  //- Stabilisation for normalisation of the interface normal
451  inline const dimensionedScalar& deltaN() const;
452 
453  //- Return the mass transfer rate for an interface
454  virtual tmp<volScalarField> dmdtf
455  (
456  const phaseInterfaceKey& key
457  ) const;
458 
459  //- Return the mass transfer rates for each phase
460  virtual PtrList<volScalarField> dmdts() const;
461 
462  //- Return the mass transfer pressure implicit coefficients
463  // for each phase
464  virtual PtrList<volScalarField> d2mdtdps() const;
465 
466  //- Return incompressibility
467  bool incompressible() const;
468 
469 
470  // Transfers
471 
472  //- Return the momentum transfer matrices for the cell-based
473  // algorithm
475 
476  //- Return the momentum transfer matrices for the face-based
477  // algorithm
479 
480  //- Return the force fluxes for the cell-based algorithm
481  virtual PtrList<surfaceScalarField> Fs() const = 0;
482 
483  //- Return the force fluxes for the face-based algorithm
484  virtual PtrList<surfaceScalarField> Ffs() const = 0;
485 
486  //- Return the inverse of the central + drag + virtual mass
487  // coefficient matrix
488  virtual void invADVs
489  (
490  const PtrList<volScalarField>& As,
494  ) const = 0;
495 
496  //- Return the inverse of the face central + drag + virtual mass
497  // coefficient matrix
499  (
500  const PtrList<surfaceScalarField>& Afs,
502  ) const = 0;
503 
504  //- Returns true if the phase pressure is treated implicitly
505  // in the phase fraction equation
506  virtual bool implicitPhasePressure(const phaseModel& phase) const;
507 
508  //- Returns true if the phase pressure is treated implicitly
509  // in the phase fraction equation for any phase
510  virtual bool implicitPhasePressure() const;
511 
512  //- Return the phase diffusivity
513  // divided by the momentum central coefficient
515  (
516  const PtrList<volScalarField>& rAs
517  ) const = 0;
518 
519  //- Return the flux corrections for the cell-based algorithm
520  virtual PtrList<surfaceScalarField> ddtCorrs() const = 0;
521 
522  //- Set the cell and faces drag correction fields
523  virtual void dragCorrs
524  (
526  PtrList<surfaceScalarField>& dragCorrf
527  ) const = 0;
528 
529  //- Re-normalise the flux of the phases
530  // around the specified mixture mean
531  void setMixturePhi
532  (
533  const PtrList<surfaceScalarField>& alphafs,
534  const surfaceScalarField& phim
535  );
536 
537  //- Return the heat transfer matrices
538  virtual autoPtr<heatTransferTable> heatTransfer() const = 0;
539 
540  //- Return the specie transfer matrices
541  virtual autoPtr<specieTransferTable> specieTransfer() const = 0;
542 
543  //- Return the surface tension force
545  (
546  const phaseModel& phase
547  ) const;
548 
549 
550  // Evolution
551 
552  //- Solve for the phase fractions
553  virtual void solve(const PtrList<volScalarField>& rAs);
554 
555  //- Correct the fluid properties other than those listed below
556  virtual void correct();
557 
558  //- Correct the continuity errors
559  virtual void correctContinuityError();
560 
561  //- Correct the kinematics
562  virtual void correctKinematics();
563 
564  //- Correct the thermodynamics
565  virtual void correctThermo();
566 
567  //- Correct the reactions
568  virtual void correctReactions();
569 
570  //- Correct the species mass fractions
571  virtual void correctSpecies();
572 
573  //- Predict the momentumTransport
574  virtual void predictMomentumTransport();
575 
576  //- Predict the energy transport e.g. alphat
577  virtual void predictThermophysicalTransport();
578 
579  //- Correct the momentumTransport
580  virtual void correctMomentumTransport();
581 
582  //- Correct the energy transport e.g. alphat
583  virtual void correctThermophysicalTransport();
584 
585  //- Update the fluid properties for mesh changes
586  virtual void meshUpdate();
587 
588  //- Correct fixed-flux BCs to be consistent with the velocity BCs
589  void correctBoundaryFlux();
590 
591  void correctPhi
592  (
593  const volScalarField& p_rgh,
594  const autoPtr<volScalarField>& divU,
597  );
598 
599 
600  // IO
601 
602  //- Read base phaseProperties dictionary
603  virtual bool read();
604 };
605 
606 
607 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
608 
611 
612 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
613 
614 } // End namespace Foam
615 
616 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
617 
618 #include "phaseSystemI.H"
619 
620 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
621 
622 #ifdef NoRepository
623  #include "phaseSystemTemplates.C"
624 #endif
625 
626 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
627 
628 #endif
629 
630 // ************************************************************************* //
Generic GeometricField class.
A HashTable specialisation for hashing pointers.
Definition: HashPtrTable.H:67
An STL-conforming hash table.
Definition: HashTable.H:127
List of MRF zones with IO functionality. MRF zones are specified by a list of dictionary entries,...
Definition: IOMRFZoneList.H:68
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:57
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: PtrList.H:75
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 fileName & name() const
Return the dictionary name.
Definition: dictionary.H:111
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Finite volume constraints.
Definition: fvConstraints.H:67
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
Finite volume models.
Definition: fvModels.H:65
Non-orthogonal solution control class. Provides non-orthogonal-loop control methods.
Word-pair based class used for keying interface models in hash tables.
Class to represent an interface between phases. Derivations can further specify the configuration of ...
Class to represent a system of phases and model interfacial transfers between them.
Definition: phaseSystem.H:73
phaseModelPartialList stationaryPhaseModels_
Stationary phase models.
Definition: phaseSystem.H:147
virtual autoPtr< specieTransferTable > specieTransfer() const =0
Return the specie transfer matrices.
virtual void correctKinematics()
Correct the kinematics.
Definition: phaseSystem.C:595
tmp< volScalarField > nearInterface() const
Indicator of the proximity of the interface.
Definition: phaseSystem.C:435
virtual PtrList< surfaceScalarField > Ffs() const =0
Return the force fluxes for the face-based algorithm.
dictionary interfacialDict(const word &name) const
Return the dictionary containing interfacial model or value.
void correctPhi(const volScalarField &p_rgh, const autoPtr< volScalarField > &divU, const pressureReference &pressureReference, nonOrthogonalSolutionControl &pimple)
Definition: phaseSystem.C:722
bool foundInterfacialModel(const phaseInterface &interface) const
Check availability of a sub model for a given interface.
virtual tmp< surfaceScalarField > alphaDByAf(const PtrList< volScalarField > &rAs) const =0
Return the phase diffusivity.
static const word propertiesName
Default name of the phase properties dictionary.
Definition: phaseSystem.H:226
void validateMassTransfer(const phaseInterface &interface) const
Check that mass transfer is supported across the given interface.
void generateInterfacialModels(const dictionary &dict, const phaseInterface &interface, PtrList< phaseInterface > &interfaces, PtrList< ModelType > &models) const
Generate interfacial-model lists.
const IOMRFZoneList & MRF() const
Return MRF zones.
Definition: phaseSystemI.H:218
const phaseModel & otherPhase(const phaseModel &phase) const
Return the phase not given as an argument in a two-phase system.
Definition: phaseSystemI.H:172
virtual void correctReactions()
Correct the reactions.
Definition: phaseSystem.C:622
const fvMesh & mesh_
Reference to the mesh.
Definition: phaseSystem.H:128
phaseModelPartialList multicomponentPhaseModels_
Multi-component phase models.
Definition: phaseSystem.H:153
phaseModelList phaseModels_
Phase models.
Definition: phaseSystem.H:141
phaseModelPartialList movingPhaseModels_
Moving phase models.
Definition: phaseSystem.H:144
void setMixtureU(const volVectorField &Um)
Re-normalise the velocity of the phases.
Definition: phaseSystem.C:80
const dimensionedScalar & deltaN() const
Stabilisation for normalisation of the interface normal.
Definition: phaseSystemI.H:248
virtual PtrList< surfaceScalarField > ddtCorrs() const =0
Return the flux corrections for the cell-based algorithm.
virtual void correctSpecies()
Correct the species mass fractions.
Definition: phaseSystem.C:631
const Foam::fvModels & fvModels() const
Access the fvModels.
Definition: phaseSystemI.H:230
const surfaceScalarField & phi() const
Return the mixture flux.
Definition: phaseSystemI.H:194
virtual void correctThermo()
Correct the thermodynamics.
Definition: phaseSystem.C:613
const pimpleNoLoopControl & pimple_
Reference to pimpleNoLoopControl.
Definition: phaseSystem.H:131
const volScalarField & dpdt() const
Return the rate of change of the pressure.
Definition: phaseSystemI.H:206
virtual void correct()
Correct the fluid properties other than those listed below.
Definition: phaseSystem.C:551
tmp< volScalarField > sigma(const phaseInterfaceKey &key) const
Return the surface tension coefficient for an interface.
Definition: phaseSystem.C:399
const phaseModelPartialList & multicomponentPhases() const
Return the models for phases that have multiple species.
Definition: phaseSystemI.H:158
virtual autoPtr< heatTransferTable > heatTransfer() const =0
Return the heat transfer matrices.
phaseModelPartialList anisothermalPhaseModels_
Anisothermal phase models.
Definition: phaseSystem.H:150
const phaseModelPartialList & movingPhases() const
Return the models for phases that are moving.
Definition: phaseSystemI.H:116
virtual void predictThermophysicalTransport()
Predict the energy transport e.g. alphat.
Definition: phaseSystem.C:649
tmp< volScalarField > K(const phaseModel &alpha1, const phaseModel &alpha2) const
Curvature of interface between two phases.
Definition: phaseSystem.C:162
virtual void solve(const PtrList< volScalarField > &rAs)
Solve for the phase fractions.
HashTable< scalar, phaseInterfaceKey, phaseInterfaceKey::hash > cAlphaTable
Definition: phaseSystem.H:122
tmp< volVectorField > U() const
Return the mixture velocity.
Definition: phaseSystem.C:371
virtual void correctMomentumTransport()
Correct the momentumTransport.
Definition: phaseSystem.C:658
void correctBoundaryFlux()
Correct fixed-flux BCs to be consistent with the velocity BCs.
Definition: phaseSystem.C:690
static autoPtr< phaseSystem > New(const fvMesh &mesh)
const phaseModelList & phases() const
Return the phase models.
Definition: phaseSystemI.H:102
PtrListDictionary< phaseModel > phaseModelList
Definition: phaseSystem.H:84
const phaseModelPartialList & anisothermalPhases() const
Return the models for phases that have variable temperature.
Definition: phaseSystemI.H:144
virtual void invADVs(const PtrList< volScalarField > &As, PtrList< volVectorField > &HVms, PtrList< PtrList< volScalarField >> &invADVs, PtrList< PtrList< surfaceScalarField >> &invADVfs) const =0
Return the inverse of the central + drag + virtual mass.
virtual tmp< volScalarField > dmdtf(const phaseInterfaceKey &key) const
Return the mass transfer rate for an interface.
Definition: phaseSystem.C:461
virtual void dragCorrs(PtrList< volVectorField > &dragCorrs, PtrList< surfaceScalarField > &dragCorrf) const =0
Set the cell and faces drag correction fields.
const dimensionedScalar deltaN_
Stabilisation for normalisation of the interface normal.
Definition: phaseSystem.H:165
interfaceSurfaceTensionModelTable interfaceSurfaceTensionModels_
Surface tension models.
Definition: phaseSystem.H:171
const pimpleNoLoopControl & pimple() const
Return pimpleNoLoopControl.
Definition: phaseSystemI.H:95
virtual autoPtr< momentumTransferTable > momentumTransferf()=0
Return the momentum transfer matrices for the face-based.
HashPtrTable< fvScalarMatrix > specieTransferTable
Definition: phaseSystem.H:82
virtual void correctThermophysicalTransport()
Correct the energy transport e.g. alphat.
Definition: phaseSystem.C:667
virtual PtrList< volScalarField > dmdts() const
Return the mass transfer rates for each phase.
Definition: phaseSystem.C:474
virtual void meshUpdate()
Update the fluid properties for mesh changes.
Definition: phaseSystem.C:676
cAlphaTable cAlphas_
Interface compression coefficients.
Definition: phaseSystem.H:162
virtual autoPtr< momentumTransferTable > momentumTransfer()=0
Return the momentum transfer matrices for the cell-based.
const phaseModelPartialList & stationaryPhases() const
Return the models for phases that are stationary.
Definition: phaseSystemI.H:130
word modelName() const
Return the model name. This is the same as the model's typename.
surfaceScalarField phi_
Total volumetric flux.
Definition: phaseSystem.H:156
HashPtrTable< fvVectorMatrix > momentumTransferTable
Definition: phaseSystem.H:78
virtual PtrList< surfaceScalarField > Fs() const =0
Return the force fluxes for the cell-based algorithm.
virtual bool implicitPhasePressure() const
Returns true if the phase pressure is treated implicitly.
Definition: phaseSystem.C:506
tmp< surfaceScalarField > surfaceTension(const phaseModel &phase) const
Return the surface tension force.
Definition: phaseSystem.C:513
tmp< volScalarField > sumAlphaMoving() const
Return the sum of the phase fractions of the moving phases.
Definition: phaseSystem.C:54
HashPtrTable< fvScalarMatrix > heatTransferTable
Definition: phaseSystem.H:80
phaseSystem(const fvMesh &mesh)
Construct from fvMesh.
Definition: phaseSystem.C:186
void setMixturePhi(const PtrList< surfaceScalarField > &alphafs, const surfaceScalarField &phim)
Re-normalise the flux of the phases.
Definition: phaseSystem.C:101
HashTable< autoPtr< interfaceSurfaceTensionModel >, phaseInterfaceKey, phaseInterfaceKey::hash > interfaceSurfaceTensionModelTable
Definition: phaseSystem.H:119
virtual ~phaseSystem()
Destructor.
Definition: phaseSystem.C:338
HashPtrTable< HashPtrTable< volScalarField >, phaseInterfaceKey, phaseInterfaceKey::hash > dmidtfTable
Definition: phaseSystem.H:104
void generateInterfacialValues(const dictionary &dict, HashTable< ValueType, phaseInterfaceKey, phaseInterfaceKey::hash > &values) const
Generate interfacial-model tables.
const fvMesh & mesh() const
Return the mesh.
Definition: phaseSystemI.H:89
TypeName("phaseSystem")
Runtime type information.
UPtrList< phaseModel > phaseModelPartialList
Definition: phaseSystem.H:86
tmp< surfaceScalarField > nHatf(const volScalarField &alpha1, const volScalarField &alpha2) const
Normal to interface between two phases dotted with face areas.
Definition: phaseSystem.C:151
tmp< volScalarField > rho() const
Return the mixture density.
Definition: phaseSystem.C:344
IOMRFZoneList MRF_
Optional MRF zones.
Definition: phaseSystem.H:134
HashPtrTable< volScalarField, phaseInterfaceKey, phaseInterfaceKey::hash > dmdtfTable
Definition: phaseSystem.H:95
const Foam::fvConstraints & fvConstraints() const
Access the fvConstraints.
Definition: phaseSystemI.H:242
word referencePhaseName_
Name of optional reference phase which is not solved for.
Definition: phaseSystem.H:138
tmp< surfaceVectorField > nHatfv(const volScalarField &alpha1, const volScalarField &alpha2) const
Normal to interface between two phases.
Definition: phaseSystem.C:125
volScalarField dpdt_
Rate of change of pressure.
Definition: phaseSystem.H:159
virtual void predictMomentumTransport()
Predict the momentumTransport.
Definition: phaseSystem.C:640
virtual PtrList< PtrList< surfaceScalarField > > invADVfs(const PtrList< surfaceScalarField > &Afs, PtrList< surfaceScalarField > &HVmfs) const =0
Return the inverse of the face central + drag + virtual mass.
static const dictionary & modelSubDict(const dictionary &dict)
Return the dictionary from which to construct a low-level.
declareRunTimeSelectionTable(autoPtr, phaseSystem, dictionary,(const fvMesh &mesh),(mesh))
virtual PtrList< volScalarField > d2mdtdps() const
Return the mass transfer pressure implicit coefficients.
Definition: phaseSystem.C:480
virtual bool read()
Read base phaseProperties dictionary.
Definition: phaseSystem.C:833
const ModelType & lookupInterfacialModel(const phaseInterface &interface) const
Return a sub model for an interface.
virtual void correctContinuityError()
Correct the continuity errors.
Definition: phaseSystem.C:560
bool incompressible() const
Return incompressibility.
Definition: phaseSystem.C:486
Pimple no-loop control class. Implements various option flags, but leaves loop controls to the deriva...
Provides controls for the pressure reference in closed-volume simulations.
A class for managing temporary objects.
Definition: tmp.H:55
A class for handling words, derived from string.
Definition: word.H:62
Forward declarations of fvMatrix specialisations.
label patchi
volScalarField sf(fieldObject, mesh)
word valueName(const direction argument)
Return the name of the value entry for the given argument.
Definition: Product2I.H:81
Namespace for OpenFOAM.
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
tmp< volScalarField > byDt(const volScalarField &vf)
Definition: phaseSystem.C:855
VolField< scalar > volScalarField
Definition: volFieldsFwd.H:64
dictionary dict
Foam::surfaceFields.