XiEqModel.H
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 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25 
26 Class
27  Foam::XiEqModel
28 
29 Description
30  Base-class for all XiEq models used by the b-XiEq combustion model.
31  The available models are :
32  \link basicXiSubXiEq.H \endlink
33  \link Gulder.H \endlink
34  \link instabilityXiEq.H \endlink
35  \link SCOPEBlendXiEq.H \endlink
36  \link SCOPEXiEq.H \endlink
37 
38 SourceFiles
39  XiEqModel.C
40 
41 \*---------------------------------------------------------------------------*/
42 
43 #ifndef XiEqModel_H
44 #define XiEqModel_H
45 
46 #include "IOdictionary.H"
47 #include "psiuReactionThermo.H"
49 #include "runTimeSelectionTables.H"
50 
51 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
52 
53 namespace Foam
54 {
55 
56 /*---------------------------------------------------------------------------*\
57  Class XiEqModel Declaration
58 \*---------------------------------------------------------------------------*/
59 
60 class XiEqModel
61 {
62 
63 protected:
64 
65  // Protected data
66 
67  //- Dictionary
69 
70  //- Thermo
72 
73  //- Turbulence
75 
76  //- Laminar burning velocity
77  const volScalarField& Su_;
78 
79 
80 private:
81 
82  // Private Member Functions
83 
84  //- No copy construct
85  XiEqModel(const XiEqModel&) = delete;
86 
87  //- No copy assignment
88  void operator=(const XiEqModel&) = delete;
89 
90 
91 public:
92 
93  //- Runtime type information
94  TypeName("XiEqModel");
95 
96 
97  // Declare run-time constructor selection table
98 
100  (
101  autoPtr,
102  XiEqModel,
103  dictionary,
104  (
105  const dictionary& XiEqProperties,
106  const psiuReactionThermo& thermo,
108  const volScalarField& Su
109  ),
110  (
111  XiEqProperties,
112  thermo,
113  turbulence,
114  Su
115  )
116  );
117 
118 
119  // Selectors
120 
121  //- Return a reference to the selected XiEq model
122  static autoPtr<XiEqModel> New
123  (
124  const dictionary& XiEqProperties,
125  const psiuReactionThermo& thermo,
127  const volScalarField& Su
128  );
129 
130 
131  // Constructors
132 
133  //- Construct from components
134  XiEqModel
135  (
136  const dictionary& XiEqProperties,
137  const psiuReactionThermo& thermo,
139  const volScalarField& Su
140  );
141 
142 
143  //- Destructor
144  virtual ~XiEqModel();
145 
146 
147  // Member Functions
148 
149  //- Return the flame-wrinkling XiEq
150  virtual tmp<volScalarField> XiEq() const
151  {
152  return turbulence_.muEff();
153  }
154 
155  //- Return the sub-grid Schelkin effect
156  tmp<volScalarField> calculateSchelkinEffect(const scalar) const;
157 
158  //- Update properties from given dictionary
159  virtual bool read(const dictionary& XiEqProperties) = 0;
160 
161  //- Write fields
162  void writeFields() const;
163 };
164 
165 
166 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
167 
168 } // End namespace Foam
169 
170 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
171 
172 #endif
173 
174 // ************************************************************************* //
Foam::XiEqModel
Base-class for all XiEq models used by the b-XiEq combustion model. The available models are : basicX...
Definition: XiEqModel.H:59
Foam::XiEqModel::calculateSchelkinEffect
tmp< volScalarField > calculateSchelkinEffect(const scalar) const
Return the sub-grid Schelkin effect.
psiuReactionThermo.H
Foam::XiEqModel::XiEqModelCoeffs_
dictionary XiEqModelCoeffs_
Dictionary.
Definition: XiEqModel.H:67
turbulence
Info<< "Reading field U\n"<< 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\n"<< endl;autoPtr< compressible::turbulenceModel > turbulence(compressible::turbulenceModel::New(rho, U, phi, thermo))
Definition: createFields.H:94
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::psiuReactionThermo
Foam::psiuReactionThermo.
Definition: psiuReactionThermo.H:55
thermo
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
Foam::XiEqModel::TypeName
TypeName("XiEqModel")
Runtime type information.
Foam::XiEqModel::XiEq
virtual tmp< volScalarField > XiEq() const
Return the flame-wrinkling XiEq.
Definition: XiEqModel.H:149
Foam::XiEqModel::declareRunTimeSelectionTable
declareRunTimeSelectionTable(autoPtr, XiEqModel, dictionary,(const dictionary &XiEqProperties, const psiuReactionThermo &thermo, const compressible::RASModel &turbulence, const volScalarField &Su),(XiEqProperties, thermo, turbulence, Su))
Su
zeroField Su
Definition: alphaSuSp.H:1
Foam::XiEqModel::read
virtual bool read(const dictionary &XiEqProperties)=0
Update properties from given dictionary.
Foam::XiEqModel::~XiEqModel
virtual ~XiEqModel()
Destructor.
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:121
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
IOdictionary.H
Foam::autoPtr
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: HashPtrTable.H:53
Foam::XiEqModel::Su_
const volScalarField & Su_
Laminar burning velocity.
Definition: XiEqModel.H:76
Foam::compressible::RASModel
RASModel< EddyDiffusivity< turbulenceModel > > RASModel
Definition: turbulentFluidThermoModel.H:66
runTimeSelectionTables.H
Macros to ease declaration of run-time selection tables.
Foam::XiEqModel::thermo_
const psiuReactionThermo & thermo_
Thermo.
Definition: XiEqModel.H:70
Foam::XiEqModel::writeFields
void writeFields() const
Write fields.
Foam::XiEqModel::New
static autoPtr< XiEqModel > New(const dictionary &XiEqProperties, const psiuReactionThermo &thermo, const compressible::RASModel &turbulence, const volScalarField &Su)
Return a reference to the selected XiEq model.
Foam::GeometricField< scalar, fvPatchField, volMesh >
Foam::XiEqModel::turbulence_
const compressible::RASModel & turbulence_
Turbulence.
Definition: XiEqModel.H:73
turbulentFluidThermoModel.H