generalizedNewtonian.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2018 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "generalizedNewtonian.H"
27 #include "volFields.H"
28 #include "surfaceFields.H"
29 #include "fvcGrad.H"
30 #include "fvcDiv.H"
31 #include "fvmLaplacian.H"
32 
33 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37 namespace laminarModels
38 {
39 
40 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
41 
42 template<class BasicTurbulenceModel>
44 (
45  const alphaField& alpha,
46  const rhoField& rho,
47  const volVectorField& U,
48  const surfaceScalarField& alphaRhoPhi,
49  const surfaceScalarField& phi,
50  const transportModel& transport,
51  const word& propertiesName
52 )
53 :
55  (
56  typeName,
57  alpha,
58  rho,
59  U,
60  alphaRhoPhi,
61  phi,
62  transport,
63  propertiesName
64  ),
65 
66  viscosityModel_
67  (
69  (
70  this->coeffDict_
71  )
72  ),
73 
74  nu_
75  (
76  IOobject
77  (
78  IOobject::groupName("generalizedNewtonian:nu", alphaRhoPhi.group()),
79  this->runTime_.timeName(),
80  this->mesh_,
83  ),
84  viscosityModel_->nu(this->nu(), strainRate())
85  )
86 {}
87 
88 
89 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
90 
91 template<class BasicTurbulenceModel>
94 {
95  return sqrt(2.0)*mag(symm(fvc::grad(this->U())));
96 }
97 
98 
99 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
100 
101 template<class BasicTurbulenceModel>
103 {
104  viscosityModel_->read(this->coeffDict_);
105 
106  return true;
107 }
108 
109 
110 template<class BasicTurbulenceModel>
113 {
114  return volScalarField::New
115  (
116  IOobject::groupName("nut", this->alphaRhoPhi_.group()),
117  this->mesh_,
119  );
120 }
121 
122 
123 template<class BasicTurbulenceModel>
126 (
127  const label patchi
128 ) const
129 {
130  return tmp<scalarField>
131  (
132  new scalarField(this->mesh_.boundary()[patchi].size(), 0.0)
133  );
134 }
135 
136 
137 template<class BasicTurbulenceModel>
140 {
141  return volScalarField::New
142  (
143  IOobject::groupName("nuEff", this->alphaRhoPhi_.group()),
144  nu_
145  );
146 }
147 
148 
149 template<class BasicTurbulenceModel>
152 (
153  const label patchi
154 ) const
155 {
156  return nu_.boundaryField()[patchi];
157 }
158 
159 
160 template<class BasicTurbulenceModel>
163 {
164  return volScalarField::New
165  (
166  IOobject::groupName("k", this->alphaRhoPhi_.group()),
167  this->mesh_,
168  dimensionedScalar(sqr(this->U_.dimensions()), 0)
169  );
170 }
171 
172 
173 template<class BasicTurbulenceModel>
176 {
177  return volScalarField::New
178  (
179  IOobject::groupName("epsilon", this->alphaRhoPhi_.group()),
180  this->mesh_,
181  dimensionedScalar(sqr(this->U_.dimensions())/dimTime, 0)
182  );
183 }
184 
185 
186 template<class BasicTurbulenceModel>
189 {
191  (
192  IOobject::groupName("R", this->alphaRhoPhi_.group()),
193  this->mesh_,
194  dimensionedSymmTensor(sqr(this->U_.dimensions()), Zero)
195  );
196 }
197 
198 
199 template<class BasicTurbulenceModel>
201 {
202  nu_ = viscosityModel_->nu(this->nu(), strainRate());
204 }
205 
206 
207 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
208 
209 } // End namespace laminarModels
210 } // End namespace Foam
211 
212 // ************************************************************************* //
Foam::surfaceFields.
static autoPtr< generalizedNewtonianViscosityModel > New(const dictionary &viscosityProperties)
Return a reference to the selected viscosity model.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
BasicTurbulenceModel::rhoField rhoField
surfaceScalarField & phi
virtual void correct()
Correct the generalizedNewtonian viscosity.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionSet dimViscosity
dimensionedScalar sqrt(const dimensionedScalar &ds)
Calculate the matrix for the laplacian of the field.
Linear viscous stress turbulence model base class.
virtual tmp< volScalarField > k() const
Return the turbulence kinetic energy.
generalizedNewtonian(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName=turbulenceModel::propertiesName)
Construct from components.
static tmp< GeometricField< scalar, fvPatchField, volMesh > > New(const word &name, const Mesh &, const dimensionSet &, const word &patchFieldType=fvPatchField< scalar >::calculatedType())
Return a temporary field constructed from name, mesh, dimensionSet.
virtual tmp< volScalarField > epsilon() const
Return the turbulence kinetic energy dissipation rate,.
Calculate the gradient of the given field.
A class for handling words, derived from string.
Definition: word.H:59
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
static word groupName(Name name, const word &group)
static const zero Zero
Definition: zero.H:97
virtual bool read()
Read turbulenceProperties dictionary.
Calculate the divergence of the given field.
virtual tmp< volScalarField > strainRate() const
virtual tmp< volScalarField > nut() const
Return the turbulence viscosity,.
dimensioned< symmTensor > dimensionedSymmTensor
Dimensioned tensor obtained from generic dimensioned type.
BasicTurbulenceModel::alphaField alphaField
label patchi
U
Definition: pEqn.H:72
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
virtual void correct()
Correct the laminar transport.
Definition: laminarModel.C:274
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
virtual tmp< volScalarField > nuEff() const
Return the effective viscosity.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
BasicTurbulenceModel::transportModel transportModel
dimensioned< scalar > mag(const dimensioned< Type > &)
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
A class for managing temporary objects.
Definition: PtrList.H:53
virtual tmp< volSymmTensorField > R() const
Return the Reynolds stress tensor.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
volScalarField & nu
Namespace for OpenFOAM.