thermoSingleLayer.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2016 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 "thermoSingleLayer.H"
27 #include "fvcDiv.H"
28 #include "fvcLaplacian.H"
29 #include "fvcFlux.H"
30 #include "fvm.H"
34 #include "mapDistribute.H"
35 #include "constants.H"
36 
37 // Sub-models
38 #include "filmThermoModel.H"
39 #include "filmViscosityModel.H"
40 #include "heatTransferModel.H"
41 #include "phaseChangeModel.H"
42 #include "filmRadiationModel.H"
43 
44 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
45 
46 namespace Foam
47 {
48 namespace regionModels
49 {
50 namespace surfaceFilmModels
51 {
52 
53 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
54 
55 defineTypeNameAndDebug(thermoSingleLayer, 0);
56 
57 addToRunTimeSelectionTable(surfaceFilmModel, thermoSingleLayer, mesh);
58 
59 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
60 
61 wordList thermoSingleLayer::hsBoundaryTypes()
62 {
63  wordList bTypes(T_.boundaryField().types());
64  forAll(bTypes, patchi)
65  {
66  if
67  (
68  T_.boundaryField()[patchi].fixesValue()
70  )
71  {
73  }
74  }
75 
76  return bTypes;
77 }
78 
79 
80 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
81 
83 {
84  // No additional properties to read
86 }
87 
88 
90 {
91  if (debug)
92  {
94  }
95 
97 
99 }
100 
101 
103 {
104  rho_ == filmThermo_->rho();
105  sigma_ == filmThermo_->sigma();
106  Cp_ == filmThermo_->Cp();
107  kappa_ == filmThermo_->kappa();
108 }
109 
110 
112 {
114 
115  volScalarField::Boundary& hsBf = hs_.boundaryFieldRef();
116 
117  forAll(hsBf, patchi)
118  {
121  {
122  hsBf[patchi] == hs(Tp, patchi);
123  }
124  }
125 }
126 
127 
129 {
131 
132  // Push boundary film temperature into wall temperature internal field
133  for (label i=0; i<intCoupledPatchIDs_.size(); i++)
134  {
136  const polyPatch& pp = regionMesh().boundaryMesh()[patchi];
139  }
141 
142  // Update film surface temperature
143  Ts_ = T_;
145 }
146 
147 
149 {
150  if (debug)
151  {
152  InfoInFunction << endl;
153  }
154 
156 
157  // Update primary region fields on local region via direct mapped (coupled)
158  // boundary conditions
160  forAll(YPrimary_, i)
161  {
162  YPrimary_[i].correctBoundaryConditions();
163  }
164 }
165 
166 
168 {
169  if (debug)
170  {
171  InfoInFunction << endl;
172  }
173 
175 
176  volScalarField::Boundary& hsSpPrimaryBf =
178 
179  // Convert accummulated source terms into per unit area per unit time
180  const scalar deltaT = time_.deltaTValue();
181  forAll(hsSpPrimaryBf, patchi)
182  {
183  const scalarField& priMagSf =
185 
186  hsSpPrimaryBf[patchi] /= priMagSf*deltaT;
187  }
188 
189  // Retrieve the source fields from the primary region via direct mapped
190  // (coupled) boundary conditions
191  // - fields require transfer of values for both patch AND to push the
192  // values into the first layer of internal cells
194 
195  // Apply enthalpy source as difference between incoming and actual states
196  hsSp_ -= rhoSp_*hs_;
197 }
198 
199 
201 {
202  if (hydrophilic_)
203  {
204  const scalar hydrophilicDry = hydrophilicDryScale_*deltaWet_;
205  const scalar hydrophilicWet = hydrophilicWetScale_*deltaWet_;
206 
207  forAll(alpha_, i)
208  {
209  if ((alpha_[i] < 0.5) && (delta_[i] > hydrophilicWet))
210  {
211  alpha_[i] = 1.0;
212  }
213  else if ((alpha_[i] > 0.5) && (delta_[i] < hydrophilicDry))
214  {
215  alpha_[i] = 0.0;
216  }
217  }
218 
220  }
221  else
222  {
223  alpha_ ==
225  }
226 }
227 
228 
230 {
231  if (debug)
232  {
233  InfoInFunction << endl;
234  }
235 
236  // Update heat transfer coefficient sub-models
237  htcs_->correct();
238  htcw_->correct();
239 
240  phaseChange_->correct
241  (
242  time_.deltaTValue(),
246  );
247 
248  // Update radiation
249  radiation_->correct();
250 
251  // Update kinematic sub-models
253 
254  // Update source fields
256  rhoSp_ += primaryMassPCTrans_/magSf()/time().deltaT();
257 
258  // Vapour recoil pressure
259  pSp_ -= sqr(primaryMassPCTrans_/magSf()/time().deltaT())/2.0/rhoPrimary_;
260 }
261 
262 
264 {
265  return
266  (
267  // Heat-transfer to the primary region
268  - fvm::Sp(htcs_->h()/Cp_, hs)
269  + htcs_->h()*(hs/Cp_ + alpha_*(TPrimary_ - T_))
270 
271  // Heat-transfer to the wall
272  - fvm::Sp(htcw_->h()/Cp_, hs)
273  + htcw_->h()*(hs/Cp_ + alpha_*(Tw_- T_))
274  );
275 }
276 
277 
279 {
280  if (debug)
281  {
282  InfoInFunction << endl;
283  }
284 
286 
287  solve
288  (
290  + fvm::div(phi_, hs_)
291  ==
292  - hsSp_
293  - rhoSp_*hs_
294  + q(hs_)
295  + radiation_->Shs()
296  );
297 
299 
300  // Evaluate viscosity from user-model
301  viscosity_->correct(pPrimary_, T_);
302 }
303 
304 
305 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
306 
307 thermoSingleLayer::thermoSingleLayer
308 (
309  const word& modelType,
310  const fvMesh& mesh,
311  const dimensionedVector& g,
312  const word& regionType,
313  const bool readFields
314 )
315 :
316  kinematicSingleLayer(modelType, mesh, g, regionType, false),
317  thermo_(mesh.lookupObject<SLGThermo>("SLGThermo")),
318  Cp_
319  (
320  IOobject
321  (
322  "Cp",
323  time().timeName(),
324  regionMesh(),
327  ),
328  regionMesh(),
330  zeroGradientFvPatchScalarField::typeName
331  ),
332  kappa_
333  (
334  IOobject
335  (
336  "kappa",
337  time().timeName(),
338  regionMesh(),
341  ),
342  regionMesh(),
344  (
345  "kappa",
347  0.0
348  ),
349  zeroGradientFvPatchScalarField::typeName
350  ),
351 
352  T_
353  (
354  IOobject
355  (
356  "Tf",
357  time().timeName(),
358  regionMesh(),
361  ),
362  regionMesh()
363  ),
364  Ts_
365  (
366  IOobject
367  (
368  "Tsf",
369  time().timeName(),
370  regionMesh(),
373  ),
374  T_,
375  zeroGradientFvPatchScalarField::typeName
376  ),
377  Tw_
378  (
379  IOobject
380  (
381  "Twf",
382  time().timeName(),
383  regionMesh(),
386  ),
387  T_,
388  zeroGradientFvPatchScalarField::typeName
389  ),
390  hs_
391  (
392  IOobject
393  (
394  "hf",
395  time().timeName(),
396  regionMesh(),
399  ),
400  regionMesh(),
401  dimensionedScalar("zero", dimEnergy/dimMass, 0.0),
402  hsBoundaryTypes()
403  ),
404 
406  (
407  IOobject
408  (
409  "primaryMassPCTrans",
410  time().timeName(),
411  regionMesh(),
414  ),
415  regionMesh(),
416  dimensionedScalar("zero", dimMass, 0),
417  zeroGradientFvPatchScalarField::typeName
418  ),
420  (
421  IOobject
422  (
423  "primaryEnergyPCTrans",
424  time().timeName(),
425  regionMesh(),
428  ),
429  regionMesh(),
430  dimensionedScalar("zero", dimEnergy, 0),
431  zeroGradientFvPatchScalarField::typeName
432  ),
433 
434  deltaWet_(readScalar(coeffs_.lookup("deltaWet"))),
435  hydrophilic_(readBool(coeffs_.lookup("hydrophilic"))),
438 
439  hsSp_
440  (
441  IOobject
442  (
443  "hsSp",
444  time().timeName(),
445  regionMesh(),
448  ),
449  regionMesh(),
451  this->mappedPushedFieldPatchTypes<scalar>()
452  ),
453 
455  (
456  IOobject
457  (
458  hsSp_.name(), // Must have same name as hSp_ to enable mapping
459  time().timeName(),
460  primaryMesh(),
463  ),
464  primaryMesh(),
465  dimensionedScalar("zero", hsSp_.dimensions(), 0.0)
466  ),
467 
468  TPrimary_
469  (
470  IOobject
471  (
472  "T", // Same name as T on primary region to enable mapping
473  time().timeName(),
474  regionMesh(),
477  ),
478  regionMesh(),
479  dimensionedScalar("zero", dimTemperature, 0.0),
480  this->mappedFieldAndInternalPatchTypes<scalar>()
481  ),
482 
483  YPrimary_(),
484 
486  htcs_
487  (
488  heatTransferModel::New(*this, coeffs().subDict("upperSurfaceModels"))
489  ),
490  htcw_
491  (
492  heatTransferModel::New(*this, coeffs().subDict("lowerSurfaceModels"))
493  ),
496  Tmin_(-VGREAT),
497  Tmax_(VGREAT)
498 {
499  if (coeffs().readIfPresent("Tmin", Tmin_))
500  {
501  Info<< " limiting minimum temperature to " << Tmin_ << endl;
502  }
503 
504  if (coeffs().readIfPresent("Tmax", Tmax_))
505  {
506  Info<< " limiting maximum temperature to " << Tmax_ << endl;
507  }
508 
510  {
511  YPrimary_.setSize(thermo_.carrier().species().size());
512 
513  forAll(thermo_.carrier().species(), i)
514  {
515  YPrimary_.set
516  (
517  i,
518  new volScalarField
519  (
520  IOobject
521  (
522  thermo_.carrier().species()[i],
523  time().timeName(),
524  regionMesh(),
527  ),
528  regionMesh(),
529  dimensionedScalar("zero", dimless, 0.0),
530  pSp_.boundaryField().types()
531  )
532  );
533  }
534  }
535 
536  if (hydrophilic_)
537  {
538  coeffs_.lookup("hydrophilicDryScale") >> hydrophilicDryScale_;
539  coeffs_.lookup("hydrophilicWetScale") >> hydrophilicWetScale_;
540  }
541 
542  if (readFields)
543  {
545 
546  correctAlpha();
547 
549 
550  // Update derived fields
551  hs_ == hs(T_);
552 
553  deltaRho_ == delta_*rho_;
554 
556  (
557  IOobject
558  (
559  "phi",
560  time().timeName(),
561  regionMesh(),
564  false
565  ),
567  );
568 
569  phi_ == phi0;
570 
571  // Evaluate viscosity from user-model
572  viscosity_->correct(pPrimary_, T_);
573  }
574 }
575 
576 
577 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
578 
580 {}
581 
582 
583 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
584 
586 (
587  const label patchi,
588  const label facei,
589  const scalar massSource,
590  const vector& momentumSource,
591  const scalar pressureSource,
592  const scalar energySource
593 )
594 {
596  (
597  patchi,
598  facei,
599  massSource,
600  momentumSource,
601  pressureSource,
602  energySource
603  );
604 
605  if (debug)
606  {
607  Info<< " energy = " << energySource << nl << endl;
608  }
609 
610  hsSpPrimary_.boundaryFieldRef()[patchi][facei] -= energySource;
611 }
612 
613 
615 {
616  if (debug)
617  {
618  InfoInFunction << endl;
619  }
620 
622 
623  // Update phase change
626 }
627 
628 
630 {
631  if (debug)
632  {
633  InfoInFunction << endl;
634  }
635 
636  // Update film coverage indicator
637  correctAlpha();
638 
639  // Update film wall and surface velocities
641 
642  // Update film wall and surface temperatures
644 
645  // Update sub-models to provide updated source contributions
646  updateSubmodels();
647 
648  // Solve continuity for deltaRho_
649  solveContinuity();
650 
651  for (int oCorr=1; oCorr<=nOuterCorr_; oCorr++)
652  {
653  // Explicit pressure source contribution
654  tmp<volScalarField> tpu(this->pu());
655 
656  // Implicit pressure source coefficient
657  tmp<volScalarField> tpp(this->pp());
658 
659  // Solve for momentum for U_
660  tmp<fvVectorMatrix> UEqn = solveMomentum(tpu(), tpp());
661 
662  // Solve energy for hs_ - also updates thermo
663  solveEnergy();
664 
665  // Film thickness correction loop
666  for (int corr=1; corr<=nCorr_; corr++)
667  {
668  // Solve thickness for delta_
669  solveThickness(tpu(), tpp(), UEqn());
670  }
671  }
672 
673  // Update deltaRho_ with new delta_
674  deltaRho_ == delta_*rho_;
675 
676  // Update temperature using latest hs_
677  T_ == T(hs_);
678 
679  // Reset source terms for next time integration
681 }
682 
683 
685 {
686  return Cp_;
687 }
688 
689 
691 {
692  return kappa_;
693 }
694 
695 
697 {
698  return T_;
699 }
700 
701 
703 {
704  return Ts_;
705 }
706 
707 
709 {
710  return Tw_;
711 }
712 
713 
715 {
716  return hs_;
717 }
718 
719 
721 {
722  return primaryMassPCTrans_;
723 }
724 
725 
727 {
729 
730  const scalarField& Tinternal = T_;
731 
732  Info<< indent << "min/mean/max(T) = "
733  << gMin(Tinternal) << ", "
734  << gAverage(Tinternal) << ", "
735  << gMax(Tinternal) << nl;
736 
737  phaseChange_->info(Info);
738 }
739 
740 
742 {
744  (
746  (
747  IOobject
748  (
749  "thermoSingleLayer::Srho",
750  time().timeName(),
751  primaryMesh(),
754  false
755  ),
756  primaryMesh(),
758  )
759  );
760 
761  scalarField& Srho = tSrho.ref();
762  const scalarField& V = primaryMesh().V();
763  const scalar dt = time_.deltaTValue();
764 
766  {
767  const label filmPatchi = intCoupledPatchIDs()[i];
768 
769  scalarField patchMass =
770  primaryMassPCTrans_.boundaryField()[filmPatchi];
771 
772  toPrimary(filmPatchi, patchMass);
773 
774  const label primaryPatchi = primaryPatchIDs()[i];
775  const unallocLabelList& cells =
776  primaryMesh().boundaryMesh()[primaryPatchi].faceCells();
777 
778  forAll(patchMass, j)
779  {
780  Srho[cells[j]] = patchMass[j]/(V[cells[j]]*dt);
781  }
782  }
783 
784  return tSrho;
785 }
786 
787 
789 (
790  const label i
791 ) const
792 {
793  const label vapId = thermo_.carrierId(filmThermo_->name());
794 
796  (
798  (
799  IOobject
800  (
801  "thermoSingleLayer::Srho(" + Foam::name(i) + ")",
802  time_.timeName(),
803  primaryMesh(),
806  false
807  ),
808  primaryMesh(),
810  )
811  );
812 
813  if (vapId == i)
814  {
815  scalarField& Srho = tSrho.ref();
816  const scalarField& V = primaryMesh().V();
817  const scalar dt = time().deltaTValue();
818 
820  {
821  const label filmPatchi = intCoupledPatchIDs_[i];
822 
823  scalarField patchMass =
824  primaryMassPCTrans_.boundaryField()[filmPatchi];
825 
826  toPrimary(filmPatchi, patchMass);
827 
828  const label primaryPatchi = primaryPatchIDs()[i];
829  const unallocLabelList& cells =
830  primaryMesh().boundaryMesh()[primaryPatchi].faceCells();
831 
832  forAll(patchMass, j)
833  {
834  Srho[cells[j]] = patchMass[j]/(V[cells[j]]*dt);
835  }
836  }
837  }
838 
839  return tSrho;
840 }
841 
842 
844 {
846  (
848  (
849  IOobject
850  (
851  "thermoSingleLayer::Sh",
852  time().timeName(),
853  primaryMesh(),
856  false
857  ),
858  primaryMesh(),
860  )
861  );
862 /*
863  phase change energy fed back into the film...
864 
865  scalarField& Sh = tSh.ref();
866  const scalarField& V = primaryMesh().V();
867  const scalar dt = time_.deltaTValue();
868 
869  forAll(intCoupledPatchIDs_, i)
870  {
871  const label filmPatchi = intCoupledPatchIDs_[i];
872 
873  scalarField patchEnergy =
874  primaryEnergyPCTrans_.boundaryField()[filmPatchi];
875 
876  toPrimary(filmPatchi, patchEnergy);
877 
878  const label primaryPatchi = primaryPatchIDs()[i];
879  const unallocLabelList& cells =
880  primaryMesh().boundaryMesh()[primaryPatchi].faceCells();
881 
882  forAll(patchEnergy, j)
883  {
884  Sh[cells[j]] += patchEnergy[j]/(V[cells[j]]*dt);
885  }
886  }
887 */
888  return tSh;
889 }
890 
891 
892 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
893 
894 } // end namespace Foam
895 } // end namespace regionModels
896 } // end namespace surfaceFilmModels
897 
898 // ************************************************************************* //
volScalarField Ts_
Temperature - surface / [K].
virtual bool read()
Read control parameters from dictionary.
volScalarField Tw_
Temperature - wall / [K].
autoPtr< filmThermoModel > filmThermo_
Film thermo model.
bool isA(const Type &t)
Check if a dynamic_cast to typeid is possible.
Definition: typeInfo.H:134
autoPtr< heatTransferModel > htcw_
Heat transfer coefficient bewteen wall and film [W/m2/K].
virtual void correctAlpha()
Correct film coverage field.
virtual const volScalarField & Cp() const
Return the film specific heat capacity [J/kg/K].
virtual const volScalarField & kappa() const
Return the film thermal conductivity [W/m/K].
virtual const volScalarField & magSf() const
Return the face area magnitudes / [m2].
autoPtr< phaseChangeModel > phaseChange_
Phase change.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
virtual void resetPrimaryRegionSourceTerms()
Reset source term fields.
Ostream & indent(Ostream &os)
Indent stream.
Definition: Ostream.H:223
const Time & time() const
Return the reference to the time database.
Definition: regionModelI.H:37
scalar Tmax_
Maximum temperature limit (optional)
virtual void transferPrimaryRegionSourceFields()
Transfer source fields from the primary region to the film region.
Type gMin(const FieldField< Field, Type > &f)
static const char *const typeName
Definition: Field.H:94
tmp< fvMatrix< Type > > Sp(const DimensionedField< scalar, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:633
virtual void transferPrimaryRegionThermoFields()
Transfer thermo fields from the primary region to the film region.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
virtual tmp< DimensionedField< scalar, volMesh > > Sh() const
Return enthalpy source - Eulerian phase only.
static autoPtr< phaseChangeModel > New(surfaceFilmModel &owner, const dictionary &dict)
Return a reference to the selected phase change model.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:65
const dimensionedScalar phi0
Magnetic flux quantum: default SI units: [Wb].
virtual void transferPrimaryRegionThermoFields()
Transfer thermo fields from the primary region to the film region.
scalarField availableMass_
Available mass for transfer via sub-models.
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:715
bool readBool(Istream &)
Definition: boolIO.C:60
autoPtr< heatTransferModel > htcs_
Heat transfer coefficient bewteen film surface and primary.
virtual tmp< DimensionedField< scalar, volMesh > > Srho() const
Return total mass source - Eulerian phase only.
virtual void updateSurfaceTemperatures()
Correct the film surface and wall temperatures.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:421
virtual void solveThickness(const volScalarField &pu, const volScalarField &pp, const fvVectorMatrix &UEqn)
Solve coupled velocity-thickness equations.
Macros for easy insertion into run-time selection tables.
virtual void solveContinuity()
Solve continuity equation.
virtual tmp< volScalarField > pp()
Implicit pressure source coefficient.
volScalarField primaryEnergyPCTrans_
Film energy evolved via phase change.
const fvMesh & regionMesh() const
Return the region mesh database.
Definition: regionModelI.H:61
volScalarField deltaRho_
Film thickness*density (helper field) / [kg/m2].
const dictionary & coeffs() const
Return the model coefficients dictionary.
Definition: regionModelI.H:96
virtual void transferPrimaryRegionSourceFields()
Transfer source fields from the primary region to the film region.
dimensionedScalar pos(const dimensionedScalar &ds)
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:58
addToRunTimeSelectionTable(surfaceFilmModel, kinematicSingleLayer, mesh)
virtual void updateSubmodels()
Update the film sub-models.
const cellShapeList & cells
virtual tmp< fvVectorMatrix > solveMomentum(const volScalarField &pu, const volScalarField &pp)
Solve for film velocity.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
scalar hydrophilicDryScale_
Length scale applied to deltaWet_ to determine when a wet.
A class for handling words, derived from string.
Definition: word.H:59
SolverPerformance< Type > solve(fvMatrix< Type > &, const dictionary &)
Solve returning the solution statistics given convergence tolerance.
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
Calculate the face-flux of the given field.
virtual tmp< fvScalarMatrix > q(volScalarField &hs) const
Return the wall/surface heat transfer term for the enthalpy equation.
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
Definition: dimensionSets.H:52
Calculate the laplacian of the given field.
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvmDdt.C:46
volScalarField kappa_
Thermal conductivity / [W/m/K].
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
word timeName
Definition: getTimeIndex.H:3
Thermo package for (S)olids (L)iquids and (G)ases Takes reference to thermo package, and provides:
Definition: SLGThermo.H:62
scalar hydrophilicWetScale_
Length scale applied to deltaWet_ to determine when a dry.
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if succesful.
Definition: doubleScalar.H:63
bool hasMultiComponentCarrier() const
Thermo database has multi-component carrier flag.
Definition: SLGThermo.C:223
scalar Tmin_
Minimum temperature limit (optional)
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
const dimensionSet & dimensions() const
Return dimensions.
virtual tmp< volScalarField > pu()
Explicit pressure source contribution.
virtual tmp< volScalarField > primaryMassTrans() const
Return mass transfer source - Eulerian phase only.
const SLGThermo & thermo_
Reference to the SLGThermo.
Calculate the divergence of the given field.
virtual const volScalarField & T() const
Return the film mean temperature [K].
virtual void evolveRegion()
Evolve the film equations.
virtual void addSources(const label patchi, const label facei, const scalar massSource, const vector &momentumSource, const scalar pressureSource, const scalar energySource)
External hook to add sources to the film.
static const char nl
Definition: Ostream.H:262
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmDiv.C:46
Type gMax(const FieldField< Field, Type > &f)
const labelList & primaryPatchIDs() const
Return the list of patch IDs on the primary region coupled.
Definition: regionModelI.H:172
static autoPtr< filmRadiationModel > New(surfaceFilmModel &owner, const dictionary &dict)
Return a reference to the selected phase change model.
virtual void resetPrimaryRegionSourceTerms()
Reset source term fields.
scalar deltaTValue() const
Return time step value.
Definition: TimeStateI.H:41
virtual void correctThermoFields()
Correct the thermo fields.
autoPtr< filmRadiationModel > radiation_
Radiation.
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
const dimensionSet dimEnergy
dictionary coeffs_
Model coefficients dictionary.
Definition: regionModel.H:105
virtual bool read()
Read control parameters from dictionary.
virtual void updateSubmodels()
Update the film sub-models.
virtual const volScalarField & Tw() const
Return the film wall temperature [K].
List< word > wordList
A List of words.
Definition: fileName.H:54
label carrierId(const word &cmptName, bool allowNotFound=false) const
Index of carrier component.
Definition: SLGThermo.C:148
label patchi
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
dimensionedScalar deltaT() const
Return time step.
Definition: TimeStateI.H:53
scalar deltaWet_
Threshold film thickness beyond which the film is considered &#39;wet&#39;.
void toPrimary(const label regionPatchi, List< Type > &regionField) const
Convert a local region field to the primary region.
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
const labelUList & faceCells() const
Return face-cell addressing.
Definition: polyPatch.C:315
volScalarField Cp_
Specific heat capacity / [J/kg/K].
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Type gAverage(const FieldField< Field, Type > &f)
const Time & time_
Reference to the time database.
Definition: regionModel.H:90
surfaceScalarField phi_
Mass flux (includes film thickness) / [kg.m/s].
static autoPtr< filmViscosityModel > New(surfaceFilmModel &owner, const dictionary &dict, volScalarField &mu)
Return a reference to the selected phase change model.
fvVectorMatrix & UEqn
Definition: UEqn.H:13
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
volScalarField alpha_
Film coverage indicator, 1 = covered, 0 = uncovered / [].
A List with indirect addressing.
Definition: fvMatrix.H:106
virtual void preEvolveRegion()
Pre-evolve film hook.
static autoPtr< heatTransferModel > New(surfaceFilmModel &owner, const dictionary &dict)
Return a reference to the selected phase change model.
virtual void updateSurfaceVelocities()
Update film surface velocities.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const basicSpecieMixture & carrier() const
Return reference to the gaseous components.
Definition: SLGThermo.C:108
virtual const volScalarField & Ts() const
Return the film surface temperature [K].
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
PtrList< volScalarField > YPrimary_
List of specie mass fractions / [0-1].
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
virtual const volScalarField & hs() const
Return the film sensible enthalpy [J/kg].
void correctBoundaryConditions()
Correct boundary field.
messageStream Info
autoPtr< filmViscosityModel > viscosity_
Viscosity model.
const labelList & intCoupledPatchIDs() const
Return the list of patch IDs internally coupled with the.
Definition: regionModelI.H:179
const speciesTable & species() const
Return the table of species.
volScalarField primaryMassPCTrans_
Film mass evolved via phase change.
virtual void addSources(const label patchi, const label facei, const scalar massSource, const vector &momentumSource, const scalar pressureSource, const scalar energySource=0)
External hook to add sources to the film.
This boundary condition provides a self-contained version of the mapped condition. It does not use information on the patch; instead it holds thr data locally.
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
labelList intCoupledPatchIDs_
List of patch IDs internally coupled with the primary region.
Definition: regionModel.H:117
A class for managing temporary objects.
Definition: PtrList.H:54
tmp< surfaceScalarField > flux(const volVectorField &vvf)
Return the face-flux field obtained from the given volVectorField.
Definition: fvcFlux.C:32
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
const fvMesh & primaryMesh() const
Return the reference to the primary mesh database.
Definition: regionModelI.H:31
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
defineTypeNameAndDebug(kinematicSingleLayer, 0)
virtual void correctHsForMappedT()
Correct sensible enthalpy for mapped temperature fields.
volScalarField hs_
Sensible enthalpy / [J/kg].
const word & name() const
Return name.
Definition: IOobject.H:260
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:57
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:451
#define InfoInFunction
Report an information message using Foam::Info.