All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
InterfaceCompositionPhaseChangePhaseSystem.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-2021 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 
28 #include "heatTransferModel.H"
30 #include "fvmSup.H"
31 
32 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
33 
34 template<class BasePhaseSystem>
37 {
39  (
40  interfaceCompositionModelTable,
41  interfaceCompositionModels_,
42  interfaceCompositionModelIter
43  )
44  {
45  const phasePair& pair =
46  this->phasePairs_[interfaceCompositionModelIter.key()];
47 
48  *dmdtfs_[pair] = Zero;
49 
50  forAllConstIter(phasePair, pair, pairIter)
51  {
52  const autoPtr<interfaceCompositionModel>& compositionModelPtr =
53  interfaceCompositionModelIter()[pairIter.index()];
54 
55  if (!compositionModelPtr.valid()) continue;
56 
57  const interfaceCompositionModel& compositionModel =
58  compositionModelPtr();
59 
60  const phaseModel& phase = *pairIter;
61 
63  (
64  hashedWordList,
65  compositionModel.species(),
66  specieIter
67  )
68  {
69  const word& specie = *specieIter;
70 
71  *dmdtfs_[pair] +=
72  (pairIter.index() == 0 ? +1 : -1)
73  *(
74  *(*dmidtfSus_[pair])[specie]
75  + *(*dmidtfSps_[pair])[specie]*phase.Y(specie)
76  );
77  }
78  }
79  }
80 }
81 
82 
83 template<class BasePhaseSystem>
86 totalDmidtfs() const
87 {
88  autoPtr<phaseSystem::dmidtfTable> totalDmidtfsPtr
89  (
90  new phaseSystem::dmidtfTable
91  );
92  phaseSystem::dmidtfTable& totalDmidtfs = totalDmidtfsPtr();
93 
95  (
96  interfaceCompositionModelTable,
97  interfaceCompositionModels_,
98  interfaceCompositionModelIter
99  )
100  {
101  const phasePair& pair =
102  this->phasePairs_[interfaceCompositionModelIter.key()];
103 
104  if (!totalDmidtfs.found(pair))
105  {
106  totalDmidtfs.insert(pair, new HashPtrTable<volScalarField>());
107  }
108 
109  forAllConstIter(phasePair, pair, pairIter)
110  {
111  const autoPtr<interfaceCompositionModel>& compositionModelPtr =
112  interfaceCompositionModelIter()[pairIter.index()];
113 
114  if (!compositionModelPtr.valid()) continue;
115 
116  const interfaceCompositionModel& compositionModel =
117  compositionModelPtr();
118 
119  const phaseModel& phase = *pairIter;
120 
122  (
123  hashedWordList,
124  compositionModel.species(),
125  specieIter
126  )
127  {
128  const word& specie = *specieIter;
129 
130 
131  tmp<volScalarField> dmidtf
132  (
133  (pairIter.index() == 0 ? +1 : -1)
134  *(
135  *(*dmidtfSus_[pair])[specie]
136  + *(*dmidtfSps_[pair])[specie]*phase.Y(specie)
137  )
138  );
139 
140  if (totalDmidtfs[pair]->found(specie))
141  {
142  *(*totalDmidtfs[pair])[specie] += dmidtf;
143  }
144  else
145  {
146  totalDmidtfs[pair]->insert(specie, dmidtf.ptr());
147  }
148  }
149  }
150  }
151 
152  return totalDmidtfsPtr;
153 }
154 
155 
156 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
157 
158 template<class BasePhaseSystem>
161 (
162  const fvMesh& mesh
163 )
164 :
165  BasePhaseSystem(mesh),
166  nInterfaceCorrectors_
167  (
168  this->template lookupOrDefault<label>("nInterfaceCorrectors", 1)
169  )
170 {
171  this->generatePairsAndSubModels
172  (
173  "interfaceComposition",
174  interfaceCompositionModels_
175  );
176 
177  this->generatePairsAndSubModels
178  (
179  "diffusiveMassTransfer",
180  diffusiveMassTransferModels_,
181  false
182  );
183 
184  // Check that models have been specified in the correct combinations
186  (
187  interfaceCompositionModelTable,
188  interfaceCompositionModels_,
189  interfaceCompositionModelIter
190  )
191  {
192  const phasePair& pair =
193  this->phasePairs_[interfaceCompositionModelIter.key()];
194 
195  this->template validateMassTransfer<interfaceCompositionModel>(pair);
196 
197  if (!this->diffusiveMassTransferModels_.found(pair))
198  {
200  << "A diffusive mass transfer model the " << pair
201  << " pair is not specified. This is required by the "
202  << "corresponding interface composition model."
203  << exit(FatalError);
204  }
205 
206  forAllConstIter(phasePair, pair, pairIter)
207  {
208  if
209  (
210  interfaceCompositionModelIter()[pairIter.index()].valid()
211  && !diffusiveMassTransferModels_[pair][pairIter.index()].valid()
212  )
213  {
215  << "A mass transfer model for the " << (*pairIter).name()
216  << " side of the " << pair << " pair is not "
217  << "specified. This is required by the corresponding "
218  << "interface composition model."
219  << exit(FatalError);
220  }
221  }
222 
223  if
224  (
225  !this->heatTransferModels_.found(pair)
226  || !this->heatTransferModels_[pair].first().valid()
227  || !this->heatTransferModels_[pair].second().valid()
228  )
229  {
231  << "A heat transfer model for both sides of the " << pair
232  << "pair is not specified. This is required by the "
233  << "corresponding interface composition model"
234  << exit(FatalError);
235  }
236  }
237 
238  // Generate mass transfer fields, initially assumed to be zero
240  (
241  interfaceCompositionModelTable,
242  interfaceCompositionModels_,
243  interfaceCompositionModelIter
244  )
245  {
246  const phasePair& pair =
247  this->phasePairs_[interfaceCompositionModelIter.key()];
248 
249  dmdtfs_.insert
250  (
251  pair,
252  new volScalarField
253  (
254  IOobject
255  (
257  (
258  "interfaceCompositionPhaseChange:dmdtf",
259  pair.name()
260  ),
261  this->mesh().time().timeName(),
262  this->mesh(),
265  ),
266  this->mesh(),
268  )
269  );
270 
271  dmidtfSus_.insert(pair, new HashPtrTable<volScalarField>());
272 
273  dmidtfSps_.insert(pair, new HashPtrTable<volScalarField>());
274 
275  Tfs_.insert
276  (
277  pair,
278  new volScalarField
279  (
280  IOobject
281  (
283  (
284  "interfaceCompositionPhaseChange:Tf",
285  pair.name()
286  ),
287  this->mesh().time().timeName(),
288  this->mesh(),
291  ),
292  (pair.phase1().thermo().T() + pair.phase2().thermo().T())/2
293  )
294  );
295 
296  forAllConstIter(phasePair, pair, pairIter)
297  {
298  const autoPtr<interfaceCompositionModel>& compositionModelPtr =
299  interfaceCompositionModelIter()[pairIter.index()];
300 
301  if (!compositionModelPtr.valid()) continue;
302 
303  const interfaceCompositionModel& compositionModel =
304  compositionModelPtr();
305 
307  (
308  hashedWordList,
309  compositionModel.species(),
310  specieIter
311  )
312  {
313  const word& specie = *specieIter;
314 
315  dmidtfSus_[pair]->insert
316  (
317  specie,
318  new volScalarField
319  (
320  IOobject
321  (
323  (
325  (
326  "interfaceCompositionPhaseChange:dmidtfSu",
327  specie
328  ),
329  pair.name()
330  ),
331  this->mesh().time().timeName(),
332  this->mesh()
333  ),
334  this->mesh(),
336  )
337  );
338 
339  dmidtfSps_[pair]->insert
340  (
341  specie,
342  new volScalarField
343  (
344  IOobject
345  (
347  (
349  (
350  "interfaceCompositionPhaseChange:dmidtfSp",
351  specie
352  ),
353  pair.name()
354  ),
355  this->mesh().time().timeName(),
356  this->mesh()
357  ),
358  this->mesh(),
360  )
361  );
362  }
363  }
364  }
365 }
366 
367 
368 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
369 
370 template<class BasePhaseSystem>
373 {}
374 
375 
376 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
377 
378 template<class BasePhaseSystem>
381 (
382  const phasePairKey& key
383 ) const
384 {
385  tmp<volScalarField> tDmdtf = BasePhaseSystem::dmdtf(key);
386 
387  if (dmdtfs_.found(key))
388  {
389  const label dmdtSign(Pair<word>::compare(this->phasePairs_[key], key));
390 
391  tDmdtf.ref() += dmdtSign**dmdtfs_[key];
392  }
393 
394  return tDmdtf;
395 }
396 
397 
398 template<class BasePhaseSystem>
401 {
402  PtrList<volScalarField> dmdts(BasePhaseSystem::dmdts());
403 
404  forAllConstIter(phaseSystem::dmdtfTable, dmdtfs_, dmdtfsIter)
405  {
406  const phasePair& pair = this->phasePairs_[dmdtfsIter.key()];
407  const phaseModel& phase = pair.phase1();
408  const phaseModel& otherPhase = pair.phase2();
409 
410  addField(phase, "dmdt", *dmdtfsIter(), dmdts);
411  addField(otherPhase, "dmdt", - *dmdtfsIter(), dmdts);
412  }
413 
414  return dmdts;
415 }
416 
417 
418 template<class BasePhaseSystem>
422 {
423  autoPtr<phaseSystem::momentumTransferTable> eqnsPtr =
425 
426  phaseSystem::momentumTransferTable& eqns = eqnsPtr();
427 
428  this->addDmdtUfs(dmdtfs_, eqns);
429 
430  return eqnsPtr;
431 }
432 
433 
434 template<class BasePhaseSystem>
438 {
439  autoPtr<phaseSystem::momentumTransferTable> eqnsPtr =
440  BasePhaseSystem::momentumTransferf();
441 
442  phaseSystem::momentumTransferTable& eqns = eqnsPtr();
443 
444  this->addDmdtUfs(dmdtfs_, eqns);
445 
446  return eqnsPtr;
447 }
448 
449 
450 template<class BasePhaseSystem>
453 heatTransfer() const
454 {
455  autoPtr<phaseSystem::heatTransferTable> eqnsPtr =
456  BasePhaseSystem::heatTransfer();
457 
458  phaseSystem::heatTransferTable& eqns = eqnsPtr();
459 
460  this->addDmidtHefs
461  (
462  totalDmidtfs(),
463  Tfs_,
465  latentHeatTransfer::mass,
466  eqns
467  );
468 
469  return eqnsPtr;
470 }
471 
472 
473 template<class BasePhaseSystem>
476 specieTransfer() const
477 {
478  autoPtr<phaseSystem::specieTransferTable> eqnsPtr =
480 
481  phaseSystem::specieTransferTable& eqns = eqnsPtr();
482 
483  // Explicit
484  /*
485  this->addDmidtYf(totalDmidtfs(), eqns);
486  */
487 
488  // Semi-implicit
490  (
491  interfaceCompositionModelTable,
492  interfaceCompositionModels_,
493  interfaceCompositionModelIter
494  )
495  {
496  const phasePair& pair =
497  this->phasePairs_[interfaceCompositionModelIter.key()];
498 
499  forAllConstIter(phasePair, pair, pairIter)
500  {
501  const autoPtr<interfaceCompositionModel>& compositionModelPtr =
502  interfaceCompositionModelIter()[pairIter.index()];
503 
504  if (!compositionModelPtr.valid()) continue;
505 
506  const interfaceCompositionModel& compositionModel =
507  compositionModelPtr();
508 
509  const phaseModel& phase = *pairIter;
510  const phaseModel& otherPhase = pairIter.otherPhase();
511 
513  (
514  hashedWordList,
515  compositionModel.species(),
516  specieIter
517  )
518  {
519  const word& specie = *specieIter;
520 
521  // Implicit transport through this phase
522  *eqns[phase.Y(specie).name()] +=
523  *(*dmidtfSus_[pair])[specie]
524  + fvm::Sp(*(*dmidtfSps_[pair])[specie], phase.Y(specie));
525 
526  // Explicit transport out of the other phase
527  if (eqns.found(IOobject::groupName(specie, otherPhase.name())))
528  {
529  *eqns[otherPhase.Y(specie).name()] -=
530  *(*dmidtfSus_[pair])[specie]
531  + *(*dmidtfSps_[pair])[specie]*phase.Y(specie);
532  }
533  }
534  }
535  }
536 
537  return eqnsPtr;
538 }
539 
540 
541 template<class BasePhaseSystem>
543 correct()
544 {
546 
547  // Sum up the contribution from each interface composition model
549  (
550  interfaceCompositionModelTable,
551  interfaceCompositionModels_,
552  interfaceCompositionModelIter
553  )
554  {
555  const phasePair& pair =
556  this->phasePairs_[interfaceCompositionModelIter.key()];
557 
558  const volScalarField& Tf(*this->Tfs_[pair]);
559 
560  forAllConstIter(phasePair, pair, pairIter)
561  {
562  const autoPtr<interfaceCompositionModel>& compositionModelPtr =
563  interfaceCompositionModelIter()[pairIter.index()];
564 
565  if (!compositionModelPtr.valid()) continue;
566 
567  const interfaceCompositionModel& compositionModel =
568  compositionModelPtr();
569 
570  const phaseModel& phase = *pairIter;
571 
572  const volScalarField K
573  (
574  diffusiveMassTransferModels_[pair][pairIter.index()]->K()
575  );
576 
578  (
579  hashedWordList,
580  compositionModel.species(),
581  specieIter
582  )
583  {
584  const word& specie = *specieIter;
585 
586  const volScalarField KD(K*compositionModel.D(specie));
587  const volScalarField Yf(compositionModel.Yf(specie, Tf));
588 
589  *(*dmidtfSus_[pair])[specie] = phase.rho()*KD*Yf;
590  *(*dmidtfSps_[pair])[specie] = - phase.rho()*KD;
591  }
592  }
593  }
594 
595  correctDmdtfs();
596 }
597 
598 
599 template<class BasePhaseSystem>
602 {
604 
605  correctDmdtfs();
606 }
607 
608 
609 template<class BasePhaseSystem>
612 {
613  // This loop solves for the interface temperatures, Tf, and updates the
614  // interface composition models.
615  //
616  // In the presence of thermally-coupled mass transfer, the relation between
617  // heat transfers across the interface between phases 1 and 2 is:
618  //
619  // Q1 + Q2 = mDot*L
620  // H1*(Tf - T1) + H2*(Tf - T1) = K*rho*(Yfi - Yi)*Li
621  //
622  // Where Q1 and Q2 are the net transfer into phases 1 and 2 respectively,
623  // H1 and H2 are the heat transfer coefficients on either side, Tf is the
624  // temperature at the interface, mDot is the mass transfer rate from phase
625  // 2 to phase 1, and L is the latent heat of phase 2 minus phase 1, K is
626  // the diffusive mass transfer coefficient, Yfi - Yi is the concentration
627  // difference of a transferring specie between the interface and the bulk
628  // driving the transfer, Li is the latent heat change of the specie, and
629  // rho is the density in the phase in which the diffusive mass transfer is
630  // being represented.
631  //
632  // Yfi is likely to be a strong non-linear (typically exponential) function
633  // of Tf, so the solution for the temperature is newton-accelerated.
634 
636  (
637  interfaceCompositionModelTable,
638  interfaceCompositionModels_,
639  interfaceCompositionModelIter
640  )
641  {
642  const phasePair& pair =
643  this->phasePairs_[interfaceCompositionModelIter.key()];
644 
645  const volScalarField H1(this->heatTransferModels_[pair].first()->K());
646  const volScalarField H2(this->heatTransferModels_[pair].second()->K());
647  const dimensionedScalar HSmall("small", heatTransferModel::dimK, small);
648 
649  const typename diffusiveMassTransferModelTable::value_type&
650  diffusiveMassTransfers = this->diffusiveMassTransferModels_[pair];
651 
652  const typename interfaceCompositionModelTable::value_type &
653  interfaceCompositions = this->interfaceCompositionModels_[pair];
654 
655  volScalarField& Tf = *this->Tfs_[pair];
656 
657  for (label i = 0; i < nInterfaceCorrectors_; ++ i)
658  {
659  tmp<volScalarField> dmdtLf =
661  (
662  IOobject::groupName("dmdtLf", pair.name()),
663  this->mesh(),
665  );
666  tmp<volScalarField> dmdtLfPrime =
668  (
669  IOobject::groupName("dmdtLfPrime", pair.name()),
670  this->mesh(),
671  dimensionedScalar(dmdtLf().dimensions()/dimTemperature, 0)
672  );
673 
674  // Add latent heats from forward and backward models
675  forAllConstIter(phasePair, pair, pairIter)
676  {
677  if (interfaceCompositions[pairIter.index()].valid())
678  {
679  const BlendedInterfacialModel<diffusiveMassTransferModel>&
680  diffusiveMassTransfer =
681  diffusiveMassTransfers[pairIter.index()];
682 
683  const interfaceCompositionModel&
684  interfaceComposition =
685  interfaceCompositions[pairIter.index()];
686 
687  const label sign = pairIter.index() == 0 ? 1 : -1;
688 
690  (
691  hashedWordList,
692  interfaceComposition.species(),
693  specieIter
694  )
695  {
696  const word& specie = *specieIter;
697 
698  const volScalarField dY
699  (
700  interfaceComposition.dY(specie, Tf)
701  );
702 
703  const volScalarField dYfPrime
704  (
705  interfaceComposition.dYfPrime(specie, Tf)
706  );
707 
708  const volScalarField rhoKDL
709  (
710  pairIter().thermo().rho()
711  *diffusiveMassTransfer.K()
712  *interfaceComposition.D(specie)
713  *this->Li
714  (
715  pair,
716  specie,
717  dY,
718  Tf,
720  )
721  );
722 
723  dmdtLf.ref() += sign*rhoKDL*dY;
724  dmdtLfPrime.ref() += sign*rhoKDL*dYfPrime;
725  }
726  }
727  }
728 
729  // Update the interface temperature by applying one step of newton's
730  // method to the interface relation
731  Tf -=
732  (
733  H1*(Tf - pair.phase1().thermo().T())
734  + H2*(Tf - pair.phase2().thermo().T())
735  - dmdtLf
736  )
737  /(
738  max(H1 + H2 - dmdtLfPrime, HSmall)
739  );
740 
741  Tf.correctBoundaryConditions();
742 
743  Info<< "Tf." << pair.name()
744  << ": min = " << min(Tf.primitiveField())
745  << ", mean = " << average(Tf.primitiveField())
746  << ", max = " << max(Tf.primitiveField())
747  << endl;
748 
749  // Update the interface compositions
750  if (this->interfaceCompositionModels_[pair].first().valid())
751  {
752  this->interfaceCompositionModels_[pair].first()->update(Tf);
753  }
754  if (this->interfaceCompositionModels_[pair].second().valid())
755  {
756  this->interfaceCompositionModels_[pair].second()->update(Tf);
757  }
758  }
759  }
760 }
761 
762 
763 template<class BasePhaseSystem>
765 {
766  if (BasePhaseSystem::read())
767  {
768  bool readOK = true;
769 
770  // Models ...
771 
772  return readOK;
773  }
774  else
775  {
776  return false;
777  }
778 }
779 
780 
781 // ************************************************************************* //
dimensionedScalar sign(const dimensionedScalar &ds)
virtual autoPtr< phaseSystem::specieTransferTable > specieTransfer() const
Return the specie transfer matrices.
virtual ~InterfaceCompositionPhaseChangePhaseSystem()
Destructor.
virtual PtrList< volScalarField > dmdts() const
Return the mass transfer rates for each phase.
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
fluidReactionThermo & thermo
Definition: createFields.H:28
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
HashPtrTable< fvVectorMatrix > momentumTransferTable
Definition: phaseSystem.H:78
InterfaceCompositionPhaseChangePhaseSystem(const fvMesh &)
Construct from fvMesh.
static int compare(const Pair< word > &a, const Pair< word > &b)
Compare Pairs.
Definition: Pair.H:153
static tmp< GeometricField< scalar, fvPatchField, volMesh > > New(const word &name, const Internal &, const PtrList< fvPatchField< scalar >> &)
Return a temporary field constructed from name,.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
HashPtrTable< fvScalarMatrix > heatTransferTable
Definition: phaseSystem.H:80
HashPtrTable< fvScalarMatrix > specieTransferTable
Definition: phaseSystem.H:82
void read(Istream &, label &, const dictionary &)
In-place read with dictionary lookup.
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
CGAL::Exact_predicates_exact_constructions_kernel K
virtual void correctInterfaceThermo()
Correct the interface temperatures.
virtual tmp< volScalarField > dmdtf(const phasePairKey &key) const
Return the mass transfer rate for an interface.
fluid correctSpecies()
static const dimensionSet dimK
Coefficient dimensions.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
const dimensionSet dimTime
dynamicFvMesh & mesh
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(rho0+psi *p, rhoMin);# 1 "/home/ubuntu/OpenFOAM-9/applications/solvers/multiphase/cavitatingFoam/alphavPsi.H" 1{ alphav=max(min((rho - rholSat)/(rhovSat - rholSat), scalar(1)), scalar(0));alphal=1.0 - alphav;Info<< "max-min alphav: "<< max(alphav).value()<< " "<< min(alphav).value()<< endl;psiModel-> correct()
Definition: pEqn.H:68
static word groupName(Name name, const word &group)
HashPtrTable< volScalarField, phasePairKey, phasePairKey::hash > dmdtfTable
Definition: phaseSystem.H:94
const dimensionSet dimDensity
static const zero Zero
Definition: zero.H:97
dimensioned< Type > average(const DimensionedField< Type, GeoMesh > &df)
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
virtual void correct()
Correct the fluid properties other than those listed below.
phaseSystem::specieTransferTable & specieTransfer(specieTransferPtr())
virtual tmp< volScalarField > Li(const phasePair &pair, const word &member, const volScalarField &dmidtf, const volScalarField &Tf, const latentHeatScheme scheme) const =0
Return the latent heat for a given pair, specie, mass transfer rate.
virtual autoPtr< phaseSystem::momentumTransferTable > momentumTransfer()
Return the momentum transfer matrices for the cell-based algorithm.
const dimensionSet dimEnergy
virtual void correctSpecies()
Correct the mass transfer rates for the new species mass fractions.
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
phaseSystem::momentumTransferTable & momentumTransfer(momentumTransferPtr())
Internal & ref()
Return a reference to the dimensioned internal field.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
const dimensionSet dimVolume
messageStream Info
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
A class for managing temporary objects.
Definition: PtrList.H:53
virtual autoPtr< phaseSystem::momentumTransferTable > momentumTransferf()
Return the momentum transfer matrices for the face-based algorithm.
virtual autoPtr< phaseSystem::heatTransferTable > heatTransfer() const
Return the heat transfer matrices.
bool found
virtual bool read()
Read base phaseProperties dictionary.
Calculate the matrix for implicit and explicit sources.
const dimensionSet dimTemperature
Class to provide interfacial heat and mass transfer between a number of phases according to a interfa...