FreeFOAM The Cross-Platform CFD Toolkit
RASModel.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
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
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
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.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "RASModel.H"
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 namespace incompressible
35 {
36 
37 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
38 
41 addToRunTimeSelectionTable(turbulenceModel, RASModel, turbulenceModel);
42 
43 // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
44 
45 void RASModel::printCoeffs()
46 {
47  if (printCoeffs_)
48  {
49  Info<< type() << "Coeffs" << coeffDict_ << endl;
50  }
51 }
52 
53 
54 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
55 
57 (
58  const word& type,
59  const volVectorField& U,
60  const surfaceScalarField& phi,
61  transportModel& lamTransportModel
62 )
63 :
64  turbulenceModel(U, phi, lamTransportModel),
65 
67  (
68  IOobject
69  (
70  "RASProperties",
71  U.time().constant(),
72  U.db(),
75  )
76  ),
77 
78  turbulence_(lookup("turbulence")),
79  printCoeffs_(lookupOrDefault<Switch>("printCoeffs", false)),
80  coeffDict_(subOrEmptyDict(type + "Coeffs")),
81 
82  k0_("k0", dimVelocity*dimVelocity, SMALL),
83  epsilon0_("epsilon0", k0_.dimensions()/dimTime, SMALL),
84  epsilonSmall_("epsilonSmall", epsilon0_.dimensions(), SMALL),
85  omega0_("omega0", dimless/dimTime, SMALL),
86  omegaSmall_("omegaSmall", omega0_.dimensions(), SMALL),
87 
88  y_(mesh_)
89 {
90  k0_.readIfPresent(*this);
91  epsilon0_.readIfPresent(*this);
92  epsilonSmall_.readIfPresent(*this);
93  omega0_.readIfPresent(*this);
94  omegaSmall_.readIfPresent(*this);
95 
96  // Force the construction of the mesh deltaCoeffs which may be needed
97  // for the construction of the derived models and BCs
98  mesh_.deltaCoeffs();
99 }
100 
101 
102 // * * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * //
103 
105 (
106  const volVectorField& U,
107  const surfaceScalarField& phi,
108  transportModel& transport
109 )
110 {
111  word modelName;
112 
113  // Enclose the creation of the dictionary to ensure it is deleted
114  // before the turbulenceModel is created otherwise the dictionary is
115  // entered in the database twice
116  {
117  IOdictionary dict
118  (
119  IOobject
120  (
121  "RASProperties",
122  U.time().constant(),
123  U.db(),
126  )
127  );
128 
129  dict.lookup("RASModel") >> modelName;
130  }
131 
132  Info<< "Selecting RAS turbulence model " << modelName << endl;
133 
134  dictionaryConstructorTable::iterator cstrIter =
135  dictionaryConstructorTablePtr_->find(modelName);
136 
137  if (cstrIter == dictionaryConstructorTablePtr_->end())
138  {
140  (
141  "RASModel::New(const volVectorField&, "
142  "const surfaceScalarField&, transportModel&)"
143  ) << "Unknown RASModel type " << modelName
144  << endl << endl
145  << "Valid RASModel types are :" << endl
146  << dictionaryConstructorTablePtr_->sortedToc()
147  << exit(FatalError);
148  }
149 
150  return autoPtr<RASModel>(cstrIter()(U, phi, transport));
151 }
152 
153 
154 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
155 
156 
157 scalar RASModel::yPlusLam(const scalar kappa, const scalar E) const
158 {
159  scalar ypl = 11.0;
160 
161  for (int i=0; i<10; i++)
162  {
163  ypl = log(max(E*ypl, 1))/kappa;
164  }
165 
166  return ypl;
167 }
168 
169 
170 tmp<scalarField> RASModel::yPlus(const label patchNo, const scalar Cmu) const
171 {
172  const fvPatch& curPatch = mesh_.boundary()[patchNo];
173 
174  tmp<scalarField> tYp(new scalarField(curPatch.size()));
175  scalarField& Yp = tYp();
176 
177  if (isA<wallFvPatch>(curPatch))
178  {
179  Yp = pow(Cmu, 0.25)
180  *y_[patchNo]
181  *sqrt(k()().boundaryField()[patchNo].patchInternalField())
182  /nu().boundaryField()[patchNo];
183  }
184  else
185  {
186  WarningIn
187  (
188  "tmp<scalarField> RASModel::yPlus(const label patchNo) const"
189  ) << "Patch " << patchNo << " is not a wall. Returning null field"
190  << nl << endl;
191 
192  Yp.setSize(0);
193  }
194 
195  return tYp;
196 }
197 
198 
200 {
202 
203  if (turbulence_ && mesh_.changing())
204  {
205  y_.correct();
206  }
207 }
208 
209 
210 bool RASModel::read()
211 {
212  if (regIOobject::read())
213  {
214  lookup("turbulence") >> turbulence_;
215 
216  if (const dictionary* dictPtr = subDictPtr(type() + "Coeffs"))
217  {
218  coeffDict_ <<= *dictPtr;
219  }
220 
221  k0_.readIfPresent(*this);
222  epsilon0_.readIfPresent(*this);
223  epsilonSmall_.readIfPresent(*this);
224  omega0_.readIfPresent(*this);
225  omegaSmall_.readIfPresent(*this);
226 
227  return true;
228  }
229  else
230  {
231  return false;
232  }
233 }
234 
235 
236 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
237 
238 } // End namespace incompressible
239 } // End namespace Foam
240 
241 // ************************ vim: set sw=4 sts=4 et: ************************ //