BUG: Modifying approach for external radiation with layers. Fixes #2476

This commit is contained in:
sergio
2022-07-05 14:42:56 -07:00
parent 5589108d73
commit 490f02fad4

View File

@ -6,7 +6,7 @@
\\/ M anipulation | \\/ M anipulation |
------------------------------------------------------------------------------- -------------------------------------------------------------------------------
Copyright (C) 2011-2017 OpenFOAM Foundation Copyright (C) 2011-2017 OpenFOAM Foundation
Copyright (C) 2015-2020 OpenCFD Ltd. Copyright (C) 2015-2022 OpenCFD Ltd.
------------------------------------------------------------------------------- -------------------------------------------------------------------------------
License License
This file is part of OpenFOAM. This file is part of OpenFOAM.
@ -376,15 +376,15 @@ void Foam::externalWallHeatFluxTemperatureFvPatchScalarField::updateCoeffs()
} }
} }
} }
scalarField hp(1/(1/htcCoeff + totalSolidRes));
const scalar Ta = const scalar Ta =
Ta_->value(this->db().time().timeOutputValue()); Ta_->value(this->db().time().timeOutputValue());
scalarField hpTa(hp*Ta); scalarField hrad(Tp.size(), Zero);
if (emissivity_ > 0) if (emissivity_ > 0)
{ {
const scalar eSig(emissivity_*sigma.value());
// Evaluate the radiative flux to the environment // Evaluate the radiative flux to the environment
// from the surface temperature ... // from the surface temperature ...
if (totalSolidRes > 0) if (totalSolidRes > 0)
@ -392,21 +392,48 @@ void Foam::externalWallHeatFluxTemperatureFvPatchScalarField::updateCoeffs()
// ... including the effect of the solid wall thermal // ... including the effect of the solid wall thermal
// resistance // resistance
scalarField TpLambda(htcCoeff/(htcCoeff + 1/totalSolidRes)); scalarField TpLambda(htcCoeff/(htcCoeff + 1/totalSolidRes));
scalarField Ts(TpLambda*Tp + (1 - TpLambda)*Ta); scalarField Ts(TpLambda*Ta + (1 - TpLambda)*Tp);
scalarField lambdaTa4(pow4((1 - TpLambda)*Ta)); hrad = eSig*((pow3(Ta) + pow3(Ts)) + Ta*Ts*(Ta + Ts));
hp += emissivity_*sigma.value()*(pow4(Ts) - lambdaTa4)/Tp; forAll(hrad, i)
hpTa += emissivity_*sigma.value()*(lambdaTa4 + pow4(Ta)); {
scalar hradTmp0 = hrad[i];
scalar TaLambda =
(htcCoeff[i] + hradTmp0)
/(htcCoeff[i] + hradTmp0 + 1/totalSolidRes);
scalar TsiNew = TaLambda*Ta + (1 - TaLambda)*Tp[i];
scalar Tsi = Ts[i];
while (mag(Tsi - TsiNew)/Tsi > 0.01)
{
Tsi = TsiNew;
scalar hradNew
(
eSig*((pow3(Ta) + pow3(Tsi)) + Ta*Tsi*(Ta + Tsi))
);
TaLambda =
(htcCoeff[i] + hradNew)
/(htcCoeff[i] + hradNew + 1/totalSolidRes);
TsiNew = TaLambda*Ta + (1 - TaLambda)*Tp[i];
};
hrad[i] =
eSig*((pow3(Ta) + pow3(Tsi)) + Ta*Tsi*(Ta + Tsi));
}
} }
else else
{ {
// ... if there is no solid wall thermal resistance use hrad = eSig*((pow3(Ta) + pow3(Tp)) + Ta*Tp*(Ta + Tp));
// the current wall temperature
hp += emissivity_*sigma.value()*pow3(Tp);
hpTa += emissivity_*sigma.value()*pow4(Ta);
} }
} }
const scalarField hp(1/(1/(htcCoeff + hrad) + totalSolidRes));
const scalarField hpTa(hp*Ta);
const scalarField kappaDeltaCoeffs const scalarField kappaDeltaCoeffs
( (
this->kappa(Tp)*patch().deltaCoeffs() this->kappa(Tp)*patch().deltaCoeffs()