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
32Info<< "Creating phaseChangeTwoPhaseMixture\n" << endl;
33autoPtr<phaseChangeTwoPhaseMixture> mixture =
34 phaseChangeTwoPhaseMixture::New(U, phi);
35
36volScalarField& alpha1(mixture->alpha1());
37volScalarField& alpha2(mixture->alpha2());
38
39const dimensionedScalar& rho1 = mixture->rho1();
40const dimensionedScalar& rho2 = mixture->rho2();
41
42
43// Need to store rho for ddt(rho, U)
44volScalarField rho
45(
46 IOobject
47 (
48 "rho",
49 runTime.timeName(),
50 mesh,
51 IOobject::READ_IF_PRESENT
52 ),
54);
55rho.oldTime();
56
57
58// Construct interface from alpha1 distribution
59interfaceProperties interface(alpha1, U, mixture());
60
61// Construct incompressible turbulence model
62autoPtr<incompressible::turbulenceModel> turbulence
63(
64 incompressible::turbulenceModel::New(U, phi, mixture())
65);
66
67
68#include "readGravitationalAcceleration.H"
69#include "readhRef.H"
70#include "gh.H"
71
72
73volScalarField p
74(
75 IOobject
76 (
77 "p",
78 runTime.timeName(),
79 mesh,
80 IOobject::NO_READ,
81 IOobject::AUTO_WRITE
82 ),
83 p_rgh + rho*gh
84);
85
86label pRefCell = 0;
87scalar pRefValue = 0.0;
89(
90 p,
91 p_rgh,
92 pimple.dict(),
95);
96
97if (p_rgh.needReference())
98{
99 p += dimensionedScalar
100 (
101 "p",
102 p.dimensions(),
103 pRefValue - getRefCellValue(p, pRefCell)
104 );
105 p_rgh = p - rho*gh;
106}
107
108mesh.setFluxRequired(p_rgh.name());
109mesh.setFluxRequired(alpha1.name());
110
111#include "createFvOptions.H"
112
114(
115 IOobject::groupName("alphaPhi0", alpha1.group()),
116 runTime.timeName(),
117 mesh,
118 IOobject::NO_READ,
119 IOobject::NO_WRITE
120);
121
122// MULES flux from previous time-step
123surfaceScalarField alphaPhi10
124(
126 phi*fvc::interpolate(alpha1)
127);
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
U
Definition: pEqn.H:72
volScalarField & p
alphaPhi10
Definition: alphaEqn.H:7
dynamicFvMesh & mesh
engineTime & runTime
compressible::turbulenceModel & turbulence
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())
Info<< "Creating temperaturePhaseChangeTwoPhaseMixture\n"<< endl;autoPtr< temperaturePhaseChangeTwoPhaseMixture > mixture
Definition: createFields.H:39
IOobject alphaPhi10Header(IOobject::groupName("alphaPhi0", alpha1.group()), runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE)