Go to the documentation of this file.
45 #include "heatTransferModel.H"
53 namespace regionModels
55 namespace surfaceFilmModels
66 wordList thermoSingleLayer::hsBoundaryTypes()
78 bTypes[patchi] = fixedValueFvPatchField<scalar>::typeName;
128 hsBf[patchi] ==
hs(Tp, patchi);
168 YPrimary_[i].correctBoundaryConditions();
182 volScalarField::Boundary& hsSpPrimaryBf =
187 forAll(hsSpPrimaryBf, patchi)
194 hsSpPrimaryBf[patchi] *= rpriMagSfdeltaT;
218 else if ((
alpha_[i] > 0.5) && (
delta_[i] < hydrophilicDry))
327 thermoSingleLayer::thermoSingleLayer
329 const word& modelType,
332 const word& regionType,
350 zeroGradientFvPatchScalarField::typeName
364 zeroGradientFvPatchScalarField::typeName
390 zeroGradientFvPatchScalarField::typeName
403 zeroGradientFvPatchScalarField::typeName
424 "primaryEnergyTrans",
432 zeroGradientFvPatchScalarField::typeName
435 deltaWet_(coeffs_.get<scalar>(
"deltaWet")),
436 hydrophilic_(coeffs_.get<
bool>(
"hydrophilic")),
437 hydrophilicDryScale_(0.0),
438 hydrophilicWetScale_(0.0),
452 this->mappedPushedFieldPatchTypes<scalar>()
481 this->mappedFieldAndInternalPatchTypes<scalar>()
500 if (coeffs().readIfPresent(
"Tmin", Tmin_))
502 Info<<
" limiting minimum temperature to " << Tmin_ <<
endl;
505 if (coeffs().readIfPresent(
"Tmax", Tmax_))
507 Info<<
" limiting maximum temperature to " << Tmax_ <<
endl;
510 if (thermo_.hasMultiComponentCarrier())
512 YPrimary_.setSize(thermo_.carrier().species().size());
514 forAll(thermo_.carrier().species(), i)
523 thermo_.carrier().species()[i],
531 pSp_.boundaryField().types()
539 coeffs_.readEntry(
"hydrophilicDryScale", hydrophilicDryScale_);
540 coeffs_.readEntry(
"hydrophilicWetScale", hydrophilicWetScale_);
545 transferPrimaryRegionThermoFields();
549 correctThermoFields();
554 deltaRho_ == delta_*rho_;
573 viscosity_->correct(pPrimary_, T_);
590 const scalar massSource,
591 const vector& momentumSource,
592 const scalar pressureSource,
593 const scalar energySource
608 Info<<
" energy = " << energySource <<
nl <<
endl;
611 hsSpPrimary_.boundaryFieldRef()[patchi][facei] -= energySource;
658 for (
int corr=1; corr<=
nCorr_; corr++)
716 <<
gMin(Tinternal) <<
", "
776 const label vapId = thermo_.carrierId(filmThermo_->name());
800 const scalar dt = time().deltaTValue();
802 forAll(intCoupledPatchIDs_, i)
804 const label filmPatchi = intCoupledPatchIDs_[i];
807 primaryMassTrans_.boundaryField()[filmPatchi];
809 toPrimary(filmPatchi, patchMass);
811 const label primaryPatchi = primaryPatchIDs()[i];
813 primaryMesh().boundaryMesh()[primaryPatchi].faceCells();
817 Srho[
cells[j]] += patchMass[j]/(V[
cells[j]]*dt);
int debug
Static debugging option.
virtual void solveEnergy()
Solve energy equation.
volScalarField Tw_
Temperature - wall [K].
virtual bool read()
Read control parameters from dictionary.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
#define InfoInFunction
Report an information message using Foam::Info.
virtual bool read()
Read control parameters from dictionary.
virtual void transferPrimaryRegionThermoFields()
Transfer thermo fields from the primary region to the film region.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
virtual void resetPrimaryRegionSourceTerms()
Reset source term fields.
virtual void info()
Provide some feedback.
virtual const volScalarField & hs() const
Return the film sensible enthalpy [J/kg].
A class for handling words, derived from Foam::string.
virtual void correctHsForMappedT()
Correct sensible enthalpy for mapped temperature fields.
virtual ~thermoSingleLayer()
Destructor.
virtual void correct(scalarField &availableMass, volScalarField &massToTransfer)
Correct kinematic transfers.
injectionModelList injection_
Cloud injection.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
virtual void preEvolveRegion()
Pre-evolve film hook.
Thermo package for (S)olids (L)iquids and (G)ases Takes reference to thermo package,...
volScalarField primaryMassTrans_
Film mass available for transfer to the primary region.
autoPtr< filmViscosityModel > viscosity_
Viscosity model.
autoPtr< heatTransferModel > htcs_
Heat transfer coefficient between film surface and primary.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero.
const dimensionSet dimEnergy
Type gAverage(const FieldField< Field, Type > &f)
autoPtr< filmTurbulenceModel > turbulence_
Turbulence model.
autoPtr< filmThermoModel > filmThermo_
Film thermo model.
virtual void solveContinuity()
Solve continuity equation.
static autoPtr< heatTransferModel > New(surfaceFilmRegionModel &film, const dictionary &dict)
Return a reference to the selected phase change model.
volScalarField primaryEnergyTrans_
Film energy transfer.
autoPtr< phaseChangeModel > phaseChange_
Phase change.
Calculate the divergence of the given field.
virtual const volScalarField & Cp() const
Return the film specific heat capacity [J/kg/K].
const Time & time() const
Return the reference to the time database.
Kinematic form of single-cell layer surface film model.
virtual void solveThickness(const volScalarField &pu, const volScalarField &pp, const fvVectorMatrix &UEqn)
Solve coupled velocity-thickness equations.
volScalarField rho_
Density [kg/m3].
volScalarField sigma_
Surface tension [m/s2].
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Ostream & endl(Ostream &os)
Add newline and flush stream.
virtual const volScalarField & kappa() const
Return the film thermal conductivity [W/m/K].
dimensionedScalar pos0(const dimensionedScalar &ds)
volScalarField T_
Temperature - mean [K].
addToRunTimeSelectionTable(surfaceFilmRegionModel, kinematicSingleLayer, mesh)
volScalarField Cp_
Specific heat capacity [J/kg/K].
scalar hydrophilicDryScale_
Length scale applied to deltaWet_ to determine when a wet.
scalarField availableMass_
Available mass for transfer via sub-models.
volScalarField pPrimary_
Pressure [Pa].
#define forAll(list, i)
Loop across all elements in list.
virtual void correctThermoFields()
Correct the thermo fields.
virtual void evolveRegion()
Evolve the film equations.
const fvMesh & primaryMesh() const
Return the reference to the primary mesh database.
List< word > wordList
A List of words.
void toPrimary(const label regionPatchi, List< Type > ®ionField) const
Convert a local region field to the primary region.
scalar deltaWet_
Threshold film thickness beyond which the film is considered 'wet'.
volScalarField cloudDiameterTrans_
Parcel diameters originating from film to cloud.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
volScalarField hsSp_
Energy [J/m2/s].
volScalarField deltaRho_
Film thickness*density (helper field) [kg/m2].
virtual void transferPrimaryRegionThermoFields()
Transfer thermo fields from the primary region to the film region.
const dimensionSet dimArea(sqr(dimLength))
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
virtual void transferPrimaryRegionSourceFields()
Transfer source fields from the primary region to the film region.
volScalarField delta_
Film thickness [m].
messageStream Info
Information stream (uses stdout - output is on the master only)
A patch is a list of labels that address the faces in the global face list.
word name(const complex &c)
Return string representation of complex.
virtual const volScalarField & magSf() const
Return the face area magnitudes / [m2].
const fvMesh & regionMesh() const
Return the region mesh database.
virtual tmp< volScalarField > pp()
Implicit pressure source coefficient.
virtual void info()
Provide some feedback.
const dimensionedScalar phi0
Magnetic flux quantum: default SI units: [Wb].
SolverPerformance< Type > solve(faMatrix< Type > &, Istream &)
Solve returning the solution statistics given convergence tolerance.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
virtual void resetPrimaryRegionSourceTerms()
Reset source term fields.
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 autoPtr< phaseChangeModel > New(surfaceFilmRegionModel &film, const dictionary &dict)
Return a reference to the selected phase change model.
volScalarField Ts_
Temperature - surface [K].
virtual tmp< fvVectorMatrix > solveMomentum(const volScalarField &pu, const volScalarField &pp)
Solve for film velocity.
virtual void updateSubmodels()
Update the film sub-models.
virtual const volScalarField & T() const
Return the film mean temperature [K].
PtrList< volScalarField > YPrimary_
List of specie mass fractions [0-1].
Macros for easy insertion into run-time selection tables.
volScalarField TPrimary_
Temperature [K].
Mesh data needed to do the Finite Volume discretisation.
const uniformDimensionedVectorField & g
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
label nCorr_
Number of PISO-like correctors.
Ostream & indent(Ostream &os)
Indent stream.
virtual tmp< fvScalarMatrix > q(volScalarField &hs) const
Return the wall/surface heat transfer term for the enthalpy equation.
virtual void preEvolveRegion()
Pre-evolve film hook.
virtual tmp< volScalarField > pu()
Explicit pressure source contribution.
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const wordHashSet &selectedFields, LIFOStack< regIOobject * > &storedObjects)
Read the selected GeometricFields of the templated type.
volScalarField pSp_
Pressure [Pa].
surfaceScalarField phi_
Mass flux (includes film thickness) [kg.m/s].
const labelList & primaryPatchIDs() const
Return the list of patch IDs on the primary region coupled.
void correctBoundaryConditions()
Correct boundary field.
volScalarField hsSpPrimary_
Energy [J/m2/s].
volScalarField hs_
Sensible enthalpy [J/kg].
defineTypeNameAndDebug(kinematicSingleLayer, 0)
const TargetType * isA(const Type &t)
Check if dynamic_cast to TargetType is possible.
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.
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const Time & time_
Reference to the time database.
Calculate the laplacian of the given field.
virtual tmp< volScalarField::Internal > Srho() const
Return total mass source - Eulerian phase only.
This boundary condition provides a self-contained version of the mapped condition....
static autoPtr< filmRadiationModel > New(surfaceFilmRegionModel &film, const dictionary &dict)
Return a reference to the selected phase change model.
virtual void correct(scalarField &availableMass, volScalarField &massToInject, volScalarField &diameterToInject)
Correct.
volScalarField alpha_
Film coverage indicator, 1 = covered, 0 = uncovered [].
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
scalar deltaTValue() const
Return time step value.
Calculate the face-flux of the given field.
volScalarField rhoPrimary_
Density [kg/m3].
const dimensionedScalar e
Elementary charge.
volScalarField kappa_
Thermal conductivity [W/m/K].
bool hydrophilic_
Activation flag.
autoPtr< filmRadiationModel > radiation_
Radiation.
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
transferModelList transfer_
Transfer with the continuous phase.
labelList intCoupledPatchIDs_
List of patch IDs internally coupled with the primary region.
A List with indirect addressing.
virtual void correctAlpha()
Correct film coverage field.
scalar hydrophilicWetScale_
Length scale applied to deltaWet_ to determine when a dry.
Calculate the first temporal derivative.
virtual tmp< volScalarField::Internal > Sh() const
Return enthalpy source - Eulerian phase only.
virtual const volScalarField & Tw() const
Return the film wall temperature [K].
autoPtr< heatTransferModel > htcw_
Heat transfer coefficient between wall and film [W/m2/K].
Type gMin(const FieldField< Field, Type > &f)
virtual const volScalarField & Ts() const
Return the film surface temperature [K].
tmp< surfaceScalarField > flux(const volVectorField &vvf)
Return the face-flux field obtained from the given volVectorField.
const dimensionSet dimVolume(pow3(dimLength))
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
const labelList & intCoupledPatchIDs() const
Return the list of patch IDs internally coupled with the.
volScalarField rhoSp_
Mass [kg/m2/s].
volScalarField cloudMassTrans_
Film mass available for transfer to cloud.
Type gMax(const FieldField< Field, Type > &f)
static autoPtr< filmViscosityModel > New(surfaceFilmRegionModel &film, const dictionary &dict, volScalarField &mu)
Return a reference to the selected phase change model.
label nOuterCorr_
Number of outer correctors.
virtual void updateSurfaceTemperatures()
Correct the film surface and wall temperatures.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
virtual void transferPrimaryRegionSourceFields()
Transfer source fields from the primary region to the film region.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.