createFields.H
Go to the documentation of this file.
1Info<< "Reading field p_rgh\n" << endl;
3(
4 IOobject
5 (
6 "p_rgh",
7 runTime.timeName(),
8 mesh,
9 IOobject::MUST_READ,
10 IOobject::AUTO_WRITE
11 ),
12 mesh
13);
14
15Info<< "Reading field U\n" << endl;
17(
18 IOobject
19 (
20 "U",
21 runTime.timeName(),
22 mesh,
23 IOobject::MUST_READ,
24 IOobject::AUTO_WRITE
25 ),
26 mesh
27);
28
29#include "createPhi.H"
30
31// Creating e based thermo
32autoPtr<twoPhaseMixtureEThermo> thermo
33(
34 new twoPhaseMixtureEThermo(U, phi)
35);
36
37// Create mixture and
38Info<< "Creating temperaturePhaseChangeTwoPhaseMixture\n" << endl;
39autoPtr<temperaturePhaseChangeTwoPhaseMixture> mixture =
40 temperaturePhaseChangeTwoPhaseMixture::New(thermo(), mesh);
41
42
43volScalarField& alpha1(thermo->alpha1());
44volScalarField& alpha2(thermo->alpha2());
45
46const dimensionedScalar& rho1 = thermo->rho1();
47const dimensionedScalar& rho2 = thermo->rho2();
48
49// Need to store rho for ddt(rho, U)
50volScalarField rho
51(
52 IOobject
53 (
54 "rho",
55 runTime.timeName(),
56 mesh,
57 IOobject::READ_IF_PRESENT,
58 IOobject::AUTO_WRITE
59 ),
61);
62rho.oldTime();
63
64// Mass flux
65surfaceScalarField rhoPhi
66(
67 IOobject
68 (
69 "rhoPhi",
70 runTime.timeName(),
71 mesh,
72 IOobject::NO_READ,
73 IOobject::NO_WRITE
74 ),
75 fvc::interpolate(rho)*phi
76);
77
78// Construct interface from alpha1 distribution
79interfaceProperties interface
80(
81 alpha1,
82 U,
83 thermo->transportPropertiesDict()
84);
85
86// Construct incompressible turbulence model
87autoPtr<incompressible::turbulenceModel> turbulence
88(
89 incompressible::turbulenceModel::New(U, phi, thermo())
90);
91
92#include "readGravitationalAcceleration.H"
93#include "readhRef.H"
94#include "gh.H"
95
96volScalarField& p = thermo->p();
97
98label pRefCell = 0;
99scalar pRefValue = 0.0;
101(
102 p,
103 p_rgh,
104 pimple.dict(),
105 pRefCell,
107);
108
109if (p_rgh.needReference())
110{
111 p += dimensionedScalar
112 (
113 "p",
114 p.dimensions(),
115 pRefValue - getRefCellValue(p, pRefCell)
116 );
117 p_rgh = p - rho*gh;
118}
119
120mesh.setFluxRequired(p_rgh.name());
121mesh.setFluxRequired(alpha1.name());
122
123#include "createMRF.H"
124#include "createFvOptions.H"
125
126// Turbulent Prandtl number
127dimensionedScalar Prt("Prt", dimless, thermo->transportPropertiesDict());
128
129volScalarField kappaEff
130(
131 IOobject
132 (
133 "kappaEff",
134 runTime.timeName(),
135 mesh,
136 IOobject::NO_READ,
137 IOobject::NO_WRITE
138 ),
139 thermo->kappa()
140);
141
142// Need to store rho for ddt(rhoCp, U)
143volScalarField rhoCp
144(
145 IOobject
146 (
147 "rhoCp",
148 runTime.timeName(),
149 mesh,
150 IOobject::NO_READ,
151 IOobject::NO_WRITE
152 ),
153 rho*thermo->Cp()
154);
155
156rhoCp.oldTime();
rhoPhi
Definition: rhoEqn.H:10
volScalarField & p_rgh
surfaceScalarField & phi
const scalar pRefValue
const label pRefCell
const volScalarField & gh
pimpleControl & pimple
const volScalarField & alpha1
volScalarField & rho2
const volScalarField & alpha2
volScalarField & rho1
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
dynamicFvMesh & mesh
engineTime & runTime
compressible::turbulenceModel & turbulence
rhoCp
Definition: TEqn.H:3
kappaEff
Definition: TEqn.H:10
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:83
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
setRefCell(p, pimple.dict(), pRefCell, pRefValue)
interfaceProperties interface(alpha1, U, thermo->transportPropertiesDict())
dimensionedScalar Prt("Prt", dimless, thermo->transportPropertiesDict())
Info<< "Creating temperaturePhaseChangeTwoPhaseMixture\n"<< endl;autoPtr< temperaturePhaseChangeTwoPhaseMixture > mixture
Definition: createFields.H:39