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
31Info<< "Reading transportProperties\n" << endl;
32incompressibleTwoPhaseMixture mixture(U, phi);
33
36
37const dimensionedScalar& rho1 = mixture.rho1();
38const dimensionedScalar& rho2 = mixture.rho2();
39
40dimensionedScalar Dab("Dab", dimViscosity, mixture);
41
42// Read the reciprocal of the turbulent Schmidt number
43dimensionedScalar alphatab("alphatab", dimless, mixture);
44
45// Need to store rho for ddt(rho, U)
47rho.oldTime();
48
49
50// Mass flux
51// Initialisation does not matter because rhoPhi is reset after the
52// alpha1 solution before it is used in the U equation.
54(
55 IOobject
56 (
57 "rhoPhi",
58 runTime.timeName(),
59 mesh,
60 IOobject::NO_READ,
61 IOobject::NO_WRITE
62 ),
63 rho1*phi
64);
65
66// Construct incompressible turbulence model
67autoPtr<incompressible::turbulenceModel> turbulence
68(
69 incompressible::turbulenceModel::New(U, phi, mixture)
70);
71
72
73#include "readGravitationalAcceleration.H"
74#include "readhRef.H"
75#include "gh.H"
76
77
78volScalarField p
79(
80 IOobject
81 (
82 "p",
83 runTime.timeName(),
84 mesh,
85 IOobject::NO_READ,
86 IOobject::AUTO_WRITE
87 ),
88 p_rgh + rho*gh
89);
90
91label pRefCell = 0;
92scalar pRefValue = 0.0;
94(
95 p,
96 p_rgh,
97 pimple.dict(),
100);
101
102if (p_rgh.needReference())
103{
104 p += dimensionedScalar
105 (
106 "p",
107 p.dimensions(),
108 pRefValue - getRefCellValue(p, pRefCell)
109 );
110 p_rgh = p - rho*gh;
111}
112
113mesh.setFluxRequired(p_rgh.name());
114mesh.setFluxRequired(alpha1.name());
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