ThermalPhaseChangePhaseSystem.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) 2015-2024 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 "heatTransferModel.H"
29 #include "fvcVolumeIntegrate.H"
30 #include "fvmSup.H"
33 
34 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
35 
36 template<class BasePhaseSystem>
38 (
39  PtrList<volScalarField>& dmdts
40 ) const
41 {
42  forAllConstIter(phaseSystem::dmdtfTable, dmdtfs_, dmdtfIter)
43  {
44  const phaseInterface interface(*this, dmdtfIter.key());
45 
46  addField(interface.phase1(), "dmdt", *dmdtfIter(), dmdts);
47  addField(interface.phase2(), "dmdt", - *dmdtfIter(), dmdts);
48  }
49 
50  forAllConstIter(phaseSystem::dmdtfTable, nDmdtfs_, nDmdtfIter)
51  {
52  const phaseInterface interface(*this, nDmdtfIter.key());
53 
54  addField(interface.phase1(), "dmdt", *nDmdtfIter(), dmdts);
55  addField(interface.phase2(), "dmdt", - *nDmdtfIter(), dmdts);
56  }
57 }
58 
59 
60 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
61 
62 template<class BasePhaseSystem>
65 (
66  const fvMesh& mesh
67 )
68 :
69  BasePhaseSystem(mesh),
70  volatile_(this->template lookupOrDefault<word>("volatile", "none")),
71  dmdt0s_(this->phases().size()),
72  pressureImplicit_
73  (
74  this->template lookupOrDefault<Switch>("pressureImplicit", true)
75  )
76 {
77  this->generateInterfacialModels(saturationModels_);
78 
79  // Check that models have been specified in the correct combinations
81  (
83  saturationModels_,
84  saturationModelIter
85  )
86  {
87  const phaseInterface& interface = saturationModelIter()->interface();
88  const phaseModel& phase1 = interface.phase1();
89  const phaseModel& phase2 = interface.phase2();
90 
91  this->template validateMassTransfer
92  <
94  >(interface);
95 
96  if
97  (
98  !this->heatTransferModels_.found(interface)
99  || !this->heatTransferModels_[interface]->haveModelInThe(phase1)
100  || !this->heatTransferModels_[interface]->haveModelInThe(phase2)
101  )
102  {
104  << "A heat transfer model for both sides of the "
105  << interface.name() << " interface is not specified. This is "
106  << "required by the corresponding saturation model"
107  << exit(FatalError);
108  }
109  }
110 
111  // Generate interfacial mass transfer fields, initially assumed to be zero
113  (
115  saturationModels_,
116  saturationModelIter
117  )
118  {
119  const phaseInterface& interface = saturationModelIter()->interface();
120 
121  dmdtfs_.insert
122  (
123  interface,
124  new volScalarField
125  (
126  IOobject
127  (
129  (
130  "thermalPhaseChange:dmdtf",
131  interface.name()
132  ),
133  this->mesh().time().name(),
134  this->mesh(),
137  ),
138  this->mesh(),
140  )
141  );
142 
143  d2mdtdpfs_.insert
144  (
145  interface,
146  new volScalarField
147  (
148  IOobject
149  (
151  (
152  "thermalPhaseChange:d2mdtdpf",
153  interface.name()
154  ),
155  this->mesh().time().name(),
156  this->mesh(),
159  ),
160  this->mesh(),
162  )
163  );
164 
165  Tfs_.insert
166  (
167  interface,
168  new volScalarField
169  (
170  IOobject
171  (
173  (
174  "thermalPhaseChange:Tf",
175  interface.name()
176  ),
177  this->mesh().time().name(),
178  this->mesh(),
181  ),
182  (
183  interface.phase1().fluidThermo().T()
184  + interface.phase2().fluidThermo().T()
185  )/2
186  )
187  );
188 
189  Tsats_.insert
190  (
191  interface,
192  new volScalarField
193  (
194  IOobject
195  (
197  (
198  "thermalPhaseChange:Tsat",
199  interface.name()
200  ),
201  this->mesh().time().name(),
202  this->mesh(),
205  ),
206  saturationModels_[interface]->Tsat
207  (
208  interface.phase1().fluidThermo().p()
209  )
210  )
211  );
212 
213  nDmdtfs_.insert
214  (
215  interface,
216  new volScalarField
217  (
218  IOobject
219  (
221  (
222  "thermalPhaseChange:nucleation:dmdtf",
223  interface.name()
224  ),
225  this->mesh().time().name(),
226  this->mesh(),
229  ),
230  this->mesh(),
232  )
233  );
234  }
235 }
236 
237 
238 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
239 
240 template<class BasePhaseSystem>
243 {}
244 
245 
246 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
247 
248 template<class BasePhaseSystem>
251 (
252  const phaseInterfaceKey& key
253 ) const
254 {
255  return saturationModels_[key];
256 }
257 
258 
259 template<class BasePhaseSystem>
262 (
263  const phaseInterfaceKey& key
264 ) const
265 {
266  tmp<volScalarField> tDmdtf = BasePhaseSystem::dmdtf(key);
267 
268  if (dmdtfs_.found(key))
269  {
270  tDmdtf.ref() += *dmdtfs_[key];
271  }
272 
273  if (nDmdtfs_.found(key))
274  {
275  tDmdtf.ref() += *nDmdtfs_[key];
276  }
277 
278  return tDmdtf;
279 }
280 
281 
282 template<class BasePhaseSystem>
285 {
286  PtrList<volScalarField> dmdts(BasePhaseSystem::dmdts());
287 
288  addDmdts(dmdts);
289 
290  return dmdts;
291 }
292 
293 
294 template<class BasePhaseSystem>
297 {
298  PtrList<volScalarField> d2mdtdps(BasePhaseSystem::d2mdtdps());
299 
300  forAllConstIter(phaseSystem::dmdtfTable, d2mdtdpfs_, d2mdtdpfIter)
301  {
302  const phaseInterface interface(*this, d2mdtdpfIter.key());
303 
304  addField(interface.phase1(), "d2mdtdp", *d2mdtdpfIter(), d2mdtdps);
305  addField(interface.phase2(), "d2mdtdp", - *d2mdtdpfIter(), d2mdtdps);
306  }
307 
308  return d2mdtdps;
309 }
310 
311 
312 template<class BasePhaseSystem>
316 {
318  BasePhaseSystem::momentumTransfer();
319 
320  phaseSystem::momentumTransferTable& eqns = eqnsPtr();
321 
322  this->addDmdtUfs(dmdtfs_, eqns);
323  this->addDmdtUfs(nDmdtfs_, eqns);
324 
325  return eqnsPtr;
326 }
327 
328 
329 template<class BasePhaseSystem>
333 {
335  BasePhaseSystem::momentumTransferf();
336 
337  phaseSystem::momentumTransferTable& eqns = eqnsPtr();
338 
339  this->addDmdtUfs(dmdtfs_, eqns);
340  this->addDmdtUfs(nDmdtfs_, eqns);
341 
342  return eqnsPtr;
343 }
344 
345 
346 template<class BasePhaseSystem>
349 {
351  BasePhaseSystem::heatTransfer();
352 
353  phaseSystem::heatTransferTable& eqns = eqnsPtr();
354 
355  // Create temperatures at which to evaluate nucleation mass transfers
357  forAllConstIter(phaseSystem::dmdtfTable, nDmdtfs_, nDmdtfIter)
358  {
359  const phaseInterface interface(*this, nDmdtfIter.key());
360  const interfaceSaturationTemperatureModel& satModel =
361  this->saturation(nDmdtfIter.key());
362 
363  Tns.insert
364  (
365  interface,
366  satModel.Tsat(interface.phase1().fluidThermo().p()).ptr()
367  );
368  }
369 
370  // Mass transfer terms
371  if (volatile_ != "none")
372  {
373  {
374  phaseSystem::dmidtfTable dmidtfs;
375 
376  forAllConstIter(phaseSystem::dmdtfTable, dmdtfs_, dmdtfIter)
377  {
378  const phaseInterface interface(*this, dmdtfIter.key());
379 
380  dmidtfs.insert(interface, new HashPtrTable<volScalarField>());
381  dmidtfs[interface]->insert
382  (
383  volatile_,
384  new volScalarField(*dmdtfIter())
385  );
386  }
387 
388  this->addDmidtHefs
389  (
390  dmidtfs,
391  //Tfs_,
392  Tsats_,
393  latentHeatScheme::upwind,
394  latentHeatTransfer::mass,
395  eqns
396  );
397  }
398  {
399  phaseSystem::dmidtfTable nDmidtfs;
400 
401  forAllConstIter(phaseSystem::dmdtfTable, nDmdtfs_, nDmdtfIter)
402  {
403  const phaseInterface interface(*this, nDmdtfIter.key());
404 
405  nDmidtfs.insert(interface, new HashPtrTable<volScalarField>());
406  nDmidtfs[interface]->insert
407  (
408  volatile_,
409  new volScalarField(*nDmdtfIter())
410  );
411  }
412 
413  this->addDmidtHefs
414  (
415  nDmidtfs,
416  Tns,
417  0,
418  latentHeatScheme::upwind,
419  eqns
420  );
421  }
422  }
423  else
424  {
425  this->addDmdtHefs
426  (
427  dmdtfs_,
428  //Tfs_,
429  Tsats_,
430  latentHeatScheme::upwind,
431  latentHeatTransfer::mass,
432  eqns
433  );
434  this->addDmdtHefs
435  (
436  nDmdtfs_,
437  Tns,
438  0,
439  latentHeatScheme::upwind,
440  eqns
441  );
442  }
443 
444  // Lagging
445  {
446  PtrList<volScalarField> dmdts(this->phases().size());
447 
448  addDmdts(dmdts);
449 
450  forAll(this->phases(), phasei)
451  {
452  const phaseModel& phase = this->phases()[phasei];
453 
454  if (dmdt0s_.set(phase.index()))
455  {
456  *eqns[phase.name()] +=
457  fvm::Sp
458  (
459  dmdt0s_[phase.index()] - dmdts[phase.index()],
460  phase.fluidThermo().he()
461  );
462  }
463  }
464  }
465 
466  return eqnsPtr;
467 }
468 
469 
470 template<class BasePhaseSystem>
473 {
475  BasePhaseSystem::specieTransfer();
476 
477  phaseSystem::specieTransferTable& eqns = eqnsPtr();
478 
479  if (volatile_ != "none")
480  {
481  {
482  phaseSystem::dmidtfTable dmidtfs;
483 
484  forAllConstIter(phaseSystem::dmdtfTable, dmdtfs_, dmdtfIter)
485  {
486  const phaseInterface interface(*this, dmdtfIter.key());
487 
488  dmidtfs.insert(interface, new HashPtrTable<volScalarField>());
489  dmidtfs[interface]->insert
490  (
491  volatile_,
492  new volScalarField(*dmdtfIter())
493  );
494  }
495 
496  this->addDmidtYf(dmidtfs, eqns);
497  }
498 
499  {
500  phaseSystem::dmidtfTable nDmidtfs;
501 
502  forAllConstIter(phaseSystem::dmdtfTable, nDmdtfs_, nDmdtfIter)
503  {
504  const phaseInterface interface(*this, nDmdtfIter.key());
505 
506  nDmidtfs.insert(interface, new HashPtrTable<volScalarField>());
507  nDmidtfs[interface]->insert
508  (
509  volatile_,
510  new volScalarField(*nDmdtfIter())
511  );
512  }
513 
514  this->addDmidtYf(nDmidtfs, eqns);
515  }
516  }
517  else
518  {
519  this->addDmdtYfs(dmdtfs_, eqns);
520  this->addDmdtYfs(nDmdtfs_, eqns);
521  }
522 
523  return eqnsPtr;
524 }
525 
526 
527 template<class BasePhaseSystem>
530 {
531  dmdt0s_ = PtrList<volScalarField>(this->phases().size());
532 
533  addDmdts(dmdt0s_);
534 
535  BasePhaseSystem::correctContinuityError();
536 }
537 
538 
539 template<class BasePhaseSystem>
540 void
542 {
543  typedef
545  wallBoilingHeatTransferModel;
546 
548  wallBoilingHeatTransferModels =
549  this->mesh().template lookupClass<wallBoilingHeatTransferModel>();
550 
551  typedef
553  alphatPhaseChangeWallFunction;
554 
556  (
558  saturationModels_,
559  saturationModelIter
560  )
561  {
562  const phaseInterface& interface = saturationModelIter()->interface();
563  const phaseModel& phase1 = interface.phase1();
564  const phaseModel& phase2 = interface.phase2();
565  const rhoFluidThermo& thermo1 = phase1.fluidThermo();
566  const rhoFluidThermo& thermo2 = phase2.fluidThermo();
567  const volScalarField& T1(thermo1.T());
568  const volScalarField& T2(thermo2.T());
569 
571  this->heatTransferModels_[interface];
572 
573  // Interfacial mass transfer update
574  {
575  volScalarField& dmdtf(*this->dmdtfs_[interface]);
576  volScalarField& Tf(*this->Tfs_[interface]);
577 
578  volScalarField& Tsat(*this->Tsats_[interface]);
579  Tsat = saturationModelIter()->Tsat(thermo1.p());
580 
581  const volScalarField L
582  (
583  volatile_ != "none"
584  ? this->Li
585  (
586  interface,
587  volatile_,
588  dmdtf,
589  Tsat,
590  latentHeatScheme::symmetric
591  )
592  : this->L
593  (
594  interface,
595  dmdtf,
596  Tsat,
597  latentHeatScheme::symmetric
598  )
599  );
600 
601  volScalarField H1(heatTransferModel.modelInThe(phase1).K(0));
602  volScalarField H2(heatTransferModel.modelInThe(phase2).K(0));
603 
604  volScalarField dmdtfNew((H1*(Tsat - T1) + H2*(Tsat - T2))/L);
605 
606  if (volatile_ != "none")
607  {
608  dmdtfNew *=
609  neg0(dmdtfNew)*phase1.Y(volatile_)
610  + pos(dmdtfNew)*phase2.Y(volatile_);
611  }
612 
613  if (pressureImplicit_)
614  {
615  volScalarField& d2mdtdpf(*this->d2mdtdpfs_[interface]);
616 
617  const dimensionedScalar dp(rootSmall*thermo1.p().average());
618 
619  const volScalarField dTsatdp
620  (
621  (
622  saturationModelIter()->Tsat(thermo1.p() + dp/2)
623  - saturationModelIter()->Tsat(thermo1.p() - dp/2)
624  )/dp
625  );
626 
627  d2mdtdpf = (H1 + H2)*dTsatdp/L;
628 
629  if (volatile_ != "none")
630  {
631  d2mdtdpf *=
632  neg0(dmdtfNew)*phase1.Y(volatile_)
633  + pos(dmdtfNew)*phase2.Y(volatile_);
634  }
635  }
636 
637  H1 = heatTransferModel.modelInThe(phase1).K();
638  H2 = heatTransferModel.modelInThe(phase2).K();
639 
640  // Limit the H[12] to avoid /0
641  H1.max(small);
642  H2.max(small);
643 
644  Tf = (H1*T1 + H2*T2 + dmdtfNew*L)/(H1 + H2);
645 
646  Info<< Tsat.name()
647  << ": min = " << gMin(Tsat.primitiveField())
648  << ", mean = " << gAverage(Tsat.primitiveField())
649  << ", max = " << gMax(Tsat.primitiveField())
650  << endl;
651 
652  Info<< Tf.name()
653  << ": min = " << gMin(Tf.primitiveField())
654  << ", mean = " << gAverage(Tf.primitiveField())
655  << ", max = " << gMax(Tf.primitiveField())
656  << endl;
657 
658  const scalar dmdtfRelax =
659  this->mesh().solution().fieldRelaxationFactor(dmdtf.member());
660 
661  dmdtf = (1 - dmdtfRelax)*dmdtf + dmdtfRelax*dmdtfNew;
662 
663  Info<< dmdtf.name()
664  << ": min = " << gMin(dmdtf.primitiveField())
665  << ", mean = " << gAverage(dmdtf.primitiveField())
666  << ", max = " << gMax(dmdtf.primitiveField())
667  << ", integral = " << fvc::domainIntegrate(dmdtf).value()
668  << endl;
669  }
670 
671  // Nucleation mass transfer update
672  {
673  volScalarField& nDmdtf(*this->nDmdtfs_[interface]);
674  nDmdtf = Zero;
675 
676  bool wallBoilingActive = false;
677 
679  (
681  wallBoilingHeatTransferModels,
682  wallBoilingHeatTransferModelIter
683  )
684  {
685  const wallBoilingHeatTransferModel& wbht =
686  *wallBoilingHeatTransferModelIter();
687 
688  if (!wbht.activePhaseInterface(interface)) continue;
689 
690  wallBoilingActive = true;
691 
692  nDmdtf += (wbht.flipSign() ? -1 : +1)*wbht.dmdtf();
693  }
694 
695  forAllConstIter(phaseInterface, interface, interfaceIter)
696  {
697  const phaseModel& phase = interfaceIter();
698 
699  const word alphatName =
700  IOobject::groupName("alphat", phase.name());
701 
702  if (!phase.mesh().foundObject<volScalarField>(alphatName))
703  continue;
704 
705  const volScalarField& alphat =
706  phase.mesh().lookupObject<volScalarField>(alphatName);
707 
708  forAll(alphat.boundaryField(), patchi)
709  {
710  const fvPatchScalarField& alphatp =
711  alphat.boundaryField()[patchi];
712 
713  if (!isA<alphatPhaseChangeWallFunction>(alphatp)) continue;
714 
715  const alphatPhaseChangeWallFunction& alphatw =
716  refCast<const alphatPhaseChangeWallFunction>(alphatp);
717 
718  if (!alphatw.activeInterface(interface)) continue;
719 
720  wallBoilingActive = true;
721 
722  UIndirectList<scalar> nDmdtfp
723  (
724  nDmdtf.primitiveFieldRef(),
725  alphatp.patch().faceCells()
726  );
727 
728  nDmdtfp =
729  scalarField(nDmdtfp)
730  - (interfaceIter.index() == 0 ? +1 : -1)
731  *alphatw.dmdtf();
732  }
733  }
734 
735  if (wallBoilingActive)
736  {
737  Info<< nDmdtf.name()
738  << ": min = " << gMin(nDmdtf.primitiveField())
739  << ", mean = " << gAverage(nDmdtf.primitiveField())
740  << ", max = " << gMax(nDmdtf.primitiveField())
741  << ", integral = " << fvc::domainIntegrate(nDmdtf).value()
742  << endl;
743  }
744  }
745  }
746 }
747 
748 
749 template<class BasePhaseSystem>
751 {
752  if (BasePhaseSystem::read())
753  {
754  bool readOK = true;
755 
756  // Models ...
757 
758  return readOK;
759  }
760  else
761  {
762  return false;
763  }
764 }
765 
766 
767 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:477
const Mesh & mesh() const
Return mesh.
dimensioned< Type > average() const
Calculate and return arithmetic average.
Generic GeometricField class.
void max(const dimensioned< Type > &)
Internal::FieldType & primitiveFieldRef()
Return a reference to the primitive field.
const Internal::FieldType & primitiveField() const
Return a const-reference to the primitive field.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
bool insert(const Key &, const T &newElmt)
Insert a new hashedEntry.
Definition: HashTableI.H:80
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
static word member(const word &name)
Return member (name without the extension)
Definition: IOobject.C:149
const word & name() const
Return name.
Definition: IOobject.H:310
static word groupName(Name name, const word &group)
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
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
Definition: Switch.H:61
Class to provide interfacial heat and mass transfer between a number of phases according the interfac...
const interfaceSaturationTemperatureModel & saturation(const phaseInterfaceKey &key) const
Return the saturation temperature model for an interface.
virtual autoPtr< phaseSystem::momentumTransferTable > momentumTransfer()
Return the momentum transfer matrices for the cell-based algorithm.
virtual autoPtr< phaseSystem::heatTransferTable > heatTransfer() const
Return the heat transfer matrices.
virtual void correctInterfaceThermo()
Correct the interface thermodynamics.
virtual tmp< volScalarField > dmdtf(const phaseInterfaceKey &key) const
Return the mass transfer rate for an interface.
virtual autoPtr< phaseSystem::specieTransferTable > specieTransfer() const
Return the specie transfer matrices.
ThermalPhaseChangePhaseSystem(const fvMesh &)
Construct from fvMesh.
virtual PtrList< volScalarField > dmdts() const
Return the mass transfer rates for each phase.
virtual autoPtr< phaseSystem::momentumTransferTable > momentumTransferf()
Return the momentum transfer matrices for the face-based algorithm.
virtual PtrList< volScalarField > d2mdtdps() const
Return the mass transfer linearisation coeffs for each phase.
virtual bool read()
Read base phaseProperties dictionary.
virtual void correctContinuityError()
Store phase dmdts at the during the continuity error update.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
virtual const volScalarField & T() const =0
Temperature [K].
virtual const volScalarField & he() const =0
Enthalpy/Internal energy [J/kg].
Abstract base-class for all alphatWallFunctions supporting phase-change.
virtual const volScalarField & p() const =0
Pressure [Pa].
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:88
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:356
virtual const labelUList & faceCells() const
Return faceCells.
Definition: fvPatch.C:99
Model for heat transfer between phases.
tmp< volScalarField > K() const
The heat transfer function K used in the enthalpy equation.
Wrapper around saturationTemperatureModel to facilitate convenient construction on interfaces.
tmp< volScalarField > Tsat(const volScalarField &p) const
Saturation temperature.
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 ...
const phaseModel & phase1() const
Return phase 1.
const phaseModel & phase2() const
Return phase 2.
virtual const PtrList< volScalarField > & Y() const =0
Return the species mass fractions.
virtual const rhoFluidThermo & fluidThermo() const =0
Return the thermophysical model.
label index() const
Return the index of the phase.
Definition: phaseModel.C:157
const word & name() const
Return the name of this phase.
Definition: phaseModel.C:145
Base-class for fluid thermodynamic properties based on density.
A class for managing temporary objects.
Definition: tmp.H:55
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
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
Volume integrate volField creating a volField.
Calculate the matrix for implicit and explicit sources.
label patchi
void read(Istream &, label &, const dictionary &)
In-place read with dictionary lookup.
dimensioned< Type > domainIntegrate(const VolField< Type > &vf)
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const VolField< Type > &)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
static const zero Zero
Definition: zero.H:97
dimensionedScalar pos(const dimensionedScalar &ds)
const dimensionSet dimPressure
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
void addField(const label phasei, const word &name, tmp< GeoField > field, PtrList< GeoField > &fieldList)
Definition: phaseSystemI.H:263
messageStream Info
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
const dimensionSet dimTime
const dimensionSet dimDensity
VolField< scalar > volScalarField
Definition: volFieldsFwd.H:64
Type gAverage(const FieldField< Field, Type > &f)
error FatalError
dimensionedScalar neg0(const dimensionedScalar &ds)
Type gMin(const FieldField< Field, Type > &f)
Type gMax(const FieldField< Field, Type > &f)