OpenFOAM logo
Open Source CFD Toolkit

wallFunctionsI.H

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-2005 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 the
00013     Free Software Foundation; either version 2 of the License, or (at your
00014     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, write to the Free Software Foundation,
00023     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
00024 
00025 Global
00026     wallFunctions
00027 
00028 Description
00029     Calculate wall dissipation from wall-functions.
00030 
00031 \*---------------------------------------------------------------------------*/
00032 
00033 {
00034     labelList cellBoundaryFaceCount(epsilon_.size(), 0);
00035 
00036     scalar Cmu25 = pow(Cmu.value(), 0.25);
00037     scalar Cmu75 = pow(Cmu.value(), 0.75);
00038 
00039     const fvPatchList& patches = mesh_.boundary();
00040 
00041     //- Initialise the near-wall G field to zero
00042     forAll(patches, patchi)
00043     {
00044         const fvPatch& curPatch = patches[patchi];
00045 
00046         if (isType<wallFvPatch>(curPatch))
00047         {
00048             forAll(curPatch, facei)
00049             {
00050                 label faceCelli = curPatch.faceCells()[facei];
00051 
00052                 epsilon_[faceCelli] = 0.0;
00053                 G[faceCelli] = 0.0;
00054             }
00055         }
00056     }
00057 
00058     //- Accumulate the wall face contributions to epsilon and G
00059     //  Increment cellBoundaryFaceCount for each face for averaging
00060     forAll(patches, patchi)
00061     {
00062         const fvPatch& curPatch = patches[patchi];
00063 
00064         if (isType<wallFvPatch>(curPatch))
00065         {
00066 #           include "checkPatchFieldTypes.H"
00067 
00068             const scalarField& rhow = rho_.boundaryField()[patchi];
00069 
00070             const scalarField muw = mu().boundaryField()[patchi];
00071             const scalarField& mutw = mut_.boundaryField()[patchi];
00072 
00073             scalarField magFaceGradU = mag(U_.boundaryField()[patchi].snGrad());
00074 
00075             forAll(curPatch, facei)
00076             {
00077                 label faceCelli = curPatch.faceCells()[facei];
00078 
00079                 scalar yPlus =
00080                     Cmu25*turbulenceModel::y_[patchi][facei]*sqrt(k_[faceCelli])
00081                     /(muw[facei]/rhow[facei]);
00082 
00083                 // For corner cells (with two boundary or more faces),
00084                 // epsilon and G in the near-wall cell are calculated
00085                 // as an average
00086 
00087                 cellBoundaryFaceCount[faceCelli]++;
00088 
00089                 epsilon_[faceCelli] +=
00090                     Cmu75*pow(k_[faceCelli], 1.5)
00091                    /(kappa_*turbulenceModel::y_[patchi][facei]);
00092 
00093                 if (yPlus > yPlusLam_)
00094                 {
00095                     G[faceCelli] +=
00096                         mutw[facei]*magFaceGradU[facei]
00097                       //- *0.5*log(2.0*yPlus/yPlusLam_)    //  UMIST correction
00098                         *Cmu25*sqrt(k_[faceCelli])
00099                         /(kappa_*turbulenceModel::y_[patchi][facei]);
00100                 }
00101             }
00102         }
00103     }
00104 
00105 
00106     // Perform the averaging
00107 
00108     forAll(patches, patchi)
00109     {
00110         const fvPatch& curPatch = patches[patchi];
00111 
00112         if (isType<wallFvPatch>(curPatch))
00113         {
00114             forAll(curPatch, facei)
00115             {
00116                 label faceCelli = curPatch.faceCells()[facei];
00117 
00118                 epsilon_[faceCelli] /= cellBoundaryFaceCount[faceCelli];
00119                 G[faceCelli] /= cellBoundaryFaceCount[faceCelli];
00120             }
00121         }
00122     }
00123 }
00124 
00125 
00126 // ************************************************************************* //
For further information go to www.openfoam.org