Limited exponent to allow single-precision operation.
[OpenFOAM-1.5.x.git] / src / turbulenceModels / LES / incompressible / wallFunc / nuSgsWallFunction / nuSgsWallFunctionFvPatchScalarField.C
blobb3c617d97df047c0705e92326292b954f47cba0d
1 /*---------------------------------------------------------------------------*\
2   =========                 |
3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
4    \\    /   O peration     |
5     \\  /    A nd           | Copyright (C) 1991-2008 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 "nuSgsWallFunctionFvPatchScalarField.H"
28 #include "LESModel.H"
29 #include "fvPatchFieldMapper.H"
30 #include "volFields.H"
31 #include "addToRunTimeSelectionTable.H"
33 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
35 namespace Foam
37 namespace incompressible
39 namespace LESModels
42 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
44 nuSgsWallFunctionFvPatchScalarField::nuSgsWallFunctionFvPatchScalarField
46     const fvPatch& p,
47     const DimensionedField<scalar, volMesh>& iF
50     fixedValueFvPatchScalarField(p, iF)
54 nuSgsWallFunctionFvPatchScalarField::nuSgsWallFunctionFvPatchScalarField
56     const nuSgsWallFunctionFvPatchScalarField& ptf,
57     const fvPatch& p,
58     const DimensionedField<scalar, volMesh>& iF,
59     const fvPatchFieldMapper& mapper
62     fixedValueFvPatchScalarField(ptf, p, iF, mapper)
66 nuSgsWallFunctionFvPatchScalarField::nuSgsWallFunctionFvPatchScalarField
68     const fvPatch& p,
69     const DimensionedField<scalar, volMesh>& iF,
70     const dictionary& dict
73     fixedValueFvPatchScalarField(p, iF, dict)
77 nuSgsWallFunctionFvPatchScalarField::nuSgsWallFunctionFvPatchScalarField
79     const nuSgsWallFunctionFvPatchScalarField& tppsf
82     fixedValueFvPatchScalarField(tppsf)
86 nuSgsWallFunctionFvPatchScalarField::nuSgsWallFunctionFvPatchScalarField
88     const nuSgsWallFunctionFvPatchScalarField& tppsf,
89     const DimensionedField<scalar, volMesh>& iF
92     fixedValueFvPatchScalarField(tppsf, iF)
96 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
98 void nuSgsWallFunctionFvPatchScalarField::evaluate
100     const Pstream::commsTypes
103     const LESModel& sgsModel
104         = db().lookupObject<LESModel>("LESProperties");
106     scalar kappa = readScalar(sgsModel.lookup("kappa"));
108     scalar E = readScalar(sgsModel.subDict("wallFunctionCoeffs").lookup("E"));
110     const scalarField& ry = patch().deltaCoeffs();
112     const fvPatchVectorField& U =
113         patch().lookupPatchField<volVectorField, vector>("U");
115     scalarField magUp = mag(U.patchInternalField() - U);
117     const scalarField& nuw =
118         patch().lookupPatchField<volScalarField, scalar>("nu");
119     scalarField& nuSgsw = *this;
122     scalarField magFaceGradU = mag(U.snGrad());
124     forAll(nuSgsw, facei)
125     {
126         scalar magUpara = magUp[facei];
128         scalar utau = sqrt((nuSgsw[facei] + nuw[facei])*magFaceGradU[facei]);
130         if(utau > VSMALL)
131         {
132             int iter = 0;
133             scalar err = GREAT;
135             do
136             {
137                 scalar kUu = min(kappa*magUpara/utau, 50);
138                 scalar fkUu = exp(kUu) - 1 - kUu*(1 + 0.5*kUu);
140                 scalar f =
141                     - utau/(ry[facei]*nuw[facei])
142                     + magUpara/utau
143                     + 1/E*(fkUu - 1.0/6.0*kUu*sqr(kUu));
145                 scalar df =
146                     - 1.0/(ry[facei]*nuw[facei])
147                     - magUpara/sqr(utau)
148                     - 1/E*kUu*fkUu/utau;
150                 scalar utauNew = utau - f/df;
151                 err = mag((utau - utauNew)/utau);
152                 utau = utauNew;
154             } while (utau > VSMALL && err > 0.01 && ++iter < 10);
156             nuSgsw[facei] =
157                 max(sqr(max(utau, 0))/magFaceGradU[facei] - nuw[facei], 0.0);
158         }
159         else
160         {
161             nuSgsw[facei] = 0;
162         }
163     }
167 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
169 makePatchTypeField(fvPatchScalarField, nuSgsWallFunctionFvPatchScalarField);
171 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
173 } // End namespace LESModels
174 } // End namespace incompressible
175 } // End namespace Foam
177 // ************************************************************************* //