Commit cf25c464 authored by Andrew Heather's avatar Andrew Heather
Browse files

added RAS variant of the Jayatilleke thermal wall function

parent 790e84a2
......@@ -15,6 +15,7 @@ wallFunctions = derivedFvPatchFields/wallFunctions
alphatWallFunctions = $(wallFunctions)/alphatWallFunctions
$(alphatWallFunctions)/alphatWallFunction/alphatWallFunctionFvPatchScalarField.C
$(alphatWallFunctions)/alphatJayatillekeWallFunction/alphatJayatillekeWallFunctionFvPatchScalarField.C
mutWallFunctions = $(wallFunctions)/mutWallFunctions
$(mutWallFunctions)/mutkRoughWallFunction/mutkRoughWallFunctionFvPatchScalarField.C
......
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2008-2009 OpenCFD Ltd.
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
\*---------------------------------------------------------------------------*/
#include "alphatJayatillekeWallFunctionFvPatchScalarField.H"
#include "RASModel.H"
#include "fvPatchFieldMapper.H"
#include "volFields.H"
#include "addToRunTimeSelectionTable.H"
#include "wallFvPatch.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace compressible
{
namespace RASModels
{
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
scalar alphatJayatillekeWallFunctionFvPatchScalarField::maxExp_ = 50.0;
scalar alphatJayatillekeWallFunctionFvPatchScalarField::tolerance_ = 0.01;
label alphatJayatillekeWallFunctionFvPatchScalarField::maxIters_ = 10;
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
void alphatJayatillekeWallFunctionFvPatchScalarField::checkType()
{
if (!isA<wallFvPatch>(patch()))
{
FatalErrorIn
(
"alphatJayatillekeWallFunctionFvPatchScalarField::checkType()"
)
<< "Patch type for patch " << patch().name() << " must be wall\n"
<< "Current patch type is " << patch().type() << nl
<< exit(FatalError);
}
}
scalar alphatJayatillekeWallFunctionFvPatchScalarField::Psmooth
(
const scalar Prat
) const
{
return 9.24*(pow(Prat, 0.75) - 1.0)*(1.0 + 0.28*exp(-0.007*Prat));
}
scalar alphatJayatillekeWallFunctionFvPatchScalarField::yPlusTherm
(
const scalar P,
const scalar Prat
) const
{
scalar ypt = 11.0;
for (int i=0; i<maxIters_; i++)
{
scalar f = ypt - (log(E_*ypt)/kappa_ + P)/Prat;
scalar df = 1.0 - 1.0/(ypt*kappa_*Prat);
scalar yptNew = ypt - f/df;
if (yptNew < VSMALL)
{
return 0;
}
else if (mag(yptNew - ypt) < tolerance_)
{
return yptNew;
}
else
{
ypt = yptNew;
}
}
return ypt;
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
alphatJayatillekeWallFunctionFvPatchScalarField::
alphatJayatillekeWallFunctionFvPatchScalarField
(
const fvPatch& p,
const DimensionedField<scalar, volMesh>& iF
)
:
fixedValueFvPatchScalarField(p, iF),
Prt_(0.85),
Cmu_(0.09),
kappa_(0.41),
E_(9.8)
{
checkType();
}
alphatJayatillekeWallFunctionFvPatchScalarField::
alphatJayatillekeWallFunctionFvPatchScalarField
(
const alphatJayatillekeWallFunctionFvPatchScalarField& ptf,
const fvPatch& p,
const DimensionedField<scalar, volMesh>& iF,
const fvPatchFieldMapper& mapper
)
:
fixedValueFvPatchScalarField(ptf, p, iF, mapper),
Prt_(ptf.Prt_),
Cmu_(ptf.Cmu_),
kappa_(ptf.kappa_),
E_(ptf.E_)
{}
alphatJayatillekeWallFunctionFvPatchScalarField::
alphatJayatillekeWallFunctionFvPatchScalarField
(
const fvPatch& p,
const DimensionedField<scalar, volMesh>& iF,
const dictionary& dict
)
:
fixedValueFvPatchScalarField(p, iF, dict),
Prt_(dict.lookupOrDefault<scalar>("Prt", 0.85)),
Cmu_(dict.lookupOrDefault<scalar>("Cmu", 0.09)),
kappa_(dict.lookupOrDefault<scalar>("kappa", 0.41)),
E_(dict.lookupOrDefault<scalar>("E", 9.8))
{
checkType();
}
alphatJayatillekeWallFunctionFvPatchScalarField::
alphatJayatillekeWallFunctionFvPatchScalarField
(
const alphatJayatillekeWallFunctionFvPatchScalarField& awfpsf
)
:
fixedValueFvPatchScalarField(awfpsf),
Prt_(awfpsf.Prt_),
Cmu_(awfpsf.Cmu_),
kappa_(awfpsf.kappa_),
E_(awfpsf.E_)
{
checkType();
}
alphatJayatillekeWallFunctionFvPatchScalarField::
alphatJayatillekeWallFunctionFvPatchScalarField
(
const alphatJayatillekeWallFunctionFvPatchScalarField& awfpsf,
const DimensionedField<scalar, volMesh>& iF
)
:
fixedValueFvPatchScalarField(awfpsf, iF),
Prt_(awfpsf.Prt_),
Cmu_(awfpsf.Cmu_),
kappa_(awfpsf.kappa_),
E_(awfpsf.E_)
{
checkType();
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
void alphatJayatillekeWallFunctionFvPatchScalarField::updateCoeffs()
{
if (updated())
{
return;
}
const label patchI = patch().index();
const RASModel& rasModel = db().lookupObject<RASModel>("RASProperties");
const scalar Cmu25 = pow(Cmu_, 0.25);
const scalarField& y = rasModel.y()[patchI];
const scalarField& muw = rasModel.mu().boundaryField()[patchI];
const scalarField& alphaw = rasModel.alpha().boundaryField()[patchI];
scalarField& alphatw = *this;
const tmp<volScalarField> tk = rasModel.k();
const volScalarField& k = tk();
const fvPatchVectorField& Uw = rasModel.U().boundaryField()[patchI];
const scalarField magUp = mag(Uw.patchInternalField() - Uw);
const scalarField magGradUw = mag(Uw.snGrad());
const scalarField& rhow = rasModel.rho().boundaryField()[patchI];
const fvPatchScalarField& hw =
patch().lookupPatchField<volScalarField, scalar>("h");
// Heat flux [W/m2] - lagging alphatw
const scalarField qDot = (alphaw + alphatw)*hw.snGrad();
// Populate boundary values
forAll(alphatw, faceI)
{
label faceCellI = patch().faceCells()[faceI];
scalar uTau = Cmu25*sqrt(k[faceCellI]);
scalar yPlus = uTau*y[faceI]/(muw[faceI]/rhow[faceI]);
// Molecular Prandtl number
scalar Pr = muw[faceI]/alphaw[faceI];
// Molecular-to-turbulenbt Prandtl number ratio
scalar Prat = Pr/Prt_;
// Thermal sublayer thickness
scalar P = Psmooth(Prat);
scalar yPlusTherm = this->yPlusTherm(P, Prat);
// Evaluate new effective thermal diffusivity
scalar alphaEff = 0.0;
if (yPlus < yPlusTherm)
{
scalar A = qDot[faceI]*rhow[faceI]*uTau*y[faceI];
scalar B = qDot[faceI]*Pr*yPlus;
scalar C = Pr*0.5*rhow[faceI]*uTau*sqr(magUp[faceI]);
alphaEff = A/(B + C + VSMALL);
}
else
{
scalar A = qDot[faceI]*rhow[faceI]*uTau*y[faceI];
scalar B = qDot[faceI]*Prt_*(1.0/kappa_*log(E_*yPlus) + P);
scalar magUc = uTau/kappa_*log(E_*yPlusTherm) - mag(Uw[faceI]);
scalar C =
0.5*rhow[faceI]*uTau
*(Prt_*sqr(magUp[faceI]) + (Pr - Prt_)*sqr(magUc));
alphaEff = A/(B + C + VSMALL);
}
// Update turbulent thermal diffusivity
alphatw[faceI] = max(0.0, alphaEff - alphaw[faceI]);
if (debug)
{
Info<< " uTau = " << uTau << nl
<< " Pr = " << Pr << nl
<< " Prt = " << Prt_ << nl
<< " qDot = " << qDot[faceI] << nl
<< " yPlus = " << yPlus << nl
<< " yPlusTherm = " << yPlusTherm << nl
<< " alphaEff = " << alphaEff << nl
<< " alphaw = " << alphaw[faceI] << nl
<< " alphatw = " << alphatw[faceI] << nl
<< endl;
}
}
fixedValueFvPatchField<scalar>::updateCoeffs();
}
void alphatJayatillekeWallFunctionFvPatchScalarField::write(Ostream& os) const
{
fvPatchField<scalar>::write(os);
os.writeKeyword("Prt") << Prt_ << token::END_STATEMENT << nl;
os.writeKeyword("Cmu") << Cmu_ << token::END_STATEMENT << nl;
os.writeKeyword("kappa") << kappa_ << token::END_STATEMENT << nl;
os.writeKeyword("E") << E_ << token::END_STATEMENT << nl;
writeEntry("value", os);
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
makePatchTypeField
(
fvPatchScalarField,
alphatJayatillekeWallFunctionFvPatchScalarField
);
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace RASModels
} // End namespace compressible
} // End namespace Foam
// ************************************************************************* //
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2008-2009 OpenCFD Ltd.
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Class
alphatJayatillekeWallFunctionFvPatchScalarField
Description
Thermal wall function for turbulent thermal diffusivity based on the
Jayatilleke thermal wall function
SourceFiles
alphatJayatillekeWallFunctionFvPatchScalarField.C
\*---------------------------------------------------------------------------*/
#ifndef alphatJayatillekeWallFunctionFvPatchScalarField_H
#define alphatJayatillekeWallFunctionFvPatchScalarField_H
#include "fixedValueFvPatchFields.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace compressible
{
namespace RASModels
{
/*---------------------------------------------------------------------------*\
Class alphatJayatillekeWallFunctionFvPatchScalarField Declaration
\*---------------------------------------------------------------------------*/
class alphatJayatillekeWallFunctionFvPatchScalarField
:
public fixedValueFvPatchScalarField
{
// Private data
//- Turbulent Prandtl number
scalar Prt_;
//- Cmu coefficient
scalar Cmu_;
//- Von Karman constant
scalar kappa_;
//- E coefficient
scalar E_;
// Solution parameters
static scalar maxExp_;
static scalar tolerance_;
static label maxIters_;
// Private member functions
//- Check the type of the patch
void checkType();
//- `P' function
scalar Psmooth(const scalar Prat) const;
//- Calculate y+ at the edge of the thermal laminar sublayer
scalar yPlusTherm
(
const scalar P,
const scalar Prat
) const;
public:
//- Runtime type information
TypeName("alphatJayatillekeWallFunction");
// Constructors
//- Construct from patch and internal field
alphatJayatillekeWallFunctionFvPatchScalarField
(
const fvPatch&,
const DimensionedField<scalar, volMesh>&
);
//- Construct from patch, internal field and dictionary
alphatJayatillekeWallFunctionFvPatchScalarField
(
const fvPatch&,
const DimensionedField<scalar, volMesh>&,
const dictionary&
);
//- Construct by mapping given an
// alphatJayatillekeWallFunctionFvPatchScalarField
// onto a new patch
alphatJayatillekeWallFunctionFvPatchScalarField
(
const alphatJayatillekeWallFunctionFvPatchScalarField&,
const fvPatch&,
const DimensionedField<scalar, volMesh>&,
const fvPatchFieldMapper&
);
//- Construct as copy
alphatJayatillekeWallFunctionFvPatchScalarField
(
const alphatJayatillekeWallFunctionFvPatchScalarField&
);
//- Construct and return a clone
virtual tmp<fvPatchScalarField> clone() const
{
return tmp<fvPatchScalarField>
(
new alphatJayatillekeWallFunctionFvPatchScalarField(*this)
);
}
//- Construct as copy setting internal field reference
alphatJayatillekeWallFunctionFvPatchScalarField
(
const alphatJayatillekeWallFunctionFvPatchScalarField&,
const DimensionedField<scalar, volMesh>&
);
//- Construct and return a clone setting internal field reference
virtual tmp<fvPatchScalarField> clone
(
const DimensionedField<scalar, volMesh>& iF
) const
{
return tmp<fvPatchScalarField>
(
new alphatJayatillekeWallFunctionFvPatchScalarField
(
*this,
iF
)
);
}
// Member functions
// Evaluation functions
//- Update the coefficients associated with the patch field
virtual void updateCoeffs();
// I-O
//- Write
void write(Ostream&) const;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace RASModels
} // End namespace compressible
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //
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