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