nutUSpaldingWallFunctionFvPatchScalarField.H 7.59 KB
Newer Older
1 2 3 4
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
OpenFOAM bot's avatar
OpenFOAM bot committed
5
    \\  /    A nd           | www.openfoam.com
OpenFOAM bot's avatar
OpenFOAM bot committed
6 7
     \\/     M anipulation  |
-------------------------------------------------------------------------------
OpenFOAM bot's avatar
OpenFOAM bot committed
8
    Copyright (C) 2011-2016 OpenFOAM Foundation
9
    Copyright (C) 2019-2020 OpenCFD Ltd.
10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33
-------------------------------------------------------------------------------
License
    This file is part of OpenFOAM.

    OpenFOAM is free software: you can redistribute it and/or modify it
    under the terms of the GNU General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.

    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
    for more details.

    You should have received a copy of the GNU General Public License
    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.

Class
    Foam::nutUSpaldingWallFunctionFvPatchScalarField

Group
    grpWallFunctions

Description
34
    This boundary condition provides a wall constraint on the turbulent
35 36
    viscosity, i.e. \c nut, based on velocity, i.e. \c U. Using Spalding's
    law gives a continuous \c nut profile to the wall.
37 38 39 40 41 42 43 44

        \f[
            y^+ = u^+ + \frac{1}{E} \left[exp(\kappa u^+) - 1 - \kappa u^+\,
                - 0.5 (\kappa u^+)^2 - \frac{1}{6} (\kappa u^+)^3\right]
        \f]

    where
    \vartable
45 46 47
        y^+     | Non-dimensional position
        u^+     | Non-dimensional velocity
        \kappa  | von Kármán constant
48 49
    \endvartable

50
Usage
51 52
    Example of the boundary condition specification:
    \verbatim
53
    <patchName>
54
    {
55
        // Mandatory entries (unmodifiable)
56
        type            nutUSpaldingWallFunction;
57

58 59 60 61 62 63
        // Optional entries (unmodifiable)
        maxIter         10;
        tolerance       0.0001;

        // Optional (inherited) entries
        ...
64 65 66
    }
    \endverbatim

67 68 69 70 71 72 73 74 75 76
    where the entries mean:
    \table
      Property  | Description                         | Type   | Req'd  | Dflt
      type      | Type name: nutUBlendedWallFunction  | word   | yes    | -
      maxIter   | Number of Newton-Raphson iterations | label  | no     | 10
      tolerance | Convergence tolerance               | scalar | no     | 0.0001
    \endtable

    The inherited entries are elaborated in:
      - \link nutWallFunctionFvPatchScalarField.H \endlink
77

78
Note
79 80
    - Suffers from non-exact restart since \c correctNut() (called through
    \c turbulence->validate) returns a slightly different value every time
81
    it is called. This is since the seed for the Newton-Raphson iteration
82 83
    uses the current value of \c *this (\c =nut ).
    - This can be avoided by overriding the tolerance. This also switches on
mattijs's avatar
mattijs committed
84 85 86
    a pre-detection whether the current nut already satisfies the turbulence
    conditions and if so does not change it at all. This means that the nut
    only changes if it 'has really changed'. This probably should be used with
87
    a tight tolerance, to make sure to kick every iteration, e.g.
mattijs's avatar
mattijs committed
88 89 90
        maxIter     100;
        tolerance   1e-7;

91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115
SourceFiles
    nutUSpaldingWallFunctionFvPatchScalarField.C

\*---------------------------------------------------------------------------*/

#ifndef nutUSpaldingWallFunctionFvPatchScalarField_H
#define nutUSpaldingWallFunctionFvPatchScalarField_H

#include "nutWallFunctionFvPatchScalarField.H"

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

namespace Foam
{

/*---------------------------------------------------------------------------*\
          Class nutUSpaldingWallFunctionFvPatchScalarField Declaration
\*---------------------------------------------------------------------------*/

class nutUSpaldingWallFunctionFvPatchScalarField
:
    public nutWallFunctionFvPatchScalarField
{
protected:

116
    // Protected Data
117

118 119 120 121 122 123 124 125 126 127 128
        //- Max iterations in calcNut
        const label maxIter_;

        //- Convergence tolerance
        const scalar tolerance_;

        //- Uncomment in case of intrumentation
        //mutable uint64_t invocations_;
        //mutable uint64_t nontrivial_;
        //mutable uint64_t nonconvergence_;
        //mutable uint64_t iterations_;
129 130


131 132
    // Protected Member Functions

133
        //- Calculate the turbulent viscosity
134 135 136 137 138
        virtual tmp<scalarField> calcNut() const;

        //- Calculate the friction velocity
        virtual tmp<scalarField> calcUTau(const scalarField& magGradU) const;

mattijs's avatar
mattijs committed
139
        //- Calculate the friction velocity and number of iterations for
140
        //- convergence
mattijs's avatar
mattijs committed
141 142 143 144 145 146 147
        virtual tmp<scalarField> calcUTau
        (
            const scalarField& magGradU,
            const label maxIter,
            scalarField& err
        ) const;

148 149 150
        //- Write local wall function variables
        virtual void writeLocalEntries(Ostream&) const;

151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175

public:

    //- Runtime type information
    TypeName("nutUSpaldingWallFunction");


    // Constructors

        //- Construct from patch and internal field
        nutUSpaldingWallFunctionFvPatchScalarField
        (
            const fvPatch&,
            const DimensionedField<scalar, volMesh>&
        );

        //- Construct from patch, internal field and dictionary
        nutUSpaldingWallFunctionFvPatchScalarField
        (
            const fvPatch&,
            const DimensionedField<scalar, volMesh>&,
            const dictionary&
        );

        //- Construct by mapping given
176 177
        //- nutUSpaldingWallFunctionFvPatchScalarField
        //- onto a new patch
178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220
        nutUSpaldingWallFunctionFvPatchScalarField
        (
            const nutUSpaldingWallFunctionFvPatchScalarField&,
            const fvPatch&,
            const DimensionedField<scalar, volMesh>&,
            const fvPatchFieldMapper&
        );

        //- Construct as copy
        nutUSpaldingWallFunctionFvPatchScalarField
        (
            const nutUSpaldingWallFunctionFvPatchScalarField&
        );

        //- Construct and return a clone
        virtual tmp<fvPatchScalarField> clone() const
        {
            return tmp<fvPatchScalarField>
            (
                new nutUSpaldingWallFunctionFvPatchScalarField(*this)
            );
        }

        //- Construct as copy setting internal field reference
        nutUSpaldingWallFunctionFvPatchScalarField
        (
            const nutUSpaldingWallFunctionFvPatchScalarField&,
            const DimensionedField<scalar, volMesh>&
        );

        //- Construct and return a clone setting internal field reference
        virtual tmp<fvPatchScalarField> clone
        (
            const DimensionedField<scalar, volMesh>& iF
        ) const
        {
            return tmp<fvPatchScalarField>
            (
                new nutUSpaldingWallFunctionFvPatchScalarField(*this, iF)
            );
        }


221 222 223 224
    //- Destructor
    virtual ~nutUSpaldingWallFunctionFvPatchScalarField();


225
    // Member Functions
226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248

        // Evaluation functions

            //- Calculate and return the yPlus at the boundary
            virtual tmp<scalarField> yPlus() const;


        // I-O

            //- Write
            virtual void write(Ostream& os) const;
};


// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

} // End namespace Foam

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

#endif

// ************************************************************************* //