00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026 #include "epsilonWallFunctionFvPatchScalarField.H"
00027 #include <compressibleRASModels/RASModel.H>
00028 #include <finiteVolume/fvPatchFieldMapper.H>
00029 #include <finiteVolume/volFields.H>
00030 #include <OpenFOAM/addToRunTimeSelectionTable.H>
00031 #include <finiteVolume/wallFvPatch.H>
00032
00033
00034
00035 namespace Foam
00036 {
00037 namespace compressible
00038 {
00039 namespace RASModels
00040 {
00041
00042
00043
00044 void epsilonWallFunctionFvPatchScalarField::checkType()
00045 {
00046 if (!isA<wallFvPatch>(patch()))
00047 {
00048 FatalErrorIn("epsilonWallFunctionFvPatchScalarField::checkType()")
00049 << "Invalid wall function specification" << nl
00050 << " Patch type for patch " << patch().name()
00051 << " must be wall" << nl
00052 << " Current patch type is " << patch().type() << nl << endl
00053 << abort(FatalError);
00054 }
00055 }
00056
00057
00058
00059
00060 epsilonWallFunctionFvPatchScalarField::epsilonWallFunctionFvPatchScalarField
00061 (
00062 const fvPatch& p,
00063 const DimensionedField<scalar, volMesh>& iF
00064 )
00065 :
00066 fixedInternalValueFvPatchField<scalar>(p, iF),
00067 UName_("U"),
00068 kName_("k"),
00069 GName_("RASModel::G"),
00070 muName_("mu"),
00071 mutName_("mut"),
00072 Cmu_(0.09),
00073 kappa_(0.41),
00074 E_(9.8)
00075 {
00076 checkType();
00077 }
00078
00079
00080 epsilonWallFunctionFvPatchScalarField::epsilonWallFunctionFvPatchScalarField
00081 (
00082 const epsilonWallFunctionFvPatchScalarField& ptf,
00083 const fvPatch& p,
00084 const DimensionedField<scalar, volMesh>& iF,
00085 const fvPatchFieldMapper& mapper
00086 )
00087 :
00088 fixedInternalValueFvPatchField<scalar>(ptf, p, iF, mapper),
00089 UName_(ptf.UName_),
00090 kName_(ptf.kName_),
00091 GName_(ptf.GName_),
00092 muName_(ptf.muName_),
00093 mutName_(ptf.mutName_),
00094 Cmu_(ptf.Cmu_),
00095 kappa_(ptf.kappa_),
00096 E_(ptf.E_)
00097 {
00098 checkType();
00099 }
00100
00101
00102 epsilonWallFunctionFvPatchScalarField::epsilonWallFunctionFvPatchScalarField
00103 (
00104 const fvPatch& p,
00105 const DimensionedField<scalar, volMesh>& iF,
00106 const dictionary& dict
00107 )
00108 :
00109 fixedInternalValueFvPatchField<scalar>(p, iF, dict),
00110 UName_(dict.lookupOrDefault<word>("U", "U")),
00111 kName_(dict.lookupOrDefault<word>("k", "k")),
00112 GName_(dict.lookupOrDefault<word>("G", "RASModel::G")),
00113 muName_(dict.lookupOrDefault<word>("mu", "mu")),
00114 mutName_(dict.lookupOrDefault<word>("mut", "mut")),
00115 Cmu_(dict.lookupOrDefault<scalar>("Cmu", 0.09)),
00116 kappa_(dict.lookupOrDefault<scalar>("kappa", 0.41)),
00117 E_(dict.lookupOrDefault<scalar>("E", 9.8))
00118 {
00119 checkType();
00120 }
00121
00122
00123 epsilonWallFunctionFvPatchScalarField::epsilonWallFunctionFvPatchScalarField
00124 (
00125 const epsilonWallFunctionFvPatchScalarField& ewfpsf
00126 )
00127 :
00128 fixedInternalValueFvPatchField<scalar>(ewfpsf),
00129 UName_(ewfpsf.UName_),
00130 kName_(ewfpsf.kName_),
00131 GName_(ewfpsf.GName_),
00132 muName_(ewfpsf.muName_),
00133 mutName_(ewfpsf.mutName_),
00134 Cmu_(ewfpsf.Cmu_),
00135 kappa_(ewfpsf.kappa_),
00136 E_(ewfpsf.E_)
00137 {
00138 checkType();
00139 }
00140
00141
00142 epsilonWallFunctionFvPatchScalarField::epsilonWallFunctionFvPatchScalarField
00143 (
00144 const epsilonWallFunctionFvPatchScalarField& ewfpsf,
00145 const DimensionedField<scalar, volMesh>& iF
00146 )
00147 :
00148 fixedInternalValueFvPatchField<scalar>(ewfpsf, iF),
00149 UName_(ewfpsf.UName_),
00150 kName_(ewfpsf.kName_),
00151 GName_(ewfpsf.GName_),
00152 muName_(ewfpsf.muName_),
00153 mutName_(ewfpsf.mutName_),
00154 Cmu_(ewfpsf.Cmu_),
00155 kappa_(ewfpsf.kappa_),
00156 E_(ewfpsf.E_)
00157 {
00158 checkType();
00159 }
00160
00161
00162
00163
00164 void epsilonWallFunctionFvPatchScalarField::updateCoeffs()
00165 {
00166 if (updated())
00167 {
00168 return;
00169 }
00170
00171 const RASModel& rasModel = db().lookupObject<RASModel>("RASProperties");
00172
00173 const scalar Cmu25 = pow(Cmu_, 0.25);
00174 const scalar Cmu75 = pow(Cmu_, 0.75);
00175
00176 const scalarField& y = rasModel.y()[patch().index()];
00177
00178 volScalarField& G = const_cast<volScalarField&>
00179 (db().lookupObject<volScalarField>(GName_));
00180
00181 volScalarField& epsilon = const_cast<volScalarField&>
00182 (db().lookupObject<volScalarField>(dimensionedInternalField().name()));
00183
00184 const volScalarField& k = db().lookupObject<volScalarField>(kName_);
00185
00186 const scalarField& muw =
00187 patch().lookupPatchField<volScalarField, scalar>(muName_);
00188
00189 const scalarField& mutw =
00190 patch().lookupPatchField<volScalarField, scalar>(mutName_);
00191
00192 const fvPatchVectorField& Uw =
00193 patch().lookupPatchField<volVectorField, vector>(UName_);
00194
00195 const scalarField magGradUw = mag(Uw.snGrad());
00196
00197
00198 forAll(mutw, faceI)
00199 {
00200 label faceCellI = patch().faceCells()[faceI];
00201
00202 epsilon[faceCellI] = Cmu75*pow(k[faceCellI], 1.5)/(kappa_*y[faceI]);
00203
00204 G[faceCellI] =
00205 (mutw[faceI] + muw[faceI])
00206 *magGradUw[faceI]
00207 *Cmu25*sqrt(k[faceCellI])
00208 /(kappa_*y[faceI]);
00209 }
00210
00211
00212
00213 fixedInternalValueFvPatchField<scalar>::updateCoeffs();
00214 }
00215
00216
00217 void epsilonWallFunctionFvPatchScalarField::evaluate
00218 (
00219 const Pstream::commsTypes commsType
00220 )
00221 {
00222 fixedInternalValueFvPatchField<scalar>::evaluate(commsType);
00223 }
00224
00225
00226 void epsilonWallFunctionFvPatchScalarField::write(Ostream& os) const
00227 {
00228 fixedInternalValueFvPatchField<scalar>::write(os);
00229 writeEntryIfDifferent<word>(os, "U", "U", UName_);
00230 writeEntryIfDifferent<word>(os, "k", "k", kName_);
00231 writeEntryIfDifferent<word>(os, "G", "RASModel::G", GName_);
00232 writeEntryIfDifferent<word>(os, "mu", "mu", muName_);
00233 writeEntryIfDifferent<word>(os, "mut", "mut", mutName_);
00234 os.writeKeyword("Cmu") << Cmu_ << token::END_STATEMENT << nl;
00235 os.writeKeyword("kappa") << kappa_ << token::END_STATEMENT << nl;
00236 os.writeKeyword("E") << E_ << token::END_STATEMENT << nl;
00237 writeEntry("value", os);
00238 }
00239
00240
00241
00242
00243 makePatchTypeField
00244 (
00245 fvPatchScalarField,
00246 epsilonWallFunctionFvPatchScalarField
00247 );
00248
00249
00250
00251 }
00252 }
00253 }
00254
00255