FreeFOAM The Cross-Platform CFD Toolkit
turbulenceModel.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 "turbulenceModel.H"
27 #include <finiteVolume/volFields.H>
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 namespace incompressible
35 {
36 
37 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
38 
39 defineTypeNameAndDebug(turbulenceModel, 0);
40 defineRunTimeSelectionTable(turbulenceModel, turbulenceModel);
41 
42 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
43 
44 turbulenceModel::turbulenceModel
45 (
46  const volVectorField& U,
47  const surfaceScalarField& phi,
48  transportModel& lamTransportModel
49 )
50 :
51  runTime_(U.time()),
52  mesh_(U.mesh()),
53 
54  U_(U),
55  phi_(phi),
56  transportModel_(lamTransportModel)
57 {}
58 
59 
60 // * * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * //
61 
62 autoPtr<turbulenceModel> turbulenceModel::New
63 (
64  const volVectorField& U,
65  const surfaceScalarField& phi,
66  transportModel& transport
67 )
68 {
69  word modelName;
70 
71  // Enclose the creation of the dictionary to ensure it is deleted
72  // before the turbulenceModel is created otherwise the dictionary is
73  // entered in the database twice
74  {
75  IOdictionary dict
76  (
77  IOobject
78  (
79  "turbulenceProperties",
80  U.time().constant(),
81  U.db(),
84  )
85  );
86 
87  dict.lookup("simulationType") >> modelName;
88  }
89 
90  Info<< "Selecting turbulence model type " << modelName << endl;
91 
92  turbulenceModelConstructorTable::iterator cstrIter =
93  turbulenceModelConstructorTablePtr_->find(modelName);
94 
95  if (cstrIter == turbulenceModelConstructorTablePtr_->end())
96  {
98  (
99  "turbulenceModel::New(const volVectorField&, "
100  "const surfaceScalarField&, transportModel&)"
101  ) << "Unknown turbulenceModel type " << modelName
102  << endl << endl
103  << "Valid turbulenceModel types are :" << endl
104  << turbulenceModelConstructorTablePtr_->sortedToc()
105  << exit(FatalError);
106  }
107 
108  return autoPtr<turbulenceModel>(cstrIter()(U, phi, transport));
109 }
110 
111 
112 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
113 
115 {
116  transportModel_.correct();
117 }
118 
119 
120 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
121 
122 } // End namespace incompressible
123 } // End namespace Foam
124 
125 // ************************ vim: set sw=4 sts=4 et: ************************ //