33 #include "phaseSystem.H"
52 alphatWallBoilingWallFunctionFvPatchScalarField::phaseTypeNames_
54 { phaseType::vaporPhase,
"vapor" },
55 { phaseType::liquidPhase,
"liquid" },
68 alphatWallBoilingWallFunctionFvPatchScalarField::
69 alphatWallBoilingWallFunctionFvPatchScalarField
76 otherPhaseName_(
"vapor"),
77 phaseType_(liquidPhase),
80 alphatConv_(
p.size(), 0),
81 dDep_(
p.size(), 1
e-5),
84 partitioningModel_(
nullptr),
85 nucleationSiteModel_(
nullptr),
86 departureDiamModel_(
nullptr),
87 departureFreqModel_(
nullptr),
88 filmBoilingModel_(
nullptr),
89 LeidenfrostModel_(
nullptr),
91 CHFSoobModel_(
nullptr),
96 AbyV_ = this->
patch().magSf();
99 const label faceCelli = this->
patch().faceCells()[facei];
100 AbyV_[facei] /= iF.
mesh().V()[faceCelli];
105 alphatWallBoilingWallFunctionFvPatchScalarField::
106 alphatWallBoilingWallFunctionFvPatchScalarField
114 otherPhaseName_(
dict.get<
word>(
"otherPhase")),
115 phaseType_(phaseTypeNames_.get(
"phaseType",
dict)),
118 alphatConv_(
p.size(), 0),
119 dDep_(
p.size(), 1
e-5),
122 partitioningModel_(
nullptr),
123 nucleationSiteModel_(
nullptr),
124 departureDiamModel_(
nullptr),
125 departureFreqModel_(
nullptr),
126 filmBoilingModel_(
nullptr),
127 LeidenfrostModel_(
nullptr),
129 CHFSoobModel_(
nullptr),
136 if (internalField().
group() == otherPhaseName_)
139 <<
"otherPhase should be the name of the vapor phase that "
140 <<
"corresponds to the liquid base of vice versa" <<
nl
141 <<
"This phase: " << internalField().group() <<
nl
142 <<
"otherPhase: " << otherPhaseName_
153 dict.subDict(
"partitioningModel")
165 dict.subDict(
"partitioningModel")
168 nucleationSiteModel_ =
171 dict.subDict(
"nucleationSiteModel")
174 departureDiamModel_ =
177 dict.subDict(
"departureDiamModel")
180 departureFreqModel_ =
183 dict.subDict(
"departureFreqModel")
187 dict.findDict(
"LeidenfrostModel");
203 const dictionary* HFSubCoolDict =
dict.findDict(
"CHFSubCoolModel");
235 if (
dict.found(
"dDep"))
240 dict.readIfPresent(
"K", K_);
242 dict.readIfPresent(
"wp", wp_);
244 if (
dict.found(
"qQuenching"))
253 if (
dict.found(
"alphatConv"))
258 AbyV_ = this->
patch().magSf();
261 const label faceCelli = this->
patch().faceCells()[facei];
262 AbyV_[facei] /= iF.
mesh().V()[faceCelli];
267 alphatWallBoilingWallFunctionFvPatchScalarField::
268 alphatWallBoilingWallFunctionFvPatchScalarField
283 otherPhaseName_(psf.otherPhaseName_),
284 phaseType_(psf.phaseType_),
285 relax_(psf.relax_.clone()),
287 alphatConv_(psf.alphatConv_, mapper),
288 dDep_(psf.dDep_, mapper),
289 qq_(psf.qq_, mapper),
291 partitioningModel_(psf.partitioningModel_),
292 nucleationSiteModel_(psf.nucleationSiteModel_),
293 departureDiamModel_(psf.departureDiamModel_),
294 filmBoilingModel_(psf.filmBoilingModel_),
295 LeidenfrostModel_(psf.LeidenfrostModel_),
296 CHFModel_(psf.CHFModel_),
297 CHFSoobModel_(psf.CHFSoobModel_),
298 MHFModel_(psf.MHFModel_),
299 TDNBModel_(psf.TDNBModel_),
304 alphatWallBoilingWallFunctionFvPatchScalarField::
305 alphatWallBoilingWallFunctionFvPatchScalarField
311 otherPhaseName_(psf.otherPhaseName_),
312 phaseType_(psf.phaseType_),
313 relax_(psf.relax_.clone()),
315 alphatConv_(psf.alphatConv_),
319 partitioningModel_(psf.partitioningModel_),
320 nucleationSiteModel_(psf.nucleationSiteModel_),
321 departureDiamModel_(psf.departureDiamModel_),
322 filmBoilingModel_(psf.filmBoilingModel_),
323 LeidenfrostModel_(psf.LeidenfrostModel_),
324 CHFModel_(psf.CHFModel_),
325 CHFSoobModel_(psf.CHFSoobModel_),
326 MHFModel_(psf.MHFModel_),
327 TDNBModel_(psf.TDNBModel_),
332 alphatWallBoilingWallFunctionFvPatchScalarField::
333 alphatWallBoilingWallFunctionFvPatchScalarField
340 otherPhaseName_(psf.otherPhaseName_),
341 phaseType_(psf.phaseType_),
342 relax_(psf.relax_.clone()),
344 alphatConv_(psf.alphatConv_),
348 partitioningModel_(psf.partitioningModel_),
349 nucleationSiteModel_(psf.nucleationSiteModel_),
350 departureDiamModel_(psf.departureDiamModel_),
351 filmBoilingModel_(psf.filmBoilingModel_),
352 LeidenfrostModel_(psf.LeidenfrostModel_),
353 CHFModel_(psf.CHFModel_),
354 CHFSoobModel_(psf.CHFSoobModel_),
355 MHFModel_(psf.MHFModel_),
356 TDNBModel_(psf.TDNBModel_),
363 bool alphatWallBoilingWallFunctionFvPatchScalarField::
376 const scalarField& alphatWallBoilingWallFunctionFvPatchScalarField::
386 <<
" dmdt requested for invalid phasePair!"
393 const scalarField& alphatWallBoilingWallFunctionFvPatchScalarField::
403 <<
" mDotL requested for invalid phasePair!"
410 void alphatWallBoilingWallFunctionFvPatchScalarField::updateCoeffs()
419 if (!partitioningModel_)
422 <<
"partitioningModel has not been constructed!"
428 refCast<const phaseSystem>
430 db().lookupObject<phaseSystem>(
"phaseProperties")
436 const label patchi =
patch().index();
438 const scalar t = this->db().time().timeOutputValue();
439 scalar
relax = relax_->value(t);
447 fluid.phases()[internalField().group()]
451 vapor.thermo().he().boundaryField()[patchi];
454 const scalarField vaporw(vapor.boundaryField()[patchi]);
460 partitioningModel_->fLiquid(1-vaporw)
473 this->operator[](i) =
475 (1 - fLiquid[i])*(alphatv[i])
476 /
max(vaporw[i], scalar(1
e-8))
484 Info<<
" alphatEffv: " <<
gMin(vaporw*(*
this + alphaw))
485 <<
" - " <<
gMax(vaporw*(*
this + alphaw)) <<
endl;
487 const scalarField qEff(vaporw*(*
this + alphaw)*hewv.snGrad());
489 scalar Qeff =
gSum(qEff*
patch().magSf());
490 Info<<
" Effective heat transfer rate to vapor:" << Qeff
498 if (!nucleationSiteModel_)
501 <<
"nucleationSiteModel has not been constructed!"
506 if (!departureDiamModel_)
509 <<
"departureDiameterModel has not been constructed!"
514 if (!departureFreqModel_)
517 <<
"departureFrequencyModel has not been constructed!"
523 fluid.phases()[internalField().group()]
550 const scalar Cmu25(
pow025(Cmu_));
565 turbModel.U().boundaryField()[patchi];
570 turbModel.rho().boundaryField()[patchi];
574 liquid.thermo().T().boundaryField()[patchi];
589 const scalarField yPlusTherm(this->yPlusTherm(P, Prat));
592 vaporTurbModel.rho().boundaryField()[patchi];
600 satModel.
Tsat(liquid.thermo().p());
604 const scalarField Tsatc(Tsatw.patchInternalField());
607 liquid.thermo().p().boundaryField()[patchi];
610 liquid.thermo().he().boundaryField()[patchi];
614 liquid.thermo().he().member() ==
"e"
621 vapor.thermo().he().member() ==
"e"
622 ? vapor.thermo().he(pw, Tsatc, patchi) + pw/rhoVaporw - hw
623 : vapor.thermo().he(pw, Tsatc, patchi) - hw
627 const scalarField liquidw(liquid.boundaryField()[patchi]);
629 const scalarField fLiquid(partitioningModel_->fLiquid(liquidw));
635 Prt_*(
log(E_*250)/kappa_ + P)
639 Tl =
max(Tc - 40, Tl);
650 && CHFSoobModel_.valid()
651 && TDNBModel_.valid()
653 && LeidenfrostModel_.valid()
654 && filmBoilingModel_.valid()
674 CHFSoobModel_->CHFSubCool
712 LeidenfrostModel_->TLeid
724 filmBoilingModel_->htcFilmBoil
743 wp_*(Tw - tDNB)/(TLeiden - tDNB),
750 Qtb = CHFtotal*(1 -
phi) +
phi*MHF;
758 nucleationSiteModel_->N
770 dDep_ = departureDiamModel_->dDeparture
783 departureFreqModel_->fDeparture
793 alphatConv_ = calcAlphat(alphatConv_);
807 if (Tw[i] > Tsatw[i])
813 regimeTypes[i] = regimeType::subcool;
817 rhow[i]*Cpw[i]*(Tsatw[i] - Tl[i])
818 /(rhoVaporw[i]*
L[i]);
821 fLiquid[i]*4.8*
exp(
min(-Ja/80,
log(VGREAT)));
825 min(
pi*
sqr(dDep_[i])*
N[i]*Al/4, scalar(1))
828 A1[i] =
max(1 - A2, scalar(1
e-4));
835 pi*
sqr(dDep_[i])*
N[i]*Al/4,
845 +
relax*(1.0/6.0)*A2E*dDep_[i]*rhoVaporw[i]
851 mDotL_[i] = dmdt_[i]*
L[i];
856 2*(alphaw[i]*Cpw[i])*fDep[i]
859 (0.8/
max(fDep[i], SMALL))/(
pi*alphaw[i]/rhow[i])
867 +
relax*A2*hQ*
max(Tw[i] - Tl[i], scalar(0))
870 this->operator[](i) =
876 (qq_[i] + mDotL_[i]/AbyV_[i])
877 /
max(hewSn[i], scalar(1
e-16))
879 /
max(liquidw[i], scalar(1
e-8)),
884 else if (Tw[i] > tDNB[i] && Tw[i] < TLeiden[i])
887 regimeTypes[i] = regimeType::transient;
890 alphatConv_[i] = 0.0;
896 relax*Qtb[i]*AbyV_[i]/
L[i]
897 + (1 -
relax)*dmdt_[i]
900 mDotL_[i] = dmdt_[i]*
L[i];
906 this->operator[](i) =
911 /
max(hewSn[i], scalar(1
e-16))
912 )/
max(liquidw[i], scalar(1
e-8)),
917 else if (Tw[i] > TLeiden[i])
919 regimeTypes[i] = regimeType::film;
922 alphatConv_[i] = 0.0;
928 relax*htcFilmBoiling[i]
929 *
max(Tw[i] - Tsatw[i], 0)
931 + (1 -
relax)*dmdt_[i]
935 mDotL_[i] = dmdt_[i]*
L[i];
942 mDotL_[i]/AbyV_[i]/
max(hewSn[i], scalar(1
e-16))
948 this->operator[](i) =
950 alphaFilm[i]/
max(liquidw[i], scalar(1
e-8))
959 regimeTypes[i] = regimeType::nonBoiling;
966 this->operator[](i) =
970 fLiquid[i]*(alphatConv_[i])
971 /
max(liquidw[i], scalar(1
e-8)),
982 liquidw*(*
this + alphaw)*hew.
snGrad()
987 Info<<
" alphatl: " <<
gMin((*
this)) <<
" - "
990 Info<<
" dmdt: " <<
gMin((dmdt_)) <<
" - "
993 Info<<
" alphatlEff: " <<
gMin(liquidw*(*
this + alphaw))
994 <<
" - " <<
gMax(liquidw*(*
this + alphaw)) <<
endl;
996 scalar Qeff =
gSum(qEff*
patch().magSf());
997 Info<<
" Effective heat transfer rate to liquid: " << Qeff
1010 switch (regimeTypes[i])
1012 case regimeType::subcool:
1016 case regimeType::transient:
1020 case regimeType::film:
1024 case regimeType::nonBoiling:
1025 nNonBoilings[i] = 1;
1030 scalar nSubCool(
gSum(nSubCools));
1031 scalar nTransient(
gSum(nTransients));
1032 scalar nFilm(
gSum(nFilms));
1033 scalar nNonBoiling(
gSum(nNonBoilings));
1037 Info<<
" sub Cool faces : " << nSubCool <<
endl;
1038 Info<<
" transient faces : " << nTransient <<
endl;
1039 Info<<
" film faces : " << nFilm <<
endl;
1040 Info<<
" non-Boiling faces : " << nNonBoiling <<
endl;
1041 Info<<
" total faces : "
1042 << nSubCool + nTransient + nFilm + nNonBoiling
1047 nNonBoilings*fLiquid*A1*(alphatConv_ + alphaw)
1052 Info<<
" Convective heat transfer: " << Qc <<
endl;
1056 relax*fLiquid*nFilms*htcFilmBoiling*(Tw - Tsatw)
1059 scalar QFilm =
gSum(qFilm*
patch().magSf());
1060 Info<<
" Film boiling heat transfer: " << QFilm <<
endl;
1062 Info<<
" Htc Film Boiling coeff: "
1063 <<
gMin(nFilms*htcFilmBoiling)
1065 <<
gMax(nFilms*htcFilmBoiling) <<
endl;
1068 gSum(fLiquid*nTransients*Qtb*
patch().magSf());
1069 Info<<
" Transient boiling heat transfer:" << Qtbtot
1080 A1*alphatConv_*hew.
snGrad()
1086 scalar QsubCool =
gSum(qSubCool*
patch().magSf());
1088 Info<<
" Sub Cool boiling heat transfer: " << QsubCool
1091 Info<<
" N: " <<
gMin(nSubCools*
N) <<
" - "
1093 Info<<
" dDep: " <<
gMin(nSubCools*dDep_) <<
" - "
1095 Info<<
" fDep: " <<
gMin(nSubCools*fDep) <<
" - "
1097 Info<<
" A1: " <<
gMin(nSubCools*A1) <<
" - "
1108 <<
"Unknown phase type. Valid types are: "
1113 fixedValueFvPatchScalarField::updateCoeffs();
1121 os.
writeEntry(
"phaseType", phaseTypeNames_[phaseType_]);
1123 relax_->writeData(os);
1130 partitioningModel_->
write(os);
1133 if (filmBoilingModel_)
1136 filmBoilingModel_->
write(os);
1140 if (LeidenfrostModel_)
1143 LeidenfrostModel_->
write(os);
1152 partitioningModel_->
write(os);
1156 nucleationSiteModel_->
write(os);
1160 departureDiamModel_->
write(os);
1164 departureFreqModel_->
write(os);
1167 if (filmBoilingModel_)
1170 filmBoilingModel_->
write(os);
1174 if (LeidenfrostModel_)
1177 LeidenfrostModel_->
write(os);
1184 CHFModel_->
write(os);
1191 CHFSoobModel_->
write(os);
1198 MHFModel_->
write(os);
1205 TDNBModel_->
write(os);
1215 os.
writeEntry(
"otherPhase", otherPhaseName_);
1221 writeEntry(
"value", os);