createFields.H
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2015 OpenFOAM Foundation
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
28Info<< "Reading field U\n" << endl;
29volVectorField U
30(
31 IOobject
32 (
33 "U",
34 runTime.timeName(),
35 mesh,
36 IOobject::MUST_READ,
37 IOobject::NO_WRITE
38 ),
39 mesh
40);
41
42Info<< "Calculating wall distance field" << endl;
43volScalarField y
44(
45 IOobject
46 (
47 "y",
48 runTime.timeName(),
49 mesh,
50 IOobject::NO_READ,
51 IOobject::NO_WRITE
52 ),
53 mesh,
54 dimensionedScalar(dimLength, Zero),
56);
57y.primitiveFieldRef() = wallDist::New(mesh).y().primitiveField();
58y.correctBoundaryConditions();
59
60
61// Set the mean boundary-layer thickness
62dimensionedScalar ybl("ybl", dimLength, Zero);
63
64if (args.found("ybl"))
65{
66 // If the boundary-layer thickness is provided use it
67 ybl.value() = args.get<scalar>("ybl");
68}
69else if (args.found("Cbl"))
70{
71 // Calculate boundary layer thickness as Cbl*mean distance to wall
72 ybl.value() = gAverage(y)*args.get<scalar>("Cbl");
73}
74
75Info<< "\nCreating boundary-layer for U of thickness "
76 << ybl.value() << " m" << nl << endl;
77
78Info<< "Creating mask field" << endl;
80(
81 IOobject
82 (
83 "mask",
84 runTime.timeName(),
85 mesh,
86 IOobject::NO_READ,
87 IOobject::NO_WRITE
88 ),
89 mesh,
90 dimensionedScalar(dimless, Zero),
92);
93
94
95// ************************************************************************* //
scalar y
T get(const label index) const
Get a value from the argument at index.
Definition: argListI.H:278
bool found(const word &optName) const
Return true if the named option is found.
Definition: argListI.H:178
U
Definition: pEqn.H:72
dynamicFvMesh & mesh
engineTime & runTime
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
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)
static const char *const typeName
The type name used in ensight case files.