transport.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 "transport.H"
27 #include "surfaceInterpolate.H"
28 #include "fvmDdt.H"
29 #include "fvcLaplacian.H"
30 #include "fvmDiv.H"
31 #include "fvmSup.H"
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38 namespace XiModels
39 {
40  defineTypeNameAndDebug(transport, 0);
41  addToRunTimeSelectionTable(XiModel, transport, dictionary);
42 }
43 }
44 
45 
46 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
47 
48 Foam::XiModels::transport::transport
49 (
50  const dictionary& XiProperties,
51  const psiuReactionThermo& thermo,
53  const volScalarField& Su,
54  const volScalarField& rho,
55  const volScalarField& b,
56  const surfaceScalarField& phi
57 )
58 :
59  XiModel(XiProperties, thermo, turbulence, Su, rho, b, phi),
60  XiShapeCoef(readScalar(XiModelCoeffs_.lookup("XiShapeCoef"))),
61  XiEqModel_(XiEqModel::New(XiProperties, thermo, turbulence, Su)),
62  XiGModel_(XiGModel::New(XiProperties, thermo, turbulence, Su))
63 {}
64 
65 
66 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
67 
69 {}
70 
71 
72 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
73 
75 {
76  return XiGModel_->Db();
77 }
78 
79 
81 (
82  const fv::convectionScheme<scalar>& mvConvection
83 )
84 {
85  volScalarField XiEqEta(XiEqModel_->XiEq());
86  volScalarField GEta(XiGModel_->G());
87 
88  volScalarField R(GEta*XiEqEta/(XiEqEta - 0.999));
89 
90  volScalarField XiEqStar(R/(R - GEta));
91 
92  volScalarField XiEq
93  (
94  1.0 + (1.0 + (2*XiShapeCoef)*(0.5 - b_))*(XiEqStar - 1.0)
95  );
96 
97  volScalarField G(R*(XiEq - 1.0)/XiEq);
98 
99  const objectRegistry& db = b_.db();
100  const volScalarField& betav = db.lookupObject<volScalarField>("betav");
101  const volScalarField& mgb = db.lookupObject<volScalarField>("mgb");
102  const surfaceScalarField& phiSt =
103  db.lookupObject<surfaceScalarField>("phiSt");
104  const volScalarField& Db = db.lookupObject<volScalarField>("Db");
105  const surfaceScalarField& nf = db.lookupObject<surfaceScalarField>("nf");
106 
107  surfaceScalarField phiXi
108  (
109  "phiXi",
110  phiSt
111  + (
112  - fvc::interpolate(fvc::laplacian(Db, b_)/mgb)*nf
114  )
115  );
116 
117  solve
118  (
119  betav*fvm::ddt(rho_, Xi_)
120  + mvConvection.fvmDiv(phi_, Xi_)
121  + fvm::div(phiXi, Xi_)
122  - fvm::Sp(fvc::div(phiXi), Xi_)
123  ==
124  betav*rho_*R
125  - fvm::Sp(betav*rho_*(R - G), Xi_)
126  );
127 
128  // Correct boundedness of Xi
129  // ~~~~~~~~~~~~~~~~~~~~~~~~~
130  Xi_.max(1.0);
131  Xi_ = min(Xi_, 2.0*XiEq);
132 }
133 
134 
135 bool Foam::XiModels::transport::read(const dictionary& XiProperties)
136 {
137  XiModel::read(XiProperties);
138 
139  XiModelCoeffs_.lookup("XiShapeCoef") >> XiShapeCoef;
140 
141  return true;
142 }
143 
144 
145 // ************************************************************************* //
autoPtr< compressible::turbulenceModel > turbulence
Definition: createFields.H:23
surfaceScalarField & phi
const volScalarField & b_
Definition: XiModel.H:121
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:221
tmp< fvMatrix< Type > > Sp(const DimensionedField< scalar, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &)
const dimensionedScalar G
Newtonian constant of gravitation.
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
volScalarField Xi_
Flame wrinking field.
Definition: XiModel.H:125
const surfaceScalarField & phi_
Definition: XiModel.H:122
tmp< GeometricField< Type, fvPatchField, volMesh > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcLaplacian.C:45
Macros for easy insertion into run-time selection tables.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
psiReactionThermo & thermo
Definition: createFields.H:31
stressControl lookup("compactNormalStress") >> compactNormalStress
virtual bool read(const dictionary &XiProperties)=0
Update properties from given dictionary.
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
Definition: createFields.H:27
const volScalarField & rho_
Definition: XiModel.H:120
SolverPerformance< Type > solve(fvMatrix< Type > &, const dictionary &)
Solve returning the solution statistics given convergence tolerance.
Calculate the laplacian of the given field.
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvmDdt.C:46
Calulate the matrix for the first temporal derivative.
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if succesful.
Definition: doubleScalar.H:63
const volScalarField & Su_
Definition: XiModel.H:119
dictionary XiModelCoeffs_
Definition: XiModel.H:115
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmDiv.C:46
defineTypeNameAndDebug(combustionModel, 0)
RASModel< EddyDiffusivity< turbulenceModel > > RASModel
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
Calculate the matrix for the divergence of the given field and flux.
virtual void correct()
Correct the flame-wrinking Xi.
Definition: transport.H:116
#define R(A, B, C, D, E, F, K, M)
virtual bool read(const dictionary &XiProperties)
Update properties from given dictionary.
void max(const dimensioned< Type > &)
virtual tmp< volScalarField > Db() const
Return the flame diffusivity.
A class for managing temporary objects.
Definition: PtrList.H:54
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
tmp< GeometricField< Type, fvPatchField, volMesh > > Su(const GeometricField< Type, fvPatchField, volMesh > &su, const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvcSup.C:44
virtual ~transport()
Destructor.
Calculate the matrix for implicit and explicit sources.
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:451