Commit d6b404db authored by Henry Weller's avatar Henry Weller
Browse files

reactingTwoPhaseEulerFoam::IATE: Added wallBoiling sub-model

to handle the size of bubbles created by boiling.  To be used in
conjunction with the alphatWallBoilingWallFunction boundary condition.

The IATE variant of the wallBoiling tutorial case is provided to
demonstrate the functionality:

tutorials/multiphase/reactingTwoPhaseEulerFoam/RAS/wallBoilingIATE
parent 0d66ffcc
......@@ -56,6 +56,8 @@ derivedFvPatchFields/wallBoilingSubModels/departureFrequencyModels/departureFreq
derivedFvPatchFields/wallBoilingSubModels/departureFrequencyModels/departureFrequencyModel/newDepartureFrequencyModel.C
derivedFvPatchFields/wallBoilingSubModels/departureFrequencyModels/Cole/Cole.C
derivedFvPatchFields/wallBoilingSubModels/IATEsource/wallBoiling.C
derivedFvPatchFields/alphatPhaseChangeJayatillekeWallFunction/alphatPhaseChangeJayatillekeWallFunctionFvPatchScalarField.C
derivedFvPatchFields/alphatFixedDmdtWallBoilingWallFunction/alphatFixedDmdtWallBoilingWallFunctionFvPatchScalarField.C
derivedFvPatchFields/alphatWallBoilingWallFunction/alphatWallBoilingWallFunctionFvPatchScalarField.C
......
......@@ -84,6 +84,7 @@ alphatWallBoilingWallFunctionFvPatchScalarField
relax_(0.5),
AbyV_(p.size(), 0),
alphatConv_(p.size(), 0),
dDep_(p.size(), 0),
partitioningModel_(nullptr),
nucleationSiteModel_(nullptr),
departureDiamModel_(nullptr),
......@@ -92,7 +93,7 @@ alphatWallBoilingWallFunctionFvPatchScalarField
AbyV_ = this->patch().magSf();
forAll(AbyV_, facei)
{
label faceCelli = this->patch().faceCells()[facei];
const label faceCelli = this->patch().faceCells()[facei];
AbyV_[facei] /= iF.mesh().V()[faceCelli];
}
}
......@@ -111,6 +112,7 @@ alphatWallBoilingWallFunctionFvPatchScalarField
relax_(dict.lookupOrDefault<scalar>("relax", 0.5)),
AbyV_(p.size(), 0),
alphatConv_(p.size(), 0),
dDep_(p.size(), 0),
partitioningModel_(nullptr),
nucleationSiteModel_(nullptr),
departureDiamModel_(nullptr),
......@@ -192,6 +194,7 @@ alphatWallBoilingWallFunctionFvPatchScalarField
relax_(psf.relax_),
AbyV_(psf.AbyV_),
alphatConv_(psf.alphatConv_, mapper),
dDep_(psf.dDep_, mapper),
partitioningModel_(psf.partitioningModel_),
nucleationSiteModel_(psf.nucleationSiteModel_),
departureDiamModel_(psf.departureDiamModel_),
......@@ -209,6 +212,7 @@ alphatWallBoilingWallFunctionFvPatchScalarField
relax_(psf.relax_),
AbyV_(psf.AbyV_),
alphatConv_(psf.alphatConv_),
dDep_(psf.dDep_),
partitioningModel_(psf.partitioningModel_),
nucleationSiteModel_(psf.nucleationSiteModel_),
departureDiamModel_(psf.departureDiamModel_),
......@@ -227,6 +231,7 @@ alphatWallBoilingWallFunctionFvPatchScalarField
relax_(psf.relax_),
AbyV_(psf.AbyV_),
alphatConv_(psf.alphatConv_),
dDep_(psf.dDep_),
partitioningModel_(psf.partitioningModel_),
nucleationSiteModel_(psf.nucleationSiteModel_),
departureDiamModel_(psf.departureDiamModel_),
......@@ -444,15 +449,12 @@ void alphatWallBoilingWallFunctionFvPatchScalarField::updateCoeffs()
);
// Bubble departure diameter:
const scalarField dDep
dDep_ = departureDiamModel_->dDeparture
(
departureDiamModel_->dDeparture
(
liquid,
vapor,
patchi,
Tsub
)
liquid,
vapor,
patchi,
Tsub
);
// Bubble departure frequency:
......@@ -463,7 +465,7 @@ void alphatWallBoilingWallFunctionFvPatchScalarField::updateCoeffs()
liquid,
vapor,
patchi,
dDep
dDep_
)
);
......@@ -473,12 +475,12 @@ void alphatWallBoilingWallFunctionFvPatchScalarField::updateCoeffs()
const scalarField Ja(rhoLiquidw*Cpw*Tsub/(rhoVaporw*L));
const scalarField Al(fLiquid*4.8*exp(-Ja/80));
const scalarField A2(min(M_PI*sqr(dDep)*N*Al/4, scalar(1)));
const scalarField A2(min(M_PI*sqr(dDep_)*N*Al/4, scalar(1)));
const scalarField A1(max(1 - A2, scalar(1e-4)));
const scalarField A2E(min(M_PI*sqr(dDep)*N*Al/4, scalar(5)));
const scalarField A2E(min(M_PI*sqr(dDep_)*N*Al/4, scalar(5)));
// Wall evaporation heat flux [kg/s3 = J/m2s]
const scalarField Qe((1.0/6.0)*A2E*dDep*rhoVaporw*fDep*L);
const scalarField Qe((1.0/6.0)*A2E*dDep_*rhoVaporw*fDep*L);
// Volumetric mass source in the near wall cell due to the
// wall boiling
......@@ -540,8 +542,8 @@ void alphatWallBoilingWallFunctionFvPatchScalarField::updateCoeffs()
Info<< " L: " << gMin(L) << " - " << gMax(L) << endl;
Info<< " Tl: " << gMin(Tl) << " - " << gMax(Tl) << endl;
Info<< " N: " << gMin(N) << " - " << gMax(N) << endl;
Info<< " dDep: " << gMin(dDep) << " - "
<< gMax(dDep) << endl;
Info<< " dDep_: " << gMin(dDep_) << " - "
<< gMax(dDep_) << endl;
Info<< " fDep: " << gMin(fDep) << " - "
<< gMax(fDep) << endl;
Info<< " Al: " << gMin(Al) << " - " << gMax(Al) << endl;
......
......@@ -178,6 +178,9 @@ private:
//- Convective turbulent thermal diffusivity
scalarField alphatConv_;
//- Departure diameter field
scalarField dDep_;
//- Run-time selected heat flux partitioning model
autoPtr<wallBoilingModels::partitioningModel>
partitioningModel_;
......@@ -266,6 +269,13 @@ public:
// Member functions
//- Calculate and return the departure diameter field
const scalarField& dDeparture() const
{
return dDep_;
}
// Evaluation functions
//- Update the coefficients associated with the patch field
......
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "wallBoiling.H"
#include "phaseCompressibleTurbulenceModel.H"
#include "alphatWallBoilingWallFunctionFvPatchScalarField.H"
#include "fvmSup.H"
#include "addToRunTimeSelectionTable.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
namespace diameterModels
{
namespace IATEsources
{
defineTypeNameAndDebug(wallBoiling, 0);
addToRunTimeSelectionTable(IATEsource, wallBoiling, dictionary);
}
}
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::diameterModels::IATEsources::wallBoiling::wallBoiling
(
const IATE& iate,
const dictionary& dict
)
:
IATEsource(iate)
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::tmp<Foam::fvScalarMatrix>
Foam::diameterModels::IATEsources::wallBoiling::R
(
volScalarField& kappai
) const
{
volScalarField::Internal R
(
IOobject
(
"wallBoiling:R",
phase().time().timeName(),
phase().mesh()
),
phase().mesh(),
dimensionedScalar("R", dimless/dimTime, 0)
);
volScalarField::Internal Rdk
(
IOobject
(
"wallBoiling:Rdk",
phase().time().timeName(),
phase().mesh()
),
phase().mesh(),
dimensionedScalar("Rdk", kappai.dimensions()/dimTime, 0)
);
const phaseCompressibleTurbulenceModel& turbulence =
const_cast<phaseCompressibleTurbulenceModel&>
(
phase().db().lookupObject<phaseCompressibleTurbulenceModel>
(
IOobject::groupName
(
turbulenceModel::propertiesName,
otherPhase().name()
)
)
);
const tmp<volScalarField> talphat(turbulence.alphat());
const volScalarField::Boundary& alphatBf = talphat().boundaryField();
const scalarField& rho = phase().rho();
typedef compressible::alphatWallBoilingWallFunctionFvPatchScalarField
alphatWallBoilingWallFunction;
forAll(alphatBf, patchi)
{
if
(
isA<alphatWallBoilingWallFunction>(alphatBf[patchi])
)
{
const alphatWallBoilingWallFunction& alphatw =
refCast<const alphatWallBoilingWallFunction>(alphatBf[patchi]);
const scalarField& dmdt = alphatw.dmdt();
const scalarField& dDep = alphatw.dDeparture();
const labelList& faceCells = alphatw.patch().faceCells();
forAll(alphatw, facei)
{
if (dmdt[facei] > SMALL)
{
const label faceCelli = faceCells[facei];
R[faceCelli] = (dmdt[facei]/rho[faceCelli]);
Rdk[faceCelli] = R[faceCelli]*(6.0/dDep[facei]);
}
}
}
}
return Rdk - fvm::Sp(R, kappai);
}
// ************************************************************************* //
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Class
Foam::diameterModels::IATEsources::wallBoiling
Description
Wall-boiling IATE source.
SourceFiles
wallBoiling.C
\*---------------------------------------------------------------------------*/
#ifndef wallBoiling_H
#define wallBoiling_H
#include "IATEsource.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace diameterModels
{
namespace IATEsources
{
/*---------------------------------------------------------------------------*\
Class wallBoiling Declaration
\*---------------------------------------------------------------------------*/
class wallBoiling
:
public IATEsource
{
public:
//- Runtime type information
TypeName("wallBoiling");
// Constructors
wallBoiling
(
const IATE& iate,
const dictionary& dict
);
//- Destructor
virtual ~wallBoiling()
{}
// Member Functions
virtual tmp<fvScalarMatrix> R(volScalarField& kappai) const;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace IATEsources
} // End namespace diameterModels
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //
......@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2013-2015 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2013-2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
......@@ -116,43 +116,33 @@ Foam::tmp<Foam::volScalarField> Foam::diameterModels::IATE::dsm() const
return max(6/max(kappai_, 6/dMax_), dMin_);
}
// Placeholder for the nucleation/condensation model
// Foam::tmp<Foam::volScalarField> Foam::diameterModels::IATE::Rph() const
// {
// const volScalarField& T = phase_thermo().T();
// const volScalarField& p = phase_.p();
//
// scalar A, B, C, sigma, vm, Rph;
//
// volScalarField ps(1e5*pow(10, A - B/(T + C)));
// volScalarField Dbc
// (
// 4*sigma*vm/(constant::physicoChemical::k*T*log(p/ps))
// );
//
// return constant::mathematical::pi*sqr(Dbc)*Rph;
// }
void Foam::diameterModels::IATE::correct()
{
// Initialise the accumulated source term to the dilatation effect
volScalarField R
fvScalarMatrix R
(
-fvm::SuSp
(
(1.0/3.0)
/max
(
fvc::average(phase_ + phase_.oldTime()),
residualAlpha_
(1.0/3.0)
/max
(
fvc::average(phase_ + phase_.oldTime()),
residualAlpha_
)
)
*(
fvc::ddt(phase_) + fvc::div(phase_.alphaPhi())
- phase_.continuityError()/phase_.rho()
),
kappai_
)
*(fvc::ddt(phase_) + fvc::div(phase_.alphaPhi()))
);
// Accumulate the run-time selectable sources
forAll(sources_, j)
{
R -= sources_[j].R();
R += sources_[j].R(kappai_);
}
fv::options& fvOptions(fv::options::New(phase_.mesh()));
......@@ -163,8 +153,7 @@ void Foam::diameterModels::IATE::correct()
fvm::ddt(kappai_) + fvm::div(phase_.phi(), kappai_)
- fvm::Sp(fvc::div(phase_.phi()), kappai_)
==
- fvm::SuSp(R, kappai_)
//+ Rph() // Omit the nucleation/condensation term
R
+ fvOptions(kappai_)
);
......
......@@ -177,7 +177,7 @@ public:
//- Return the bubble Webber number
tmp<volScalarField> We() const;
virtual tmp<volScalarField> R() const = 0;
virtual tmp<fvScalarMatrix> R(volScalarField& kappai) const = 0;
};
......
......@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2013-2015 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2013-2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
......@@ -24,6 +24,7 @@ License
\*---------------------------------------------------------------------------*/
#include "randomCoalescence.H"
#include "fvmSup.H"
#include "addToRunTimeSelectionTable.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
......@@ -59,42 +60,39 @@ randomCoalescence
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::tmp<Foam::volScalarField>
Foam::diameterModels::IATEsources::randomCoalescence::R() const
Foam::tmp<Foam::fvScalarMatrix>
Foam::diameterModels::IATEsources::randomCoalescence::R
(
volScalarField& kappai
) const
{
tmp<volScalarField> tR
volScalarField::Internal R
(
new volScalarField
IOobject
(
IOobject
(
"R",
iate_.phase().time().timeName(),
iate_.phase().mesh()
),
iate_.phase().mesh(),
dimensionedScalar("R", dimless/dimTime, 0)
)
"randomCoalescence:R",
iate_.phase().time().timeName(),
iate_.phase().mesh()
),
iate_.phase().mesh(),
dimensionedScalar("R", dimless/dimTime, 0)
);
volScalarField R = tR();
scalar Crc = Crc_.value();
scalar C = C_.value();
scalar alphaMax = alphaMax_.value();
volScalarField Ut(this->Ut());
const scalar Crc = Crc_.value();
const scalar C = C_.value();
const scalar alphaMax = alphaMax_.value();
const volScalarField Ut(this->Ut());