createFields.H
Go to the documentation of this file.
1#include "createRDeltaT.H"
2
3Info<< "Reading field p_rgh\n" << endl;
5(
6 IOobject
7 (
8 "p_rgh",
9 runTime.timeName(),
10 mesh,
11 IOobject::MUST_READ,
12 IOobject::AUTO_WRITE
13 ),
14 mesh
15);
16
17Info<< "Reading field U\n" << endl;
19(
20 IOobject
21 (
22 "U",
23 runTime.timeName(),
24 mesh,
25 IOobject::MUST_READ,
26 IOobject::AUTO_WRITE
27 ),
28 mesh
29);
30
31#include "createPhi.H"
32
33immiscibleIncompressibleThreePhaseMixture mixture(U, phi);
34
37volScalarField& alpha3(mixture.alpha3());
38
39const dimensionedScalar& rho1 = mixture.rho1();
40const dimensionedScalar& rho2 = mixture.rho2();
41const dimensionedScalar& rho3 = mixture.rho3();
42
43dimensionedScalar D23("D23", dimViscosity, mixture);
44
45// Need to store rho for ddt(rho, U)
47(
48 IOobject
49 (
50 "rho",
51 runTime.timeName(),
52 mesh,
53 IOobject::READ_IF_PRESENT
54 ),
55 alpha1*rho1 + alpha2*rho2 + alpha3*rho3
56);
57rho.oldTime();
58
59
60// Mass flux
61// Initialisation does not matter because rhoPhi is reset after the
62// alpha solution before it is used in the U equation.
64(
65 IOobject
66 (
67 "rhoPhi",
68 runTime.timeName(),
69 mesh,
70 IOobject::NO_READ,
71 IOobject::NO_WRITE
72 ),
73 rho1*phi
74);
75
76
77// Construct incompressible turbulence model
78autoPtr<incompressible::turbulenceModel> turbulence
79(
80 incompressible::turbulenceModel::New(U, phi, mixture)
81);
82
83
84#include "readGravitationalAcceleration.H"
85#include "readhRef.H"
86#include "gh.H"
87
88
89volScalarField p
90(
91 IOobject
92 (
93 "p",
94 runTime.timeName(),
95 mesh,
96 IOobject::NO_READ,
97 IOobject::AUTO_WRITE
98 ),
99 p_rgh + rho*gh
100);
101
102label pRefCell = 0;
103scalar pRefValue = 0.0;
105(
106 p,
107 p_rgh,
108 pimple.dict(),
109 pRefCell,
111);
112
113if (p_rgh.needReference())
114{
115 p += dimensionedScalar
116 (
117 "p",
118 p.dimensions(),
119 pRefValue - getRefCellValue(p, pRefCell)
120 );
121 p_rgh = p - rho*gh;
122}
123
124mesh.setFluxRequired(p_rgh.name());
125mesh.setFluxRequired(alpha2.name());
126
127#include "createMRF.H"
128#include "createFvOptions.H"
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
U
Definition: pEqn.H:72
volScalarField & p
dynamicFvMesh & mesh
engineTime & runTime
compressible::turbulenceModel & turbulence
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
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)
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
setRefCell(p, pimple.dict(), pRefCell, pRefValue)
Info<< "Creating temperaturePhaseChangeTwoPhaseMixture\n"<< endl;autoPtr< temperaturePhaseChangeTwoPhaseMixture > mixture
Definition: createFields.H:39