createFields.H
Go to the documentation of this file.
1Info<< "Reading thermophysical properties\n" << endl;
2
3Info<< "Reading field T\n" << endl;
5(
6 IOobject
7 (
8 "T",
9 runTime.timeName(),
10 mesh,
11 IOobject::MUST_READ,
12 IOobject::AUTO_WRITE
13 ),
14 mesh
15);
16
17Info<< "Reading field p_rgh\n" << endl;
19(
20 IOobject
21 (
22 "p_rgh",
23 runTime.timeName(),
24 mesh,
25 IOobject::MUST_READ,
26 IOobject::AUTO_WRITE
27 ),
28 mesh
29);
30
31Info<< "Reading field U\n" << endl;
33(
34 IOobject
35 (
36 "U",
37 runTime.timeName(),
38 mesh,
39 IOobject::MUST_READ,
40 IOobject::AUTO_WRITE
41 ),
42 mesh
43);
44
45#include "createPhi.H"
46
48
49Info<< "Creating turbulence model\n" << endl;
50autoPtr<incompressible::turbulenceModel> turbulence
51(
52 incompressible::turbulenceModel::New(U, phi, laminarTransport)
53);
54
55// Kinematic density for buoyancy force
56volScalarField rhok
57(
58 IOobject
59 (
60 "rhok",
61 runTime.timeName(),
62 mesh
63 ),
64 1.0 - beta*(T - TRef)
65);
66
67// kinematic turbulent thermal thermal conductivity m2/s
68Info<< "Reading field alphat\n" << endl;
69volScalarField alphat
70(
71 IOobject
72 (
73 "alphat",
74 runTime.timeName(),
75 mesh,
76 IOobject::MUST_READ,
77 IOobject::AUTO_WRITE
78 ),
79 mesh
80);
81
82
83#include "readGravitationalAcceleration.H"
84#include "readhRef.H"
85#include "gh.H"
86
87
88volScalarField p
89(
90 IOobject
91 (
92 "p",
93 runTime.timeName(),
94 mesh,
95 IOobject::NO_READ,
96 IOobject::AUTO_WRITE
97 ),
98 p_rgh + rhok*gh
99);
100
101label pRefCell = 0;
102scalar pRefValue = 0.0;
104(
105 p,
106 p_rgh,
107 simple.dict(),
108 pRefCell,
110);
111
112if (p_rgh.needReference())
113{
114 p += dimensionedScalar
115 (
116 "p",
117 p.dimensions(),
118 pRefValue - getRefCellValue(p, pRefCell)
119 );
120}
121
122mesh.setFluxRequired(p_rgh.name());
123
124#include "createMRF.H"
126#include "createFvOptions.H"
volScalarField & p_rgh
surfaceScalarField & phi
const scalar pRefValue
const label pRefCell
const volScalarField & gh
U
Definition: pEqn.H:72
volScalarField & p
const volScalarField & T
dynamicFvMesh & mesh
engineTime & runTime
compressible::turbulenceModel & turbulence
rhok
Definition: TEqn.H:27
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
const dictionary & simple
dimensionedScalar TRef("TRef", dimTemperature, laminarTransport)
dimensionedScalar beta("beta", dimless/dimTemperature, laminarTransport)
singlePhaseTransportModel laminarTransport(U, phi)
setRefCell(p, pimple.dict(), pRefCell, pRefValue)