createFields.H
Go to the documentation of this file.
1#include "createRDeltaT.H"
2
3#include "readGravitationalAcceleration.H"
4
5Info<< "Reading thermophysical properties\n" << endl;
6autoPtr<psiReactionThermo> pThermo(psiReactionThermo::New(mesh));
7psiReactionThermo& thermo = pThermo();
8thermo.validate(args.executable(), "h", "e");
9
10SLGThermo slgThermo(mesh, thermo);
11
12basicSpecieMixture& composition = thermo.composition();
13PtrList<volScalarField>& Y = composition.Y();
14
15const word inertSpecie(thermo.get<word>("inertSpecie"));
16if (!composition.species().found(inertSpecie))
17{
19 << "Inert specie " << inertSpecie << " not found in available species "
20 << composition.species()
21 << exit(FatalIOError);
22}
23
24volScalarField& p = thermo.p();
25
26multivariateSurfaceInterpolationScheme<scalar>::fieldTable fields;
27
29{
30 fields.add(Y[i]);
31}
32fields.add(thermo.he());
33
34volScalarField rho
35(
36 IOobject
37 (
38 "rho",
39 runTime.timeName(),
40 mesh,
41 IOobject::NO_READ,
42 IOobject::AUTO_WRITE
43 ),
44 thermo.rho()
45);
46
47// lagrangian effective density field - used externally (optional)
48volScalarField rhoEffLagrangian
49(
50 IOobject
51 (
52 "rhoEffLagrangian",
53 runTime.timeName(),
54 mesh,
55 IOobject::NO_READ,
56 IOobject::AUTO_WRITE
57 ),
58 mesh,
59 dimensionedScalar(dimDensity, Zero)
60);
61
62// dynamic pressure field - used externally (optional)
63volScalarField pDyn
64(
65 IOobject
66 (
67 "pDyn",
68 runTime.timeName(),
69 mesh,
70 IOobject::NO_READ,
71 IOobject::AUTO_WRITE
72 ),
73 mesh,
74 dimensionedScalar(dimPressure, Zero)
75);
76
77
78Info<< "\nReading field U\n" << endl;
79volVectorField U
80(
81 IOobject
82 (
83 "U",
84 runTime.timeName(),
85 mesh,
86 IOobject::MUST_READ,
87 IOobject::AUTO_WRITE
88 ),
89 mesh
90);
91
93
94mesh.setFluxRequired(p.name());
95
96Info<< "Creating turbulence model\n" << endl;
97autoPtr<compressible::turbulenceModel> turbulence
98(
99 compressible::turbulenceModel::New
100 (
101 rho,
102 U,
103 phi,
104 thermo
105 )
106);
107
108Info<< "Creating combustion model\n" << endl;
109autoPtr<CombustionModel<psiReactionThermo>> combustion
110(
111 CombustionModel<psiReactionThermo>::New(thermo, turbulence())
112);
113
114volScalarField Qdot
115(
116 IOobject
117 (
118 "Qdot",
119 runTime.timeName(),
120 mesh,
121 IOobject::READ_IF_PRESENT,
122 IOobject::AUTO_WRITE
123 ),
124 mesh,
125 dimensionedScalar(dimEnergy/dimVolume/dimTime, Zero)
126);
127
128#include "createDpdt.H"
129
130#include "createK.H"
131
132#include "createMRF.H"
133#include "createClouds.H"
134#include "createRadiationModel.H"
135#include "createFvOptions.H"
surfaceScalarField & phi
const word & executable() const noexcept
Name of executable without the path.
Definition: argListI.H:51
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
U
Definition: pEqn.H:72
volScalarField & p
Creates and initialises the face-flux field phi.
dynamicFvMesh & mesh
engineTime & runTime
#define FatalIOErrorIn(functionName, ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:467
compressible::turbulenceModel & turbulence
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
Foam::argList args(argc, argv)
scalar Qdot
Definition: solveChemistry.H:2
Info<< "Reading thermophysical properties\n"<< endl;autoPtr< psiReactionThermo > pThermo(psiReactionThermo::New(mesh))
const word inertSpecie(thermo.get< word >("inertSpecie"))
basicSpecieMixture & composition
Definition: createFields.H:8
SLGThermo slgThermo(mesh, thermo)
Info<< "Creating combustion model\n"<< endl;autoPtr< CombustionModel< psiReactionThermo > > combustion(CombustionModel< psiReactionThermo >::New(thermo, turbulence()))
PtrList< volScalarField > & Y
Definition: createFields.H:9
multivariateSurfaceInterpolationScheme< scalar >::fieldTable fields
Definition: createFields.H:97
volScalarField rhoEffLagrangian(IOobject("rhoEffLagrangian", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), mesh, dimensionedScalar(dimDensity, Zero))
volScalarField pDyn(IOobject("pDyn", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), mesh, dimensionedScalar(dimPressure, Zero))
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333