LESdelta.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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2015 OpenFOAM Foundation
9  Copyright (C) 2019-2021 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "LESdelta.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
38 }
39 
40 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
41 
43 (
44  const word& name,
46 )
47 :
48  turbulenceModel_(turbulence),
49  delta_
50  (
51  IOobject
52  (
53  name,
54  turbulence.mesh().time().timeName(),
55  turbulence.mesh(),
56  IOobject::NO_READ,
57  IOobject::NO_WRITE
58  ),
59  turbulence.mesh(),
61  calculatedFvPatchScalarField::typeName
62  )
63 {}
64 
65 
66 // * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //
67 
69 (
70  const word& name,
72  const dictionary& dict,
73  const word& lookupName
74 )
75 {
76  const word deltaType(dict.get<word>(lookupName));
77 
78  Info<< "Selecting LES " << lookupName << " type " << deltaType << endl;
79 
80  auto* ctorPtr = dictionaryConstructorTable(deltaType);
81 
82  if (!ctorPtr)
83  {
85  (
86  dict,
87  "LESdelta",
88  deltaType,
89  *dictionaryConstructorTablePtr_
90  ) << exit(FatalIOError);
91  }
92 
93  return autoPtr<LESdelta>(ctorPtr(name, turbulence, dict));
94 }
95 
96 
98 (
99  const word& name,
101  const dictionary& dict,
102  const dictionaryConstructorTableType& additionalConstructors,
103  const word& lookupName
104 )
105 {
106  const word deltaType(dict.get<word>(lookupName));
107 
108  Info<< "Selecting LES " << lookupName << " type " << deltaType << endl;
109 
110  // First any additional ones
111  {
112  auto ctorIter = additionalConstructors.cfind(deltaType);
113 
114  if (ctorIter.found())
115  {
116  return autoPtr<LESdelta>(ctorIter.val()(name, turbulence, dict));
117  }
118  }
119 
120  auto* ctorPtr = dictionaryConstructorTable(deltaType);
121 
122  if (!ctorPtr)
123  {
125  (
126  dict,
127  "LESdelta",
128  deltaType,
129  *dictionaryConstructorTablePtr_
130  );
131 
132  if (!additionalConstructors.empty())
133  {
135  << " and " << additionalConstructors.sortedToc() << nl;
136  }
137 
139  << exit(FatalIOError);
140  }
141 
142  return autoPtr<LESdelta>(ctorPtr(name, turbulence, dict));
143 }
144 
145 
146 // ************************************************************************* //
dictionary dict
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:125
LESdelta(const LESdelta &)=delete
No copy construct.
A list of keyword definitions, which are a keyword followed by a number of values (eg...
Definition: dictionary.H:120
constexpr char nl
The newline &#39;\n&#39; character (0x0a)
Definition: Ostream.H:49
Abstract base class for LES deltas.
Definition: LESdelta.H:49
ThermalDiffusivity< CompressibleTurbulenceModel< fluidThermo > > turbulenceModel
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:487
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a T. FatalIOError if not found, or if the number of tokens is incorrect.
Abstract base class for turbulence models (RAS, LES and laminar).
static autoPtr< LESdelta > New(const word &name, const turbulenceModel &turbulence, const dictionary &dict, const word &lookupName="delta")
Return a reference to the selected LES delta.
Definition: LESdelta.C:62
word timeName
Definition: getTimeIndex.H:3
dynamicFvMesh & mesh
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:52
A class for handling words, derived from Foam::string.
Definition: word.H:63
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
Info<< "Reading field U\"<< endl;volVectorField U(IOobject("U", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);volScalarField rho(IOobject("rho", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), thermo.rho());volVectorField rhoU(IOobject("rhoU", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), rho *U);volScalarField rhoE(IOobject("rhoE", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), rho *(e+0.5 *magSqr(U)));surfaceScalarField pos(IOobject("pos", runTime.timeName(), mesh), mesh, dimensionedScalar("pos", dimless, 1.0));surfaceScalarField neg(IOobject("neg", runTime.timeName(), mesh), mesh, dimensionedScalar("neg", dimless, -1.0));surfaceScalarField phi("phi", fvc::flux(rhoU));Info<< "Creating turbulence model\"<< endl;autoPtr< compressible::turbulenceModel > turbulence(compressible::turbulenceModel::New(rho, U, phi, thermo))
Definition: createFields.H:94
defineTypeNameAndDebug(combustionModel, 0)
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
messageStream Info
Information stream (stdout output on master, null elsewhere)
Defines the attributes of an object for which implicit objectRegistry management is supported...
Definition: IOobject.H:166
#define FatalIOErrorInLookup(ios, lookupTag, lookupName, lookupTable)
Report an error message using Foam::FatalIOError.
Definition: error.H:615
Namespace for OpenFOAM.
IOerror FatalIOError
Error stream (stdout output on all processes), with additional &#39;FOAM FATAL IO ERROR&#39; header text and ...