basic.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) 2011-2016 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
28#include "basic.H"
30
31// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32
33namespace Foam
34{
35namespace PDRDragModels
36{
37 defineTypeNameAndDebug(basic, 0);
38 addToRunTimeSelectionTable(PDRDragModel, basic, dictionary);
39}
40}
41
42
43// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
44
46(
47 const dictionary& PDRProperties,
48 const compressible::RASModel& turbulence,
49 const volScalarField& rho,
50 const volVectorField& U,
51 const surfaceScalarField& phi
52)
53:
54 PDRDragModel(PDRProperties, turbulence, rho, U, phi),
55 Csu("Csu", dimless, PDRDragModelCoeffs_),
56 Csk("Csk", dimless, PDRDragModelCoeffs_),
57
58 Aw_
59 (
60 IOobject
61 (
62 "Aw",
63 U_.mesh().facesInstance(),
64 U_.mesh(),
65 IOobject::MUST_READ,
66 IOobject::NO_WRITE
67 ),
68 U_.mesh()
69 ),
70
71 CR_
72 (
73 IOobject
74 (
75 "CR",
76 U_.mesh().facesInstance(),
77 U_.mesh(),
78 IOobject::MUST_READ,
79 IOobject::NO_WRITE
80 ),
81 U_.mesh()
82 )
83{}
84
85
86// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
87
89{}
90
91
92// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
93
95{
96 tmp<volSymmTensorField> tDragDcu
97 (
99 (
100 IOobject
101 (
102 "tDragDcu",
103 U_.mesh().time().constant(),
104 U_.mesh(),
107 ),
108 U_.mesh(),
110 )
111 );
112
113 volSymmTensorField& DragDcu = tDragDcu.ref();
114
115 if (on_)
116 {
117 const volScalarField& betav =
118 U_.db().lookupObject<volScalarField>("betav");
119
120 DragDcu =
121 (0.5*rho_)*CR_*mag(U_) + (Csu*I)*betav*turbulence_.muEff()*sqr(Aw_);
122 }
123
124 return tDragDcu;
125}
126
127
129{
130 tmp<volScalarField> tGk
131 (
133 (
134 IOobject
135 (
136 "tGk",
137 U_.mesh().time().constant(),
138 U_.mesh(),
141 ),
142 U_.mesh(),
144 )
145 );
146
147 volScalarField& Gk = tGk.ref();
148
149 if (on_)
150 {
151 const volScalarField& betav =
152 U_.db().lookupObject<volScalarField>("betav");
153
154 const volSymmTensorField& CT =
156
157 Gk =
158 (0.5*rho_)*mag(U_)*(U_ & CT & U_)
159 + Csk*betav*turbulence_.muEff()*sqr(Aw_)*magSqr(U_);
160 }
161
162 return tGk;
163}
164
165
166bool Foam::PDRDragModels::basic::read(const dictionary& PDRProperties)
167{
168 PDRDragModel::read(PDRProperties);
169
170 PDRDragModelCoeffs_.readEntry("Csu", Csu.value());
171 PDRDragModelCoeffs_.readEntry("Csk", Csk.value());
172
173 return true;
174}
175
176
178{
179 Aw_.write();
180 CR_.write();
181}
182
183// ************************************************************************* //
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
Internal & ref(const bool updateAccessTime=true)
Return a reference to the dimensioned internal field.
const objectRegistry & db() const noexcept
Return the local objectRegistry.
Definition: IOobject.C:500
virtual bool read()
Inherit read from regIOobject.
Basic sub-grid obstacle drag model. Details supplied by J Puttock 2/7/06.
Definition: basic.H:100
virtual ~basic()
Destructor.
virtual tmp< volSymmTensorField > Dcu() const
Return the momentum drag coefficient.
virtual tmp< volScalarField > Gk() const
Return the momentum drag turbulence generation rate.
void writeFields() const
Write fields.
virtual bool read()
Re-read model coefficients if they have changed.
const Type & lookupObject(const word &name, const bool recursive=false) const
A class for managing temporary objects.
Definition: tmp.H:65
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
U
Definition: pEqn.H:72
dynamicFvMesh & mesh
compressible::turbulenceModel & turbulence
Namespace for OpenFOAM.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const dimensionSet dimless
Dimensionless.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:52
dimensionedScalar pow3(const dimensionedScalar &ds)
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
static const Identity< scalar > I
Definition: Identity.H:94
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
dimensioned< symmTensor > dimensionedSymmTensor
Dimensioned tensor obtained from generic dimensioned type.
GeometricField< symmTensor, fvPatchField, volMesh > volSymmTensorField
Definition: volFieldsFwd.H:86
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:60
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:51
const volScalarField & betav