Skip to content
GitLab
Commit f83bfb62 authored by Andrew Heather's avatar Andrew Heather
Browse files

ENH: Updates to surface film models

parent 97ee7cfb
Hide whitespace changes
Inline Side-by-side
src/surfaceFilmModels/submodels/thermo/phaseChangeModel/standardPhaseChange/standardPhaseChange.C
View file @ f83bfb62
......@@ -74,6 +74,7 @@ Foam::surfaceFilmModels::standardPhaseChange::standardPhaseChange
)
:
phaseChangeModel(typeName, owner, dict),
Tb_(readScalar(coeffs_.lookup("Tb"))),
deltaMin_(readScalar(coeffs_.lookup("deltaMin"))),
L_(readScalar(coeffs_.lookup("L")))
{}
......@@ -116,48 +117,63 @@ void Foam::surfaceFilmModels::standardPhaseChange::correct
const scalarField hInf = film.htcs().h();
const scalarField hFilm = film.htcw().h();
const vectorField dU = film.UPrimary() - film.Us();
const scalarField availableMass = delta*rho*magSf;
const scalarField availableMass = (delta - deltaMin_)*rho*magSf;
// Reynolds number
const scalarField Re = rho*mag(dU)*L_/mu;
// molecular weight of vapour
const scalar Wvap = thermo.carrier().W(vapId);
// molecular weight of liquid
const scalar Wliq = liq.W();
forAll(dMass, cellI)
{
if (delta[cellI] > deltaMin_)
{
// cell pressure
// cell pressure [Pa]
const scalar pc = pInf[cellI];
// saturation pressure
const scalar pSat = liq.pv(pc, Ts[cellI]);
// local temperature - impose lower limit of 200 K for stability
const scalar Tloc = min(liq.Tc(), max(200.0, Ts[cellI]));
// saturation pressure [Pa]
const scalar pSat = liq.pv(pc, Tloc);
// latent heat
const scalar hVap = liq.hl(pc, Ts[cellI]);
// latent heat [J/kg]
const scalar hVap = liq.hl(pc, Tloc);
// calculate mass transfer
if (pSat > pc)
if (pSat >= pc)
{
// boiling
const scalar qDotInf = hInf[cellI]*(TInf[cellI] - T[cellI]);
const scalar qDotFilm = hFilm[cellI]*(T[cellI] - Tw[cellI]);
dMass[cellI] = dt*magSf[cellI]/hVap*(qDotInf - qDotFilm);
const scalar cp = liq.cp(pc, Tloc);
const scalar qCorr = availableMass[cellI]*cp*(T[cellI] - Tb_);
dMass[cellI] =
dt*magSf[cellI]/hVap*(qDotInf + qDotFilm)
+ qCorr/hVap;
}
else
{
// Reynolds number
const scalarField Re = rho*mag(dU)*L_/mu;
// molecular weight of vapour [kg/kmol]
const scalar Wvap = thermo.carrier().W(vapId);
// molecular weight of liquid [kg/kmol]
const scalar Wliq = liq.W();
// vapour mass fraction at interface
const scalar Ys = Wliq*pSat/(Wliq*pSat + Wvap*(pc - pSat));
// bulk gas average density
const scalar rhoAve = pc/(specie::RR*Ts[cellI]);
// bulk gas average density [kg/m3]
scalar Winf = 0;
forAll(thermo.carrier().Y(), i)
{
Winf += film.YPrimary()[i][cellI]*thermo.carrier().W(i);
}
const scalar rhoInf = Winf*pc/(specie::RR*Tloc);
// vapour diffusivity [m2/s]
const scalar Dab = liq.D(pc, Ts[cellI]);
const scalar Dab = liq.D(pc, Tloc);
// Schmidt number
const scalar Sc = mu[cellI]/(rho[cellI]*(Dab + ROOTVSMALL));
......@@ -170,7 +186,7 @@ void Foam::surfaceFilmModels::standardPhaseChange::correct
// add mass contribution to source
dMass[cellI] =
dt*magSf[cellI]*rhoAve*hm*(Ys - YInf[cellI])/(1.0 - Ys);
dt*magSf[cellI]*rhoInf*hm*(Ys - YInf[cellI])/(1.0 - Ys);
}
dMass[cellI] = min(availableMass[cellI], max(0.0, dMass[cellI]));
......
src/surfaceFilmModels/submodels/thermo/phaseChangeModel/standardPhaseChange/standardPhaseChange.H
View file @ f83bfb62
......@@ -67,10 +67,13 @@ protected:
// Protected data
//- Boiling temperature / [K]
const scalar Tb_;
//- Minimum film height for model to be active
const scalar deltaMin_;
//- Length scalae / [m]
//- Length scale / [m]
const scalar L_;
......
src/surfaceFilmModels/surfaceFilmModel/kinematicSingleLayer/kinematicSingleLayer.C
View file @ f83bfb62
......@@ -222,6 +222,10 @@ transferPrimaryRegionFields()
rhoSp_.field() /= magSf_*deltaT;
USp_.field() /= magSf_*deltaT;
pSp_.field() /= magSf_*deltaT;
// reset transfer to primary fields
massForPrimary_ == dimensionedScalar("zero", dimMass, 0.0);
diametersForPrimary_ == dimensionedScalar("zero", dimLength, -1.0);
}
......@@ -271,9 +275,6 @@ Foam::surfaceFilmModels::kinematicSingleLayer::pp()
void Foam::surfaceFilmModels::kinematicSingleLayer::correctDetachedFilm()
{
massForPrimary_ == dimensionedScalar("zero", dimMass, 0.0);
diametersForPrimary_ == dimensionedScalar("zero", dimLength, -1.0);
const scalarField gNorm = this->gNorm();
forAll(gNorm, i)
......@@ -308,8 +309,7 @@ void Foam::surfaceFilmModels::kinematicSingleLayer::updateSubmodels()
// Update source fields
const dimensionedScalar deltaT = time_.deltaT();
rhoSp_ -= massForPrimary_/magSf_/deltaT;
USp_ -= massForPrimary_*U_/magSf_/deltaT;
rhoSp_ -= (massForPrimary_ + massPhaseChangeForPrimary_)/magSf_/deltaT;
}
......@@ -457,6 +457,12 @@ Foam::surfaceFilmModels::kinematicSingleLayer::solveMomentum
USp_
+ tau(U_)
+ fvc::grad(sigma_)
- fvm::Sp
(
(massForPrimary_ + massPhaseChangeForPrimary_)
/magSf_/time_.deltaT(),
U_
)
);
fvVectorMatrix& UEqn = tUEqn();
......@@ -1239,7 +1245,7 @@ Foam::surfaceFilmModels::kinematicSingleLayer::Srho() const
const label filmPatchI = filmBottomPatchIDs_[i];
scalarField patchMass =
massPhaseChangeForPrimary().boundaryField()[filmPatchI];
massPhaseChangeForPrimary_.boundaryField()[filmPatchI];
distMap.distribute(patchMass);
......
src/surfaceFilmModels/surfaceFilmModel/noFilm/noFilm.C
View file @ f83bfb62
......@@ -138,10 +138,8 @@ void Foam::surfaceFilmModels::noFilm::addSources
const Foam::volVectorField& Foam::surfaceFilmModels::noFilm::U() const
{
FatalErrorIn
(
"const volScalarField& noFilm::U() const"
) << "U field not available for " << type() << abort(FatalError);
FatalErrorIn("const volScalarField& noFilm::U() const")
<< "U field not available for " << type() << abort(FatalError);
return volVectorField::null();
}
......@@ -149,10 +147,8 @@ const Foam::volVectorField& Foam::surfaceFilmModels::noFilm::U() const
const Foam::volVectorField& Foam::surfaceFilmModels::noFilm::Us() const
{
FatalErrorIn
(
"const volScalarField& noFilm::Us() const"
) << "Us field not available for " << type() << abort(FatalError);
FatalErrorIn("const volScalarField& noFilm::Us() const")
<< "Us field not available for " << type() << abort(FatalError);
return volVectorField::null();
}
......@@ -160,10 +156,8 @@ const Foam::volVectorField& Foam::surfaceFilmModels::noFilm::Us() const
const Foam::volVectorField& Foam::surfaceFilmModels::noFilm::Uw() const
{
FatalErrorIn
(
"const volScalarField& noFilm::Uw() const"
) << "Uw field not available for " << type() << abort(FatalError);
FatalErrorIn("const volScalarField& noFilm::Uw() const")
<< "Uw field not available for " << type() << abort(FatalError);
return volVectorField::null();
}
......@@ -171,10 +165,8 @@ const Foam::volVectorField& Foam::surfaceFilmModels::noFilm::Uw() const
const Foam::volScalarField& Foam::surfaceFilmModels::noFilm::rho() const
{
FatalErrorIn
(
"const volScalarField& noFilm::rho() const"
) << "rho field not available for " << type() << abort(FatalError);
FatalErrorIn("const volScalarField& noFilm::rho() const")
<< "rho field not available for " << type() << abort(FatalError);
return volScalarField::null();
}
......@@ -182,10 +174,8 @@ const Foam::volScalarField& Foam::surfaceFilmModels::noFilm::rho() const
const Foam::volScalarField& Foam::surfaceFilmModels::noFilm::T() const
{
FatalErrorIn
(
"const Foam::volScalarField& Foam::noFilm::T() const"
) << "T field not available for " << type() << abort(FatalError);
FatalErrorIn("const Foam::volScalarField& Foam::noFilm::T() const")
<< "T field not available for " << type() << abort(FatalError);
return volScalarField::null();
}
......@@ -193,10 +183,8 @@ const Foam::volScalarField& Foam::surfaceFilmModels::noFilm::T() const
const Foam::volScalarField& Foam::surfaceFilmModels::noFilm::Ts() const
{
FatalErrorIn
(
"const Foam::volScalarField& Foam::noFilm::Ts() const"
) << "Ts field not available for " << type() << abort(FatalError);
FatalErrorIn("const Foam::volScalarField& Foam::noFilm::Ts() const")
<< "Ts field not available for " << type() << abort(FatalError);
return volScalarField::null();
}
......@@ -204,10 +192,8 @@ const Foam::volScalarField& Foam::surfaceFilmModels::noFilm::Ts() const
const Foam::volScalarField& Foam::surfaceFilmModels::noFilm::Tw() const
{
FatalErrorIn
(
"const Foam::volScalarField& Foam::noFilm::Tw() const"
) << "Tw field not available for " << type() << abort(FatalError);
FatalErrorIn("const Foam::volScalarField& Foam::noFilm::Tw() const")
<< "Tw field not available for " << type() << abort(FatalError);
return volScalarField::null();
}
......@@ -215,10 +201,8 @@ const Foam::volScalarField& Foam::surfaceFilmModels::noFilm::Tw() const
const Foam::volScalarField& Foam::surfaceFilmModels::noFilm::cp() const
{
FatalErrorIn
(
"const volScalarField& noFilm::cp() const"
) << "cp field not available for " << type() << abort(FatalError);
FatalErrorIn("const volScalarField& noFilm::cp() const")
<< "cp field not available for " << type() << abort(FatalError);
return volScalarField::null();
}
......@@ -226,10 +210,8 @@ const Foam::volScalarField& Foam::surfaceFilmModels::noFilm::cp() const
const Foam::volScalarField& Foam::surfaceFilmModels::noFilm::kappa() const
{
FatalErrorIn
(
"const volScalarField& noFilm::kappa() const"
) << "kappa field not available for " << type() << abort(FatalError);
FatalErrorIn("const volScalarField& noFilm::kappa() const")
<< "kappa field not available for " << type() << abort(FatalError);
return volScalarField::null();
}
......@@ -238,10 +220,8 @@ const Foam::volScalarField& Foam::surfaceFilmModels::noFilm::kappa() const
const Foam::volScalarField&
Foam::surfaceFilmModels::noFilm::massForPrimary() const
{
FatalErrorIn
(
"const volScalarField& noFilm::massForPrimary() const"
) << "massForPrimary field not available for " << type()
FatalErrorIn("const volScalarField& noFilm::massForPrimary() const")
<< "massForPrimary field not available for " << type()
<< abort(FatalError);
return volScalarField::null();
......@@ -251,10 +231,8 @@ Foam::surfaceFilmModels::noFilm::massForPrimary() const
const Foam::volScalarField&
Foam::surfaceFilmModels::noFilm::diametersForPrimary() const
{
FatalErrorIn
(
"const volScalarField& noFilm::diametersForPrimary() const"
) << "diametersForPrimary field not available for " << type()
FatalErrorIn("const volScalarField& noFilm::diametersForPrimary() const")
<< "diametersForPrimary field not available for " << type()
<< abort(FatalError);
return volScalarField::null();
......@@ -267,7 +245,7 @@ Foam::surfaceFilmModels::noFilm::massPhaseChangeForPrimary() const
FatalErrorIn
(
"const volScalarField& noFilm::massPhaseChangeForPrimary() const"
) << "massPhaseChangeForPrimary field not available for " << type()
) << "massPhaseChange field not available for " << type()
<< abort(FatalError);
return volScalarField::null();
......
src/surfaceFilmModels/surfaceFilmModel/thermoSingleLayer/thermoSingleLayer.C
View file @ f83bfb62
......@@ -175,8 +175,8 @@ void Foam::surfaceFilmModels::thermoSingleLayer::updateSubmodels()
htcw_->correct();
// Update phase change
massPhaseChangeForPrimary_ == dimensionedScalar("zero", dimMass, 0.0);
energyPhaseChangeForPrimary_ == dimensionedScalar("zero", dimEnergy, 0.0);
massPhaseChangeForPrimary_.internalField() = 0.0;
energyPhaseChangeForPrimary_.internalField() = 0.0;
phaseChange_->correct
(
......@@ -191,11 +191,7 @@ void Foam::surfaceFilmModels::thermoSingleLayer::updateSubmodels()
kinematicSingleLayer::updateSubmodels();
// Update source fields
const dimensionedScalar deltaT = time_.deltaT();
rhoSp_ -= massPhaseChangeForPrimary_/magSf_/deltaT;
USp_ -= massPhaseChangeForPrimary_*U_/magSf_/deltaT;
hsSp_ -=
(massForPrimary_*hs_ + energyPhaseChangeForPrimary_)/magSf_/deltaT;
hsSp_ -= energyPhaseChangeForPrimary_/magSf_/time_.deltaT();
}
......@@ -208,7 +204,7 @@ Foam::tmp<Foam::fvScalarMatrix> Foam::surfaceFilmModels::thermoSingleLayer::q
return
(
- fvm::Sp(htcs_->h()/cp_, hs) - htcs_->h()*(Tstd - Ts_)
- fvm::Sp(htcs_->h()/cp_, hs) - htcs_->h()*(Tstd - TPrimary_)
- fvm::Sp(htcw_->h()/cp_, hs) - htcw_->h()*(Tstd - Tw_)
);
}
......@@ -228,8 +224,9 @@ void Foam::surfaceFilmModels::thermoSingleLayer::solveEnergy()
fvm::ddt(deltaRho_, hs_)
+ fvm::div(phi_, hs_)
==
hsSp_
fvm::Sp(hsSp_/hs_, hs_)
+ q(hs_)
- fvm::Sp(massForPrimary_/magSf_/time_.deltaT(), hs_)
);
correctThermoFields();
......@@ -627,7 +624,7 @@ Foam::surfaceFilmModels::thermoSingleLayer::Srho(const label i) const
const directMappedWallPolyPatch& wpp =
refCast<const directMappedWallPolyPatch>
(
mesh_.boundaryMesh()[primaryPatchI]
mesh_.boundaryMesh()[primaryPatchI]
);
const mapDistribute& distMap = wpp.map();
......@@ -635,7 +632,7 @@ Foam::surfaceFilmModels::thermoSingleLayer::Srho(const label i) const
const label filmPatchI = filmBottomPatchIDs_[i];
scalarField patchMass =
massPhaseChangeForPrimary().boundaryField()[filmPatchI];
massPhaseChangeForPrimary_.boundaryField()[filmPatchI];
distMap.distribute(patchMass);
......@@ -672,7 +669,7 @@ Foam::surfaceFilmModels::thermoSingleLayer::Sh() const
dimensionedScalar("zero", dimEnergy/dimVolume/dimTime, 0.0)
)
);
/*
scalarField& Sh = tSh();
const scalarField& V = mesh_.V();
const scalar dt = time_.deltaTValue();
......@@ -683,28 +680,25 @@ Foam::surfaceFilmModels::thermoSingleLayer::Sh() const
const directMappedWallPolyPatch& wpp =
refCast<const directMappedWallPolyPatch>
(
mesh_.boundaryMesh()[primaryPatchI]
mesh_.boundaryMesh()[primaryPatchI]
);
const mapDistribute& distMap = wpp.map();
const label filmPatchI = filmBottomPatchIDs_[i];
scalarField patchMass =
massPhaseChangeForPrimary().boundaryField()[filmPatchI];
distMap.distribute(patchMass);
scalarField patchEnthalpy = hs_.boundaryField()[filmPatchI];
distMap.distribute(patchEnthalpy);
scalarField patchEnergy =
energyPhaseChangeForPrimary_.boundaryField()[filmPatchI];
distMap.distribute(patchEnergy);
const unallocLabelList& cells = wpp.faceCells();
forAll(patchMass, j)
{
Sh[cells[j]] += patchMass[j]*patchEnthalpy[j]/(V[cells[j]]*dt);
Sh[cells[j]] += patchEnergy[j]/(V[cells[j]]*dt);
}
}
*/
return tSh;
}
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment