New BCs: nutkWallFunction (equivalent to nutWallFunction)
[OpenFOAM-1.6.x.git] / src / turbulenceModels / compressible / RAS / derivedFvPatchFields / wallFunctions / mutWallFunctions / mutkWallFunction / mutkWallFunctionFvPatchScalarField.C
blob5ab8324cda85beb987ea456744ebee6ae8e4356c
1 /*---------------------------------------------------------------------------*\
2   =========                 |
3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
4    \\    /   O peration     |
5     \\  /    A nd           | Copyright (C) 2008-2009 OpenCFD Ltd.
6      \\/     M anipulation  |
7 -------------------------------------------------------------------------------
8 License
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
19     for more details.
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"
28 #include "RASModel.H"
29 #include "fvPatchFieldMapper.H"
30 #include "volFields.H"
31 #include "wallFvPatch.H"
32 #include "addToRunTimeSelectionTable.H"
34 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
36 namespace Foam
38 namespace compressible
40 namespace RASModels
43 // * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //
45 void mutkWallFunctionFvPatchScalarField::checkType()
47     if (!isA<wallFvPatch>(patch()))
48     {
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
54             << abort(FatalError);
55     }
59 scalar mutkWallFunctionFvPatchScalarField::calcYPlusLam
61     const scalar kappa,
62     const scalar E
63 ) const
65     scalar ypl = 11.0;
67     for (int i=0; i<10; i++)
68     {
69         ypl = log(E*ypl)/kappa;
70     }
72     return ypl;
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();
91     forAll(mutw, faceI)
92     {
93         label faceCellI = patch().faceCells()[faceI];
95         scalar yPlus =
96             Cmu25*y[faceI]*sqrt(k[faceCellI])/(muw[faceI]/rhow[faceI]);
98         if (yPlus > yPlusLam_)
99         {
100             mutw[faceI] = muw[faceI]*(yPlus*kappa_/log(E_*yPlus) - 1);
101         }
102     }
104     return tmutw;
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
120     const fvPatch& p,
121     const DimensionedField<scalar, volMesh>& iF
124     fixedValueFvPatchScalarField(p, iF),
125     Cmu_(0.09),
126     kappa_(0.41),
127     E_(9.8),
128     yPlusLam_(calcYPlusLam(kappa_, E_))
132 mutkWallFunctionFvPatchScalarField::mutkWallFunctionFvPatchScalarField
134     const mutkWallFunctionFvPatchScalarField& ptf,
135     const fvPatch& p,
136     const DimensionedField<scalar, volMesh>& iF,
137     const fvPatchFieldMapper& mapper
140     fixedValueFvPatchScalarField(ptf, p, iF, mapper),
141     Cmu_(ptf.Cmu_),
142     kappa_(ptf.kappa_),
143     E_(ptf.E_),
144     yPlusLam_(ptf.yPlusLam_)
148 mutkWallFunctionFvPatchScalarField::mutkWallFunctionFvPatchScalarField
150     const fvPatch& p,
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),
169     Cmu_(wfpsf.Cmu_),
170     kappa_(wfpsf.kappa_),
171     E_(wfpsf.E_),
172     yPlusLam_(wfpsf.yPlusLam_)
176 mutkWallFunctionFvPatchScalarField::mutkWallFunctionFvPatchScalarField
178     const mutkWallFunctionFvPatchScalarField& wfpsf,
179     const DimensionedField<scalar, volMesh>& iF
182     fixedValueFvPatchScalarField(wfpsf, iF),
183     Cmu_(wfpsf.Cmu_),
184     kappa_(wfpsf.kappa_),
185     E_(wfpsf.E_),
186     yPlusLam_(wfpsf.yPlusLam_)
190 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
192 void mutkWallFunctionFvPatchScalarField::updateCoeffs()
194     if (updated())
195     {
196         return;
197     }
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 // ************************************************************************* //