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 "nutWallFunctionFvPatchScalarField.H"
00027 #include <incompressibleRASModels/RASModel.H>
00028 #include <finiteVolume/fvPatchFieldMapper.H>
00029 #include <finiteVolume/volFields.H>
00030 #include <finiteVolume/wallFvPatch.H>
00031 #include <OpenFOAM/addToRunTimeSelectionTable.H>
00032
00033
00034
00035 namespace Foam
00036 {
00037 namespace incompressible
00038 {
00039 namespace RASModels
00040 {
00041
00042
00043
00044 void nutWallFunctionFvPatchScalarField::checkType()
00045 {
00046 if (!isA<wallFvPatch>(patch()))
00047 {
00048 FatalErrorIn("nutWallFunctionFvPatchScalarField::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 scalar nutWallFunctionFvPatchScalarField::calcYPlusLam
00059 (
00060 const scalar kappa,
00061 const scalar E
00062 ) const
00063 {
00064 scalar ypl = 11.0;
00065
00066 for (int i=0; i<10; i++)
00067 {
00068 ypl = log(E*ypl)/kappa;
00069 }
00070
00071 return ypl;
00072 }
00073
00074
00075 tmp<scalarField> nutWallFunctionFvPatchScalarField::calcNut() const
00076 {
00077 const label patchI = patch().index();
00078
00079 const RASModel& rasModel = db().lookupObject<RASModel>("RASProperties");
00080 const scalarField& y = rasModel.y()[patchI];
00081 const tmp<volScalarField> tk = rasModel.k();
00082 const volScalarField& k = tk();
00083 const scalarField& nuw = rasModel.nu().boundaryField()[patchI];
00084
00085 const scalar Cmu25 = pow(Cmu_, 0.25);
00086
00087 tmp<scalarField> tnutw(new scalarField(patch().size(), 0.0));
00088 scalarField& nutw = tnutw();
00089
00090 forAll(nutw, faceI)
00091 {
00092 label faceCellI = patch().faceCells()[faceI];
00093
00094 scalar yPlus = Cmu25*y[faceI]*sqrt(k[faceCellI])/nuw[faceI];
00095
00096 if (yPlus > yPlusLam_)
00097 {
00098 nutw[faceI] = nuw[faceI]*(yPlus*kappa_/log(E_*yPlus) - 1.0);
00099 }
00100 }
00101
00102 return tnutw;
00103 }
00104
00105
00106 void nutWallFunctionFvPatchScalarField::writeLocalEntries(Ostream& os) const
00107 {
00108 os.writeKeyword("Cmu") << Cmu_ << token::END_STATEMENT << nl;
00109 os.writeKeyword("kappa") << kappa_ << token::END_STATEMENT << nl;
00110 os.writeKeyword("E") << E_ << token::END_STATEMENT << nl;
00111 }
00112
00113
00114
00115
00116 nutWallFunctionFvPatchScalarField::nutWallFunctionFvPatchScalarField
00117 (
00118 const fvPatch& p,
00119 const DimensionedField<scalar, volMesh>& iF
00120 )
00121 :
00122 fixedValueFvPatchScalarField(p, iF),
00123 Cmu_(0.09),
00124 kappa_(0.41),
00125 E_(9.8),
00126 yPlusLam_(calcYPlusLam(kappa_, E_))
00127 {
00128 checkType();
00129 }
00130
00131
00132 nutWallFunctionFvPatchScalarField::nutWallFunctionFvPatchScalarField
00133 (
00134 const nutWallFunctionFvPatchScalarField& ptf,
00135 const fvPatch& p,
00136 const DimensionedField<scalar, volMesh>& iF,
00137 const fvPatchFieldMapper& mapper
00138 )
00139 :
00140 fixedValueFvPatchScalarField(ptf, p, iF, mapper),
00141 Cmu_(ptf.Cmu_),
00142 kappa_(ptf.kappa_),
00143 E_(ptf.E_),
00144 yPlusLam_(ptf.yPlusLam_)
00145 {
00146 checkType();
00147 }
00148
00149
00150 nutWallFunctionFvPatchScalarField::nutWallFunctionFvPatchScalarField
00151 (
00152 const fvPatch& p,
00153 const DimensionedField<scalar, volMesh>& iF,
00154 const dictionary& dict
00155 )
00156 :
00157 fixedValueFvPatchScalarField(p, iF, dict),
00158 Cmu_(dict.lookupOrDefault<scalar>("Cmu", 0.09)),
00159 kappa_(dict.lookupOrDefault<scalar>("kappa", 0.41)),
00160 E_(dict.lookupOrDefault<scalar>("E", 9.8)),
00161 yPlusLam_(calcYPlusLam(kappa_, E_))
00162 {
00163 checkType();
00164 }
00165
00166
00167 nutWallFunctionFvPatchScalarField::nutWallFunctionFvPatchScalarField
00168 (
00169 const nutWallFunctionFvPatchScalarField& wfpsf
00170 )
00171 :
00172 fixedValueFvPatchScalarField(wfpsf),
00173 Cmu_(wfpsf.Cmu_),
00174 kappa_(wfpsf.kappa_),
00175 E_(wfpsf.E_),
00176 yPlusLam_(wfpsf.yPlusLam_)
00177 {
00178 checkType();
00179 }
00180
00181
00182 nutWallFunctionFvPatchScalarField::nutWallFunctionFvPatchScalarField
00183 (
00184 const nutWallFunctionFvPatchScalarField& wfpsf,
00185 const DimensionedField<scalar, volMesh>& iF
00186 )
00187 :
00188 fixedValueFvPatchScalarField(wfpsf, iF),
00189 Cmu_(wfpsf.Cmu_),
00190 kappa_(wfpsf.kappa_),
00191 E_(wfpsf.E_),
00192 yPlusLam_(wfpsf.yPlusLam_)
00193 {
00194 checkType();
00195 }
00196
00197
00198
00199
00200 void nutWallFunctionFvPatchScalarField::updateCoeffs()
00201 {
00202 operator==(calcNut());
00203
00204 fixedValueFvPatchScalarField::updateCoeffs();
00205 }
00206
00207
00208 tmp<scalarField> nutWallFunctionFvPatchScalarField::yPlus() const
00209 {
00210 const label patchI = patch().index();
00211
00212 const RASModel& rasModel = db().lookupObject<RASModel>("RASProperties");
00213 const scalarField& y = rasModel.y()[patchI];
00214
00215 const tmp<volScalarField> tk = rasModel.k();
00216 const volScalarField& k = tk();
00217 const scalarField kwc = k.boundaryField()[patchI].patchInternalField();
00218 const scalarField& nuw = rasModel.nu().boundaryField()[patchI];
00219
00220 return pow(Cmu_, 0.25)*y*sqrt(kwc)/nuw;
00221 }
00222
00223
00224 void nutWallFunctionFvPatchScalarField::write(Ostream& os) const
00225 {
00226 fvPatchField<scalar>::write(os);
00227 writeLocalEntries(os);
00228 writeEntry("value", os);
00229 }
00230
00231
00232
00233
00234 makePatchTypeField(fvPatchScalarField, nutWallFunctionFvPatchScalarField);
00235
00236
00237
00238 }
00239 }
00240 }
00241
00242