![]() |
|
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 GenSGSStress 00027 00028 Description 00029 General base class for all models that directly solve for the SGS stress 00030 tensor B. Contains tensor fields B (the SGS stress tensor) as well as 00031 scalar fields for k (SGS turbulent energy) gamma (SGS viscosity) and 00032 epsilon (SGS dissipation). 00033 00034 SourceFiles 00035 GenSGSStress.C 00036 00037 \*---------------------------------------------------------------------------*/ 00038 00039 #ifndef GenSGSStress_H 00040 #define GenSGSStress_H 00041 00042 #include "LESmodel.H" 00043 00044 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // 00045 00046 namespace Foam 00047 { 00048 namespace LESmodels 00049 { 00050 00051 /*---------------------------------------------------------------------------*\ 00052 Class GenSGSStress Declaration 00053 \*---------------------------------------------------------------------------*/ 00054 00055 class GenSGSStress 00056 : 00057 virtual public LESmodel 00058 { 00059 // Private Member Functions 00060 00061 // Disallow default bitwise copy construct and assignment 00062 GenSGSStress(const GenSGSStress&); 00063 GenSGSStress& operator=(const GenSGSStress&); 00064 00065 00066 protected: 00067 00068 dimensionedScalar ce_; 00069 00070 scalar couplingFactor_; 00071 00072 volTensorField B_; 00073 volScalarField nuSgs_; 00074 00075 00076 public: 00077 00078 // Constructors 00079 00080 //- Constructor from components 00081 GenSGSStress 00082 ( 00083 const volVectorField& U, 00084 const surfaceScalarField& phi, 00085 transportModel& transport 00086 ); 00087 00088 00089 // Destructor 00090 00091 virtual ~GenSGSStress() 00092 {} 00093 00094 00095 // Member Functions 00096 00097 //- Return the SGS stress. 00098 virtual tmp<volTensorField> B() const 00099 { 00100 return B_; 00101 } 00102 00103 //- Return the SGS turbulent kinetic energy. 00104 virtual tmp<volScalarField> k() const 00105 { 00106 return 0.5*tr(B_); 00107 } 00108 00109 //- Return the SGS turbulent dissipation. 00110 virtual tmp<volScalarField> epsilon() const 00111 { 00112 volScalarField K = k(); 00113 return ce_*K*sqrt(K)/delta(); 00114 } 00115 00116 //- Return the SGS viscosity. 00117 virtual tmp<volScalarField> nuSgs() const 00118 { 00119 return nuSgs_; 00120 } 00121 00122 //- Returns div(B). 00123 // This is the additional term due to the filtering of the NSE. 00124 virtual tmp<fvVectorMatrix> divB(volVectorField& U) const; 00125 00126 //- Read turbulenceProperties dictionary 00127 bool read(); 00128 }; 00129 00130 00131 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // 00132 00133 } // End namespace LESmodels 00134 } // End namespace Foam 00135 00136 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // 00137 00138 #endif 00139 00140 // ************************************************************************* //