XiEqModel.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) 2011-2015 OpenFOAM Foundation
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 "XiEqModel.H"
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33  defineTypeNameAndDebug(XiEqModel, 0);
34  defineRunTimeSelectionTable(XiEqModel, dictionary);
35 }
36 
37 
38 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
39 
40 Foam::XiEqModel::XiEqModel
41 (
42  const dictionary& XiEqProperties,
43  const psiuReactionThermo& thermo,
45  const volScalarField& Su
46 )
47 :
48  XiEqModelCoeffs_
49  (
50  XiEqProperties.subDict
51  (
52  word(XiEqProperties.lookup("XiEqModel")) + "Coeffs"
53  )
54  ),
55  thermo_(thermo),
56  turbulence_(turbulence),
57  Su_(Su)
58 {}
59 
60 
61 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
62 
64 {}
65 
66 
67 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
68 
69 bool Foam::XiEqModel::read(const dictionary& XiEqProperties)
70 {
71  XiEqModelCoeffs_ = XiEqProperties.subDict(type() + "Coeffs");
72 
73  return true;
74 }
75 
76 
78 {
79  //***HGW It is not clear why B is written here
80  if (Su_.mesh().foundObject<volSymmTensorField>("B"))
81  {
82  const volSymmTensorField& B =
83  Su_.mesh().lookupObject<volSymmTensorField>("B");
84  B.write();
85  }
86 }
87 
88 
90 Foam::XiEqModel::calculateSchelkinEffect(const scalar uPrimeCoef) const
91 {
92  const fvMesh& mesh = Su_.mesh();
93 
94  const volVectorField& U = mesh.lookupObject<volVectorField>("U");
95  const volSymmTensorField& CT = mesh.lookupObject<volSymmTensorField>("CT");
96  const volScalarField& Nv = mesh.lookupObject<volScalarField>("Nv");
97  const volSymmTensorField& nsv =
98  mesh.lookupObject<volSymmTensorField>("nsv");
99 
100  tmp<volScalarField> tN
101  (
102  new volScalarField
103  (
104  IOobject
105  (
106  "tN",
107  mesh.time().timeName(),
108  mesh,
111  false
112  ),
113  mesh,
114  dimensionedScalar("zero", Nv.dimensions(), 0.0),
115  zeroGradientFvPatchVectorField::typeName
116  )
117  );
118 
119  volScalarField& N = tN();
120 
121  N.internalField() = Nv.internalField()*pow(mesh.V(), 2.0/3.0);
122 
123  tmp<volSymmTensorField> tns
124  (
126  (
127  IOobject
128  (
129  "tns",
130  mesh.time().timeName(),
131  mesh,
134  ),
135  mesh,
137  (
138  "zero",
139  nsv.dimensions(),
140  pTraits<symmTensor>::zero
141  )
142  )
143  );
144 
145  volSymmTensorField& ns = tns();
146 
147  ns.internalField() = nsv.internalField()*pow(mesh.V(), 2.0/3.0);
148 
149  const volVectorField Uhat
150  (
151  U/(mag(U) + dimensionedScalar("Usmall", U.dimensions(), 1e-4))
152  );
153 
154  const volScalarField nr(sqrt(max(N - (Uhat & ns & Uhat), scalar(1e-4))));
155 
156  const scalarField cellWidth(pow(mesh.V(), 1.0/3.0));
157 
158  const scalarField upLocal(uPrimeCoef*sqrt((U & CT & U)*cellWidth));
159 
160  const scalarField deltaUp(upLocal*(max(scalar(1.0), pow(nr, 0.5)) - 1.0));
161 
162  // Re use tN
163  N.internalField() = upLocal*(max(scalar(1.0), pow(nr, 0.5)) - 1.0);
164 
165  return tN;
166 }
167 
168 
169 // ************************************************************************* //
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
void writeFields() const
Write fields.
dimensioned< scalar > mag(const dimensioned< Type > &)
dictionary XiEqModelCoeffs_
Dictionary.
Definition: XiEqModel.H:65
autoPtr< compressible::turbulenceModel > turbulence
Definition: createFields.H:23
const Mesh & mesh() const
Return mesh.
dynamicFvMesh & mesh
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:638
const volScalarField & Su_
Laminar burning velocity.
Definition: XiEqModel.H:74
RASModel< EddyDiffusivity< turbulenceModel > > RASModel
Namespace for OpenFOAM.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
tmp< GeometricField< Type, fvPatchField, volMesh > > Su(const GeometricField< Type, fvPatchField, volMesh > &su, const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvcSup.C:44
virtual bool write() const
Write using setting from DB.
virtual ~XiEqModel()
Destructor.
const double e
Elementary charge.
Definition: doubleFloat.H:78
stressControl lookup("compactNormalStress") >> compactNormalStress
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
GeometricField< symmTensor, fvPatchField, volMesh > volSymmTensorField
Definition: volFieldsFwd.H:58
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
virtual bool read(const dictionary &XiEqProperties)=0
Update properties from given dictionary.
tmp< volScalarField > calculateSchelkinEffect(const scalar) const
Return the sub-grid Schelkin effect.
dimensioned< symmTensor > dimensionedSymmTensor
Dimensioned tensor obtained from generic dimensioned type.
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
psiReactionThermo & thermo
Definition: createFields.H:32
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A class for managing temporary objects.
Definition: PtrList.H:118
defineTypeNameAndDebug(combustionModel, 0)