FreeFOAM The Cross-Platform CFD Toolkit
Hosted by SourceForge:
Get FreeFOAM at SourceForge.net.
            Fast, secure and Free Open Source software downloads

mutSpalartAllmarasStandardRoughWallFunctionFvPatchScalarField.C

Go to the documentation of this file.
00001 /*---------------------------------------------------------------------------*\
00002   =========                 |
00003   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
00004    \\    /   O peration     |
00005     \\  /    A nd           | Copyright (C) 1991-2010 OpenCFD Ltd.
00006      \\/     M anipulation  |
00007 -------------------------------------------------------------------------------
00008 License
00009     This file is part of OpenFOAM.
00010 
00011     OpenFOAM is free software: you can redistribute it and/or modify it
00012     under the terms of the GNU General Public License as published by
00013     the Free Software Foundation, either version 3 of the License, or
00014     (at your option) any later version.
00015 
00016     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
00017     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
00018     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
00019     for more details.
00020 
00021     You should have received a copy of the GNU General Public License
00022     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
00023 
00024 \*---------------------------------------------------------------------------*/
00025 
00026 #include "mutSpalartAllmarasStandardRoughWallFunctionFvPatchScalarField.H"
00027 #include <finiteVolume/fvPatchFieldMapper.H>
00028 #include <finiteVolume/volFields.H>
00029 #include <compressibleRASModels/RASModel.H>
00030 #include <OpenFOAM/addToRunTimeSelectionTable.H>
00031 
00032 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00033 
00034 namespace Foam
00035 {
00036 namespace compressible
00037 {
00038 namespace RASModels
00039 {
00040 
00041 // * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //
00042 
00043 tmp<scalarField>
00044 mutSpalartAllmarasStandardRoughWallFunctionFvPatchScalarField::calcYPlus
00045 (
00046     const scalarField& magUp
00047 ) const
00048 {
00049     const label patchI = patch().index();
00050 
00051     const RASModel& rasModel = db().lookupObject<RASModel>("RASProperties");
00052     const scalarField& y = rasModel.y()[patchI];
00053     const scalarField& muw = rasModel.mu().boundaryField()[patchI];
00054     const fvPatchScalarField& rho = rasModel.rho().boundaryField()[patchI];
00055 
00056     tmp<scalarField> tyPlus(new scalarField(patch().size(), 0.0));
00057     scalarField& yPlus = tyPlus();
00058 
00059     if (roughnessHeight_ > 0.0)
00060     {
00061         // Rough Walls
00062         const scalar c_1 = 1/(90 - 2.25) + roughnessConstant_;
00063         static const scalar c_2 = 2.25/(90 - 2.25);
00064         static const scalar c_3 = 2.0*atan(1.0)/log(90/2.25);
00065         static const scalar c_4 = c_3*log(2.25);
00066 
00067         //if (KsPlusBasedOnYPlus_)
00068         {
00069             // If KsPlus is based on YPlus the extra term added to the law
00070             // of the wall will depend on yPlus
00071             forAll(yPlus, facei)
00072             {
00073                 const scalar magUpara = magUp[facei];
00074                 const scalar Re = rho[facei]*magUpara*y[facei]/muw[facei];
00075                 const scalar kappaRe = kappa_*Re;
00076 
00077                 scalar yp = yPlusLam_;
00078                 const scalar ryPlusLam = 1.0/yp;
00079 
00080                 int iter = 0;
00081                 scalar yPlusLast = 0.0;
00082                 scalar dKsPlusdYPlus = roughnessHeight_/y[facei];
00083 
00084                 // Enforce the roughnessHeight to be less than the distance to
00085                 // the first cell centre.
00086                 if (dKsPlusdYPlus > 1)
00087                 {
00088                     dKsPlusdYPlus = 1;
00089                 }
00090 
00091                 // Additional tuning parameter (fudge factor) - nominally = 1
00092                 dKsPlusdYPlus *= roughnessFudgeFactor_;
00093 
00094                 do
00095                 {
00096                     yPlusLast = yp;
00097 
00098                     // The non-dimensional roughness height
00099                     scalar KsPlus = yp*dKsPlusdYPlus;
00100 
00101                     // The extra term in the law-of-the-wall
00102                     scalar G = 0.0;
00103 
00104                     scalar yPlusGPrime = 0.0;
00105 
00106                     if (KsPlus >= 90)
00107                     {
00108                         const scalar t_1 = 1 + roughnessConstant_*KsPlus;
00109                         G = log(t_1);
00110                         yPlusGPrime = roughnessConstant_*KsPlus/t_1;
00111                     }
00112                     else if (KsPlus > 2.25)
00113                     {
00114                         const scalar t_1 = c_1*KsPlus - c_2;
00115                         const scalar t_2 = c_3*log(KsPlus) - c_4;
00116                         const scalar sint_2 = sin(t_2);
00117                         const scalar logt_1 = log(t_1);
00118                         G = logt_1*sint_2;
00119                         yPlusGPrime =
00120                             (c_1*sint_2*KsPlus/t_1) + (c_3*logt_1*cos(t_2));
00121                     }
00122 
00123                     scalar denom = 1.0 + log(E_*yp) - G - yPlusGPrime;
00124                     if (mag(denom) > VSMALL)
00125                     {
00126                         yp = (kappaRe + yp*(1 - yPlusGPrime))/denom;
00127                     }
00128 
00129                 } while
00130                 (
00131                     mag(ryPlusLam*(yp - yPlusLast)) > 0.0001
00132                  && ++iter < 10
00133                  && yp > VSMALL
00134                 );
00135 
00136                 yPlus[facei] = max(0.0, yp);
00137             }
00138         }
00139     }
00140     else
00141     {
00142         // Smooth Walls
00143         forAll(yPlus, facei)
00144         {
00145             const scalar magUpara = magUp[facei];
00146             const scalar Re = rho[facei]*magUpara*y[facei]/muw[facei];
00147             const scalar kappaRe = kappa_*Re;
00148 
00149             scalar yp = yPlusLam_;
00150             const scalar ryPlusLam = 1.0/yp;
00151 
00152             int iter = 0;
00153             scalar yPlusLast = 0.0;
00154 
00155             do
00156             {
00157                 yPlusLast = yp;
00158                 yp = (kappaRe + yp)/(1.0 + log(E_*yp));
00159 
00160             } while
00161             (
00162                 mag(ryPlusLam*(yp - yPlusLast)) > 0.0001
00163              && ++iter < 10
00164             );
00165 
00166             yPlus[facei] = max(0.0, yp);
00167         }
00168     }
00169 
00170     return tyPlus;
00171 }
00172 
00173 
00174 tmp<scalarField>
00175 mutSpalartAllmarasStandardRoughWallFunctionFvPatchScalarField::calcMut() const
00176 {
00177     const label patchI = patch().index();
00178 
00179     const RASModel& rasModel = db().lookupObject<RASModel>("RASProperties");
00180     const scalarField& y = rasModel.y()[patchI];
00181     const fvPatchVectorField& Uw = rasModel.U().boundaryField()[patchI];
00182     const scalarField& muw = rasModel.mu().boundaryField()[patchI];
00183     const fvPatchScalarField& rho = rasModel.rho().boundaryField()[patchI];
00184 
00185     scalarField magUp = mag(Uw.patchInternalField() - Uw);
00186 
00187     tmp<scalarField> tyPlus = calcYPlus(magUp);
00188     scalarField& yPlus = tyPlus();
00189 
00190     tmp<scalarField> tmutw(new scalarField(patch().size(), 0.0));
00191     scalarField& mutw = tmutw();
00192 
00193     forAll(yPlus, facei)
00194     {
00195         if (yPlus[facei] > yPlusLam_)
00196         {
00197             const scalar Re = rho[facei]*magUp[facei]*y[facei]/muw[facei];
00198             mutw[facei] = muw[facei]*(sqr(yPlus[facei])/Re - 1);
00199         }
00200     }
00201 
00202     return tmutw;
00203 }
00204 
00205 
00206 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
00207 
00208 mutSpalartAllmarasStandardRoughWallFunctionFvPatchScalarField::
00209 mutSpalartAllmarasStandardRoughWallFunctionFvPatchScalarField
00210 (
00211     const fvPatch& p,
00212     const DimensionedField<scalar, volMesh>& iF
00213 )
00214 :
00215     mutWallFunctionFvPatchScalarField(p, iF),
00216     roughnessHeight_(pTraits<scalar>::zero),
00217     roughnessConstant_(pTraits<scalar>::zero),
00218     roughnessFudgeFactor_(pTraits<scalar>::zero)
00219 {}
00220 
00221 
00222 mutSpalartAllmarasStandardRoughWallFunctionFvPatchScalarField::
00223 mutSpalartAllmarasStandardRoughWallFunctionFvPatchScalarField
00224 (
00225     const mutSpalartAllmarasStandardRoughWallFunctionFvPatchScalarField& ptf,
00226     const fvPatch& p,
00227     const DimensionedField<scalar, volMesh>& iF,
00228     const fvPatchFieldMapper& mapper
00229 )
00230 :
00231     mutWallFunctionFvPatchScalarField(ptf, p, iF, mapper),
00232     roughnessHeight_(ptf.roughnessHeight_),
00233     roughnessConstant_(ptf.roughnessConstant_),
00234     roughnessFudgeFactor_(ptf.roughnessFudgeFactor_)
00235 {}
00236 
00237 
00238 mutSpalartAllmarasStandardRoughWallFunctionFvPatchScalarField::
00239 mutSpalartAllmarasStandardRoughWallFunctionFvPatchScalarField
00240 (
00241     const fvPatch& p,
00242     const DimensionedField<scalar, volMesh>& iF,
00243     const dictionary& dict
00244 )
00245 :
00246     mutWallFunctionFvPatchScalarField(p, iF, dict),
00247     roughnessHeight_(readScalar(dict.lookup("roughnessHeight"))),
00248     roughnessConstant_(readScalar(dict.lookup("roughnessConstant"))),
00249     roughnessFudgeFactor_(readScalar(dict.lookup("roughnessFudgeFactor")))
00250 {}
00251 
00252 
00253 mutSpalartAllmarasStandardRoughWallFunctionFvPatchScalarField::
00254 mutSpalartAllmarasStandardRoughWallFunctionFvPatchScalarField
00255 (
00256     const mutSpalartAllmarasStandardRoughWallFunctionFvPatchScalarField& rwfpsf
00257 )
00258 :
00259     mutWallFunctionFvPatchScalarField(rwfpsf),
00260     roughnessHeight_(rwfpsf.roughnessHeight_),
00261     roughnessConstant_(rwfpsf.roughnessConstant_),
00262     roughnessFudgeFactor_(rwfpsf.roughnessFudgeFactor_)
00263 {}
00264 
00265 
00266 mutSpalartAllmarasStandardRoughWallFunctionFvPatchScalarField::
00267 mutSpalartAllmarasStandardRoughWallFunctionFvPatchScalarField
00268 (
00269     const mutSpalartAllmarasStandardRoughWallFunctionFvPatchScalarField& rwfpsf,
00270     const DimensionedField<scalar, volMesh>& iF
00271 )
00272 :
00273     mutWallFunctionFvPatchScalarField(rwfpsf, iF),
00274     roughnessHeight_(rwfpsf.roughnessHeight_),
00275     roughnessConstant_(rwfpsf.roughnessConstant_),
00276     roughnessFudgeFactor_(rwfpsf.roughnessFudgeFactor_)
00277 {}
00278 
00279 
00280 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
00281 
00282 tmp<scalarField>
00283 mutSpalartAllmarasStandardRoughWallFunctionFvPatchScalarField::yPlus() const
00284 {
00285     const label patchI = patch().index();
00286 
00287     const RASModel& rasModel = db().lookupObject<RASModel>("RASProperties");
00288     const fvPatchVectorField& Uw = rasModel.U().boundaryField()[patchI];
00289     const scalarField magUp = mag(Uw.patchInternalField() - Uw);
00290 
00291     return calcYPlus(magUp);
00292 }
00293 
00294 
00295 void mutSpalartAllmarasStandardRoughWallFunctionFvPatchScalarField::write
00296 (
00297     Ostream& os
00298 ) const
00299 {
00300     fixedValueFvPatchScalarField::write(os);
00301     writeLocalEntries(os);
00302     os.writeKeyword("roughnessHeight")
00303         << roughnessHeight_ << token::END_STATEMENT << nl;
00304     os.writeKeyword("roughnessConstant")
00305         << roughnessConstant_ << token::END_STATEMENT << nl;
00306     os.writeKeyword("roughnessFudgeFactor")
00307         << roughnessFudgeFactor_ << token::END_STATEMENT << nl;
00308 }
00309 
00310 
00311 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00312 
00313 makePatchTypeField
00314 (
00315     fvPatchScalarField,
00316     mutSpalartAllmarasStandardRoughWallFunctionFvPatchScalarField
00317 );
00318 
00319 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00320 
00321 } // End namespace RASModels
00322 } // End namespace compressible
00323 } // End namespace Foam
00324 
00325 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines