41 turbulentTemperatureTwoPhaseRadCoupledMixedFvPatchScalarField::regionType
44turbulentTemperatureTwoPhaseRadCoupledMixedFvPatchScalarField::regionTypeNames_
46 { regionType::solid,
"solid" },
47 { regionType::fluid,
"fluid" },
54 turbulentTemperatureTwoPhaseRadCoupledMixedFvPatchScalarField::KMethodType
57turbulentTemperatureTwoPhaseRadCoupledMixedFvPatchScalarField::KMethodTypeNames_
59 { KMethodType::mtSolidThermo,
"solidThermo" },
60 { KMethodType::mtLookup,
"lookup" },
61 { KMethodType::mtPhaseSystem,
"phaseSystem" }
75 const polyMesh&
mesh =
patch().boundaryMesh().mesh();
76 const label patchi =
patch().index();
82 const solidThermo&
thermo =
83 mesh.lookupObject<solidThermo>(basicThermo::dictName);
85 return thermo.kappa(patchi);
91 if (
mesh.foundObject<volScalarField>(kappaName_))
98 else if (
mesh.foundObject<volSymmTensorField>(kappaName_))
108 return n & KWall &
n;
113 <<
"Did not find field " << kappaName_
114 <<
" on mesh " <<
mesh.name() <<
" patch " <<
patch().name()
116 <<
" Please set 'kappa' to the name of a volScalarField"
117 <<
" or volSymmTensorField."
129 const phaseSystem&
fluid =
131 mesh.lookupObject<phaseSystem>(
"phaseProperties")
156 <<
"Unimplemented method " << KMethodTypeNames_[method_] <<
nl
157 <<
"Please set 'kappaMethod' to one of "
159 <<
"and 'kappa' to the name of the volScalar"
177turbulentTemperatureTwoPhaseRadCoupledMixedFvPatchScalarField::
178turbulentTemperatureTwoPhaseRadCoupledMixedFvPatchScalarField
181 const DimensionedField<scalar, volMesh>& iF
184 mixedFvPatchScalarField(
p, iF),
188 otherPhaseName_(
"vapor"),
189 TnbrName_(
"undefined-Tnbr"),
190 qrNbrName_(
"undefined-qrNbr"),
191 qrName_(
"undefined-qr")
193 this->refValue() = 0.0;
194 this->refGrad() = 0.0;
195 this->valueFraction() = 1.0;
199turbulentTemperatureTwoPhaseRadCoupledMixedFvPatchScalarField::
200turbulentTemperatureTwoPhaseRadCoupledMixedFvPatchScalarField
202 const turbulentTemperatureTwoPhaseRadCoupledMixedFvPatchScalarField& psf,
204 const DimensionedField<scalar, volMesh>& iF,
205 const fvPatchFieldMapper& mapper
208 mixedFvPatchScalarField(psf,
p, iF, mapper),
209 regionType_(psf.regionType_),
210 method_(psf.method_),
211 kappaName_(psf.kappaName_),
212 otherPhaseName_(psf.otherPhaseName_),
213 TnbrName_(psf.TnbrName_),
214 qrNbrName_(psf.qrNbrName_),
219turbulentTemperatureTwoPhaseRadCoupledMixedFvPatchScalarField::
220turbulentTemperatureTwoPhaseRadCoupledMixedFvPatchScalarField
223 const DimensionedField<scalar, volMesh>& iF,
224 const dictionary&
dict
227 mixedFvPatchScalarField(
p, iF),
228 regionType_(regionTypeNames_.
get(
"region",
dict)),
229 method_(KMethodTypeNames_.
get(
"kappaMethod",
dict)),
230 kappaName_(
dict.getOrDefault<word>(
"kappa",
"none")),
231 otherPhaseName_(
dict.
get<word>(
"otherPhase")),
232 TnbrName_(
dict.getOrDefault<word>(
"Tnbr",
"T")),
233 qrNbrName_(
dict.getOrDefault<word>(
"qrNbr",
"none")),
234 qrName_(
dict.getOrDefault<word>(
"qr",
"none"))
236 if (!isA<mappedPatchBase>(this->
patch().
patch()))
239 <<
"' not type '" << mappedPatchBase::typeName <<
"'"
240 <<
"\n for patch " <<
p.name()
241 <<
" of field " << internalField().name()
242 <<
" in file " << internalField().objectPath()
248 if (
dict.found(
"refValue"))
260 valueFraction() = 1.0;
265turbulentTemperatureTwoPhaseRadCoupledMixedFvPatchScalarField::
266turbulentTemperatureTwoPhaseRadCoupledMixedFvPatchScalarField
268 const turbulentTemperatureTwoPhaseRadCoupledMixedFvPatchScalarField& psf,
269 const DimensionedField<scalar, volMesh>& iF
272 mixedFvPatchScalarField(psf, iF),
273 regionType_(psf.regionType_),
274 method_(psf.method_),
275 kappaName_(psf.kappaName_),
276 otherPhaseName_(psf.otherPhaseName_),
277 TnbrName_(psf.TnbrName_),
278 qrNbrName_(psf.qrNbrName_),
285void turbulentTemperatureTwoPhaseRadCoupledMixedFvPatchScalarField::
293 const polyMesh&
mesh =
patch().boundaryMesh().mesh();
297 int oldTag = UPstream::msgType();
298 UPstream::msgType() = oldTag+1;
301 const label patchi =
patch().index();
302 const mappedPatchBase& mpp =
303 refCast<const mappedPatchBase>(
patch().
patch());
304 const polyMesh& nbrMesh = mpp.sampleMesh();
305 const label samplePatchi = mpp.samplePolyPatch().index();
306 const fvPatch& nbrPatch =
307 refCast<const fvMesh>(nbrMesh).boundary()[samplePatchi];
311 const turbulentTemperatureTwoPhaseRadCoupledMixedFvPatchScalarField&
313 <
const turbulentTemperatureTwoPhaseRadCoupledMixedFvPatchScalarField>
320 mpp.distribute(TcNbr);
325 KDeltaNbr = nbrField.kappa(nbrField)*nbrPatch.deltaCoeffs();
326 mpp.distribute(KDeltaNbr);
331 if (qrName_ !=
"none")
337 if (qrNbrName_ !=
"none")
339 qrNbr = nbrPatch.lookupPatchField<
volScalarField, scalar>(qrNbrName_);
340 mpp.distribute(qrNbr);
344 if (regionType_ == solid)
347 const phaseSystem&
fluid =
349 nbrMesh.lookupObject<phaseSystem>(
"phaseProperties")
353 const phaseModel& liquid
355 fluid.phases()[nbrField.internalField().group()]
358 const phaseModel& vapor(
fluid.phases()[otherPhaseName_]);
364 alphal*(liquid.kappaEff(samplePatchi))*nbrPatch.deltaCoeffs();
365 mpp.distribute(KDeltaLiqNbr);
370 alphav*(vapor.kappaEff(samplePatchi))*nbrPatch.deltaCoeffs();
371 mpp.distribute(KDeltaVapNbr);
375 vapor.thermo().
T().boundaryField()[samplePatchi];
376 TvNbr = Tv.patchInternalField();
377 mpp.distribute(TvNbr);
385 scalarField KDeltaLiqVapNbr(KDeltaLiqNbr + KDeltaVapNbr);
386 valueFraction() = KDeltaLiqVapNbr/(KDeltaLiqVapNbr + KDelta);
387 refValue() =
c/KDeltaLiqVapNbr;
388 refGrad() = (qr + qrNbr)/
kappa(Tp);
397 <<
" heat transfer rate from solid:" << Q
398 <<
" walltemperature "
399 <<
" min:" <<
gMin(Tp)
400 <<
" max:" <<
gMax(Tp)
405 else if (regionType_ ==
fluid)
407 const phaseSystem&
fluid =
409 mesh.lookupObject<phaseSystem>(
"phaseProperties")
412 const phaseModel& liquid
417 const phaseModel& vapor(
fluid.phases()[otherPhaseName_]);
420 vapor.thermo().
T().boundaryField()[patchi];
426 alphav*(vapor.kappaEff(patchi))*
patch().deltaCoeffs()
433 alphal*(liquid.kappaEff(patchi))*
patch().deltaCoeffs()
438 const scalarField c(TcNbr*KDeltaNbr + Tv.patchInternalField()*KdeltaVap);
442 valueFraction() = a/(a + KdeltaLiq);
444 refGrad() = (qr + qrNbr)/
kappa(Tp);
450 scalarField qVap((Tp - Tv.patchInternalField())*KdeltaVap);
457 scalar QLiq =
gSum(qLiq*
patch().magSf());
458 scalar QVap =
gSum(qVap*
patch().magSf());
460 Info<<
" Heat transfer to Liq: " << QLiq <<
endl;
461 Info<<
" Heat transfer to Vap: " << QVap <<
endl;
463 Info<<
" walltemperature "
464 <<
" min:" <<
gMin(Tp)
465 <<
" max:" <<
gMax(Tp)
473 <<
"Unknown phase type. Valid types are: "
474 << regionTypeNames_ <<
nl <<
exit(FatalError);
480 UPstream::msgType() = oldTag;
484void turbulentTemperatureTwoPhaseRadCoupledMixedFvPatchScalarField::write
490 os.writeEntry(
"kappaMethod", KMethodTypeNames_[method_]);
491 os.writeEntryIfDifferent<word>(
"kappa",
"none", kappaName_);
493 os.writeEntry(
"Tnbr", TnbrName_);
495 os.writeEntryIfDifferent<word>(
"qrNbr",
"none", qrNbrName_);
496 os.writeEntryIfDifferent<word>(
"qr",
"none", qrName_);
498 os.writeEntry(
"region", regionTypeNames_[regionType_]);
499 os.writeEntry(
"otherPhase", otherPhaseName_);
508 turbulentTemperatureTwoPhaseRadCoupledMixedFvPatchScalarField
Macros for easy insertion into run-time selection tables.
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
tmp< Field< Type > > T() const
Return the field transpose (only defined for second rank tensors)
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
virtual bool write()
Write the output fields.
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
A class for managing temporary objects.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
OBJstream os(runTime.globalPath()/outputName)
#define makePatchTypeField(PatchTypeField, typePatchTypeField)
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
constexpr const char *const group
Group name for atomic constants.
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
const dimensionedScalar c
Speed of light in a vacuum.
const std::string patch
OpenFOAM patch number as a std::string.
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Type gSum(const FieldField< Field, Type > &f)
To & refCast(From &r)
Reference type cast template function.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
messageStream Info
Information stream (stdout output on master, null elsewhere)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Ostream & endl(Ostream &os)
Add newline and flush stream.
FlatOutput::OutputAdaptor< Container, Delimiters > flatOutput(const Container &obj, Delimiters delim)
Global flatOutput() function with specified output delimiters.
Field< vector > vectorField
Specialisation of Field<T> for vector.
GeometricField< symmTensor, fvPatchField, volMesh > volSymmTensorField
Field< symmTensor > symmTensorField
Specialisation of Field<T> for symmTensor.
Type gAverage(const FieldField< Field, Type > &f)
Type gMin(const FieldField< Field, Type > &f)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Type gMax(const FieldField< Field, Type > &f)
fvPatchField< scalar > fvPatchScalarField
constexpr char nl
The newline '\n' character (0x0a)
#define forAll(list, i)
Loop across all elements in list.