1 /*---------------------------------------------------------------------------*\
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
5 \\ / A nd | Copyright (C) 2008-2009 OpenCFD Ltd.
7 -------------------------------------------------------------------------------
9 This file is part of OpenFOAM.
11 OpenFOAM is free software; you can redistribute it and/or modify it
12 under the terms of the GNU General Public License as published by the
13 Free Software Foundation; either version 2 of the License, or (at your
14 option) any later version.
16 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 You should have received a copy of the GNU General Public License
22 along with OpenFOAM; if not, write to the Free Software Foundation,
23 Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
25 \*---------------------------------------------------------------------------*/
27 #include "mutkWallFunctionFvPatchScalarField.H"
29 #include "fvPatchFieldMapper.H"
30 #include "volFields.H"
31 #include "wallFvPatch.H"
32 #include "addToRunTimeSelectionTable.H"
34 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
38 namespace compressible
43 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
45 void mutkWallFunctionFvPatchScalarField::checkType()
47 if (!isA<wallFvPatch>(patch()))
49 FatalErrorIn("mutkWallFunctionFvPatchScalarField::checkType()")
50 << "Invalid wall function specification" << nl
51 << " Patch type for patch " << patch().name()
52 << " must be wall" << nl
53 << " Current patch type is " << patch().type() << nl << endl
59 scalar mutkWallFunctionFvPatchScalarField::calcYPlusLam
67 for (int i=0; i<10; i++)
69 ypl = log(E*ypl)/kappa;
76 tmp<scalarField> mutkWallFunctionFvPatchScalarField::calcMut() const
78 const label patchI = patch().index();
79 const RASModel& rasModel = db().lookupObject<RASModel>("RASProperties");
80 const scalarField& y = rasModel.y()[patchI];
81 const scalarField& rhow = rasModel.rho().boundaryField()[patchI];
82 const tmp<volScalarField> tk = rasModel.k();
83 const volScalarField& k = tk();
84 const scalarField& muw = rasModel.mu().boundaryField()[patchI];
86 const scalar Cmu25 = pow(Cmu_, 0.25);
88 tmp<scalarField> tmutw(new scalarField(patch().size(), 0.0));
89 scalarField& mutw = tmutw();
93 label faceCellI = patch().faceCells()[faceI];
96 Cmu25*y[faceI]*sqrt(k[faceCellI])/(muw[faceI]/rhow[faceI]);
98 if (yPlus > yPlusLam_)
100 mutw[faceI] = muw[faceI]*(yPlus*kappa_/log(E_*yPlus) - 1);
108 void mutkWallFunctionFvPatchScalarField::writeLocalEntries(Ostream& os) const
110 os.writeKeyword("Cmu") << Cmu_ << token::END_STATEMENT << nl;
111 os.writeKeyword("kappa") << kappa_ << token::END_STATEMENT << nl;
112 os.writeKeyword("E") << E_ << token::END_STATEMENT << nl;
116 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
118 mutkWallFunctionFvPatchScalarField::mutkWallFunctionFvPatchScalarField
121 const DimensionedField<scalar, volMesh>& iF
124 fixedValueFvPatchScalarField(p, iF),
128 yPlusLam_(calcYPlusLam(kappa_, E_))
132 mutkWallFunctionFvPatchScalarField::mutkWallFunctionFvPatchScalarField
134 const mutkWallFunctionFvPatchScalarField& ptf,
136 const DimensionedField<scalar, volMesh>& iF,
137 const fvPatchFieldMapper& mapper
140 fixedValueFvPatchScalarField(ptf, p, iF, mapper),
144 yPlusLam_(ptf.yPlusLam_)
148 mutkWallFunctionFvPatchScalarField::mutkWallFunctionFvPatchScalarField
151 const DimensionedField<scalar, volMesh>& iF,
152 const dictionary& dict
155 fixedValueFvPatchScalarField(p, iF, dict),
156 Cmu_(dict.lookupOrDefault<scalar>("Cmu", 0.09)),
157 kappa_(dict.lookupOrDefault<scalar>("kappa", 0.41)),
158 E_(dict.lookupOrDefault<scalar>("E", 9.8)),
159 yPlusLam_(calcYPlusLam(kappa_, E_))
163 mutkWallFunctionFvPatchScalarField::mutkWallFunctionFvPatchScalarField
165 const mutkWallFunctionFvPatchScalarField& wfpsf
168 fixedValueFvPatchScalarField(wfpsf),
170 kappa_(wfpsf.kappa_),
172 yPlusLam_(wfpsf.yPlusLam_)
176 mutkWallFunctionFvPatchScalarField::mutkWallFunctionFvPatchScalarField
178 const mutkWallFunctionFvPatchScalarField& wfpsf,
179 const DimensionedField<scalar, volMesh>& iF
182 fixedValueFvPatchScalarField(wfpsf, iF),
184 kappa_(wfpsf.kappa_),
186 yPlusLam_(wfpsf.yPlusLam_)
190 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
192 void mutkWallFunctionFvPatchScalarField::updateCoeffs()
199 operator==(calcMut());
201 fixedValueFvPatchScalarField::updateCoeffs();
205 tmp<scalarField> mutkWallFunctionFvPatchScalarField::yPlus() const
207 const label patchI = patch().index();
209 const RASModel& rasModel = db().lookupObject<RASModel>("RASProperties");
210 const scalarField& y = rasModel.y()[patchI];
212 const tmp<volScalarField> tk = rasModel.k();
213 const volScalarField& k = tk();
214 const scalarField kwc = k.boundaryField()[patchI].patchInternalField();
215 const scalarField& muw = rasModel.mu().boundaryField()[patchI];
216 const scalarField& rhow = rasModel.rho().boundaryField()[patchI];
218 return pow(Cmu_, 0.25)*y*sqrt(kwc)/(muw/rhow);
222 void mutkWallFunctionFvPatchScalarField::write(Ostream& os) const
224 fvPatchField<scalar>::write(os);
225 writeLocalEntries(os);
226 writeEntry("value", os);
230 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
232 makePatchTypeField(fvPatchScalarField, mutkWallFunctionFvPatchScalarField);
234 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
236 } // End namespace RASModels
237 } // End namespace compressible
238 } // End namespace Foam
240 // ************************************************************************* //