CourantNo.C
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) 2013-2016 OpenFOAM Foundation
9 Copyright (C) 2020 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29#include "CourantNo.H"
30#include "surfaceFields.H"
31#include "fvcSurfaceIntegrate.H"
34
35// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
36
37namespace Foam
38{
39namespace functionObjects
40{
43}
44}
45
46
47// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
48
50Foam::functionObjects::CourantNo::byRho
51(
52 const tmp<volScalarField::Internal>& Co
53) const
54{
55 if (Co().dimensions() == dimDensity)
56 {
57 return Co/obr_.lookupObject<volScalarField>(rhoName_);
58 }
59
60 return Co;
61}
62
63
64bool Foam::functionObjects::CourantNo::calc()
65{
66 if (foundObject<surfaceScalarField>(fieldName_))
67 {
68 const surfaceScalarField& phi =
69 lookupObject<surfaceScalarField>(fieldName_);
70
71 tmp<volScalarField::Internal> Coi
72 (
73 byRho
74 (
75 (0.5*mesh_.time().deltaT())
77 /mesh_.V()
78 )
79 );
80
81 if (foundObject<volScalarField>(resultName_, false))
82 {
83 volScalarField& Co =
84 lookupObjectRef<volScalarField>(resultName_);
85
86 Co.ref() = Coi();
87 Co.correctBoundaryConditions();
88 }
89 else
90 {
92 (
93 IOobject
94 (
95 resultName_,
96 mesh_.time().timeName(),
97 mesh_,
100 ),
101 mesh_,
103 zeroGradientFvPatchScalarField::typeName
104 );
105 tCo.ref().ref() = Coi();
106 tCo.ref().correctBoundaryConditions();
107 mesh_.objectRegistry::store(tCo.ptr());
108 }
109
110 return true;
111 }
112
113 return false;
114}
115
116
117// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
118
120(
121 const word& name,
122 const Time& runTime,
123 const dictionary& dict
124)
125:
127 rhoName_("rho")
128{
129 setResultName("Co", "phi");
130 read(dict);
131}
132
133
134// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
135
137{
139
140 rhoName_ = dict.getOrDefault<word>("rho", "rho");
141
142 return true;
143}
144
145
146// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
surfaceScalarField & phi
virtual bool read()
Re-read model coefficients if they have changed.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:80
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
Abstract base-class for Time/database function objects.
Computes the Courant number field for time-variant simulations.
Definition: CourantNo.H:148
virtual bool read(const dictionary &)
Read the CourantNo data.
Definition: CourantNo.C:136
Intermediate class for handling field expression function objects (e.g. blendingFactor etc....
void setResultName(const word &typeName, const word &defaultArg)
Set the name of result field.
const objectRegistry & obr_
Reference to the region objectRegistry.
const Type & lookupObject(const word &name, const bool recursive=false) const
A class for managing temporary objects.
Definition: tmp.H:65
A class for handling words, derived from Foam::string.
Definition: word.H:68
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
engineTime & runTime
Surface integrate surfaceField creating a volField. Surface sum a surfaceField creating a volField.
tmp< GeometricField< Type, fvPatchField, volMesh > > surfaceSum(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Namespace for OpenFOAM.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const dimensionSet dimless
Dimensionless.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
const dimensionSet dimDensity
dictionary dict
Foam::surfaceFields.