createFields.H
Go to the documentation of this file.
1Info<< "Reading velocity field U\n" << endl;
2volVectorField U
3(
4 IOobject
5 (
6 "U",
7 runTime.timeName(),
8 mesh,
9 IOobject::MUST_READ,
10 IOobject::AUTO_WRITE
11 ),
12 mesh
13);
14
15// Initialise the velocity internal field to zero
16U = dimensionedVector(U.dimensions(), Zero);
17
18surfaceScalarField phi
19(
20 IOobject
21 (
22 "phi",
23 runTime.timeName(),
24 mesh,
25 IOobject::NO_READ,
26 IOobject::AUTO_WRITE
27 ),
28 fvc::flux(U)
29);
30
31if (args.found("initialiseUBCs"))
32{
33 U.correctBoundaryConditions();
34 phi = fvc::flux(U);
35}
36
37
38// Construct a pressure field
39// If it is available read it otherwise construct from the velocity BCs
40// converting fixed-value BCs to zero-gradient and vice versa.
41
42// Allow override from command-line -pName option
43const word pName = args.getOrDefault<word>("pName", "p");
44
45// Infer the pressure BCs from the velocity
46wordList pBCTypes
47(
48 U.boundaryField().size(),
50);
51
52forAll(U.boundaryField(), patchi)
53{
54 if (U.boundaryField()[patchi].fixesValue())
55 {
57 }
58}
59
60Info<< "Constructing pressure field " << pName << nl << endl;
62(
63 IOobject
64 (
65 pName,
66 runTime.timeName(),
67 mesh,
68 IOobject::READ_IF_PRESENT,
69 IOobject::NO_WRITE
70 ),
71 mesh,
72 dimensionedScalar(sqr(dimVelocity), Zero),
74);
75
76// Infer the velocity potential BCs from the pressure
77wordList PhiBCTypes
78(
79 p.boundaryField().size(),
81);
82
83forAll(p.boundaryField(), patchi)
84{
85 if (p.boundaryField()[patchi].fixesValue())
86 {
87 PhiBCTypes[patchi] = fixedValueFvPatchScalarField::typeName;
88 }
89}
90
91Info<< "Constructing velocity potential field Phi\n" << endl;
93(
94 IOobject
95 (
96 "Phi",
97 runTime.timeName(),
98 mesh,
99 IOobject::READ_IF_PRESENT,
100 IOobject::NO_WRITE
101 ),
102 mesh,
103 dimensionedScalar(dimLength*dimVelocity, Zero),
104 PhiBCTypes
105);
106
107label PhiRefCell = 0;
108scalar PhiRefValue = 0;
110(
111 Phi,
112 potentialFlow.dict(),
113 PhiRefCell,
114 PhiRefValue
115);
116mesh.setFluxRequired(Phi.name());
117
118#include "createMRF.H"
surfaceScalarField & phi
bool found(const word &optName) const
Return true if the named option is found.
Definition: argListI.H:178
T getOrDefault(const word &optName, const T &deflt) const
Get a value from the named option if present, or return default.
Definition: argListI.H:307
U
Definition: pEqn.H:72
volScalarField & p
dynamicFvMesh & mesh
engineTime & runTime
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
List< word > wordList
A List of words.
Definition: fileName.H:63
dimensionedSymmTensor sqr(const dimensionedVector &dv)
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
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
Foam::argList args(argc, argv)
const dictionary & potentialFlow(mesh.solutionDict().subDict("potentialFlow"))
wordList pBCTypes(U.boundaryField().size(), fixedValueFvPatchScalarField::typeName)
setRefCell(p, pimple.dict(), pRefCell, pRefValue)
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
static const char *const typeName
The type name used in ensight case files.