![]() |
|
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 Class 00026 locDynOneEqEddy 00027 00028 Description 00029 <pre> 00030 Localised Dynamic One Equation Eddy Viscosity Model 00031 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 00032 Eddy viscosity SGS model using a modeled balance equation to simulate 00033 the behaviour of k, hence 00034 00035 d/dt(k) + div(U*k) - div(nuSgs*grad(k)) 00036 = 00037 -B*L - ce*rho*k^3/2/delta 00038 00039 and 00040 00041 B = 2/3*k*I - 2*nuEff*dev(D) 00042 00043 where 00044 00045 nuEff = cD*delta^2*||D|| + nu 00046 00047 A dynamic procedure is here applied to evaluate ck and ce 00048 00049 ck=<L.M>/<M.M> 00050 00051 and 00052 00053 ce=<e*m>/<m*m> 00054 00055 where 00056 00057 K = 0.5*(F(U.U) - F(U).F(U)) 00058 L = (F(U*U) - F(U)*F(U) - 0.33*K*I) 00059 M = delta*(F(sqrt(k)*D) - 2*sqrt(K + filter(k))*F(D)) 00060 m = pow(K + F(k), 3.0/2.0)/(2*delta) - F(pow(k, 3.0/2.0))/delta 00061 e = 2*delta*ck*(F(sqrt(k)*(D && D)) - 2*sqrt(K + F(k))*(F(D) && F(D)))/ 00062 </pre> 00063 00064 SourceFiles 00065 locDynOneEqEddy.C 00066 00067 \*---------------------------------------------------------------------------*/ 00068 00069 #ifndef locDynOneEqEddy_H 00070 #define locDynOneEqEddy_H 00071 00072 #include "GenEddyVisc.H" 00073 #include "simpleFilter.H" 00074 #include "LESfilter.H" 00075 00076 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // 00077 00078 namespace Foam 00079 { 00080 namespace LESmodels 00081 { 00082 00083 /*---------------------------------------------------------------------------*\ 00084 Class locDynOneEqEddy Declaration 00085 \*---------------------------------------------------------------------------*/ 00086 00087 class locDynOneEqEddy 00088 : 00089 public GenEddyVisc 00090 { 00091 // Private data 00092 00093 volScalarField k_; 00094 00095 simpleFilter simpleFilter_; 00096 autoPtr<LESfilter> filterPtr_; 00097 LESfilter& filter_; 00098 00099 00100 // Private Member Functions 00101 00102 //- Calculate ck, ce by filtering the velocity field U. 00103 volScalarField ck(const volTensorField&, const volScalarField&) const; 00104 volScalarField ce(const volTensorField&, const volScalarField&) const; 00105 00106 // Disallow default bitwise copy construct and assignment 00107 locDynOneEqEddy(const locDynOneEqEddy&); 00108 locDynOneEqEddy& operator=(const locDynOneEqEddy&); 00109 00110 00111 public: 00112 00113 //- Runtime type information 00114 TypeName("locDynOneEqEddy"); 00115 00116 // Constructors 00117 00118 //- from components 00119 locDynOneEqEddy 00120 ( 00121 const volVectorField& U, 00122 const surfaceScalarField& phi, 00123 transportModel& transport 00124 ); 00125 00126 00127 // Destructor 00128 00129 ~locDynOneEqEddy(); 00130 00131 00132 // Member Functions 00133 00134 //- Return SGS kinetic energy 00135 tmp<volScalarField> k() const 00136 { 00137 return k_; 00138 } 00139 00140 //- Return the effective diffusivity for k 00141 tmp<volScalarField> DkEff() const 00142 { 00143 return tmp<volScalarField> 00144 ( 00145 new volScalarField("DkEff", nuSgs_ + nu()) 00146 ); 00147 } 00148 00149 //- Correct Eddy-Viscosity and related properties 00150 void correct(const tmp<volTensorField>& gradU); 00151 00152 //- Read turbulenceProperties dictionary 00153 bool read(); 00154 }; 00155 00156 00157 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // 00158 00159 } // End namespace LESmodels 00160 } // End namespace Foam 00161 00162 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // 00163 00164 #endif 00165 00166 // ************************************************************************* //