generalizedNewtonian.H
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 Class
25  Foam::laminarModels::generalizedNewtonian
26 
27 Description
28  Turbulence model for shear-dependent Non-Newtonian flow.
29 
30 SourceFiles
31  generalizedNewtonian.C
32 
33 \*---------------------------------------------------------------------------*/
34 
35 #ifndef generalizedNewtonian_H
36 #define generalizedNewtonian_H
37 
38 #include "laminarModel.H"
39 #include "linearViscousStress.H"
41 
42 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
43 
44 namespace Foam
45 {
46 namespace laminarModels
47 {
48 
49 /*---------------------------------------------------------------------------*\
50  Class generalizedNewtonian Declaration
51 \*---------------------------------------------------------------------------*/
52 
53 template<class BasicTurbulenceModel>
55 :
56  public linearViscousStress<laminarModel<BasicTurbulenceModel>>
57 {
58 
59 protected:
60 
61  // Protected data
62 
63  //- Run-time selectable non-Newtonian viscosity model
65 
66  //- The non-Newtonian viscosity field
68 
69 
70  // Protected Member Functions
71 
72  virtual tmp<volScalarField> strainRate() const;
73 
74 
75 public:
76 
77  typedef typename BasicTurbulenceModel::alphaField alphaField;
78  typedef typename BasicTurbulenceModel::rhoField rhoField;
79  typedef typename BasicTurbulenceModel::transportModel transportModel;
80 
81 
82  //- Runtime type information
83  TypeName("generalizedNewtonian");
84 
85 
86  // Constructors
87 
88  //- Construct from components
90  (
91  const alphaField& alpha,
92  const rhoField& rho,
93  const volVectorField& U,
94  const surfaceScalarField& alphaRhoPhi,
95  const surfaceScalarField& phi,
96  const transportModel& transport,
97  const word& propertiesName = turbulenceModel::propertiesName
98  );
99 
100 
101  // Selectors
102 
103  //- Return a reference to the selected turbulence model
105  (
106  const alphaField& alpha,
107  const rhoField& rho,
108  const volVectorField& U,
109  const surfaceScalarField& alphaRhoPhi,
110  const surfaceScalarField& phi,
111  const transportModel& transport,
112  const word& propertiesName = turbulenceModel::propertiesName
113  );
114 
115 
116  //- Destructor
117  virtual ~generalizedNewtonian()
118  {}
119 
120 
121  // Member Functions
122 
123  //- Read turbulenceProperties dictionary
124  virtual bool read();
125 
126  //- Return the turbulence viscosity,
127  // i.e. 0 for generalized Newtonian flow
128  virtual tmp<volScalarField> nut() const;
129 
130  //- Return the turbulence viscosity on patch
131  virtual tmp<scalarField> nut(const label patchi) const;
132 
133  //- Return the effective viscosity
134  // i.e. the generalizedNewtonian viscosity
135  virtual tmp<volScalarField> nuEff() const;
136 
137  //- Return the effective viscosity on patch
138  virtual tmp<scalarField> nuEff(const label patchi) const;
139 
140  //- Return the turbulence kinetic energy
141  // i.e. 0 for generalizedNewtonian flow
142  virtual tmp<volScalarField> k() const;
143 
144  //- Return the turbulence kinetic energy dissipation rate,
145  // i.e. 0 for generalizedNewtonian flow
146  virtual tmp<volScalarField> epsilon() const;
147 
148  //- Return the Reynolds stress tensor
149  // i.e. 0 for generalizedNewtonian flow
150  virtual tmp<volSymmTensorField> R() const;
151 
152  //- Correct the generalizedNewtonian viscosity
153  virtual void correct();
154 };
155 
156 
157 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
158 
159 } // End namespace laminarModels
160 } // End namespace Foam
161 
162 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
163 
164 #ifdef NoRepository
165  #include "generalizedNewtonian.C"
166 #endif
167 
168 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
169 
170 #endif
171 
172 // ************************************************************************* //
volScalarField nu_
The non-Newtonian viscosity field.
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.
Linear viscous stress turbulence model base class.
Turbulence model for shear-dependent Non-Newtonian flow.
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.
virtual tmp< volScalarField > epsilon() const
Return the turbulence kinetic energy dissipation rate,.
static const word propertiesName
Default name of the turbulence properties dictionary.
A class for handling words, derived from string.
Definition: word.H:59
static autoPtr< generalizedNewtonian > New(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName=turbulenceModel::propertiesName)
Return a reference to the selected turbulence model.
virtual bool read()
Read turbulenceProperties dictionary.
TypeName("generalizedNewtonian")
Runtime type information.
virtual tmp< volScalarField > strainRate() const
virtual tmp< volScalarField > nut() const
Return the turbulence viscosity,.
BasicTurbulenceModel::alphaField alphaField
label patchi
U
Definition: pEqn.H:72
virtual tmp< volScalarField > nuEff() const
Return the effective viscosity.
BasicTurbulenceModel::transportModel transportModel
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
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.
autoPtr< generalizedNewtonianViscosityModel > viscosityModel_
Run-time selectable non-Newtonian viscosity model.
Namespace for OpenFOAM.