linearViscousStress.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) 2013-2023 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 "linearViscousStress.H"
27 #include "fvcGrad.H"
28 #include "fvcDiv.H"
29 #include "fvmLaplacian.H"
30 
31 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
32 
33 template<class BasicMomentumTransportModel>
35 (
36  const word& modelName,
37  const alphaField& alpha,
38  const rhoField& rho,
39  const volVectorField& U,
40  const surfaceScalarField& alphaRhoPhi,
41  const surfaceScalarField& phi,
42  const viscosity& viscosity
43 )
44 :
45  BasicMomentumTransportModel
46  (
47  modelName,
48  alpha,
49  rho,
50  U,
51  alphaRhoPhi,
52  phi,
53  viscosity
54  )
55 {}
56 
57 
58 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
59 
60 template<class BasicMomentumTransportModel>
62 {
64 }
65 
66 
67 template<class BasicMomentumTransportModel>
70 {
72  (
73  this->groupName("devTau"),
74  (-(this->alpha_*this->rho_*this->nuEff()))
75  *dev(twoSymm(fvc::grad(this->U_)))
76  );
77 }
78 
79 
80 template<class BasicMomentumTransportModel>
83 (
85 ) const
86 {
87  const fvVectorMatrix divDevTauCorr
88  (
89  this->divDevTauCorr
90  (
91  -(this->alpha_*this->rho_*this->nuEff())*dev2(T(fvc::grad(U))),
92  U
93  )
94  );
95 
96  return
97  (
98  divDevTauCorr
99  - fvm::laplacian(this->alpha_*this->rho_*this->nuEff(), U)
100  );
101 }
102 
103 
104 template<class BasicMomentumTransportModel>
107 (
108  const volScalarField& rho,
110 ) const
111 {
112  return
113  (
114  - fvm::laplacian(this->alpha_*rho*this->nuEff(), U)
115  + this->divDevTauCorr
116  (
117  -(this->alpha_*rho*this->nuEff())*dev2(T(fvc::grad(U))),
118  U
119  )
120  );
121 }
122 
123 
124 template<class BasicMomentumTransportModel>
126 {
128 }
129 
130 
131 // ************************************************************************* //
Generic GeometricField class.
static tmp< GeometricField< Type, PatchField, GeoMesh > > New(const word &name, const Internal &, const PtrList< PatchField< Type >> &, const HashPtrTable< Source > &=HashPtrTable< Source >())
Return a temporary field constructed from name,.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:118
BasicMomentumTransportModel::alphaField alphaField
linearViscousStress(const word &modelName, const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const viscosity &viscosity)
Construct from components.
virtual tmp< volSymmTensorField > devTau() const
Return the effective stress tensor.
virtual bool read()=0
Re-read model coefficients if they have changed.
virtual void correct()=0
Solve the turbulence equations and correct the turbulence viscosity.
virtual tmp< fvVectorMatrix > divDevTau(volVectorField &U) const
Return the source term for the momentum equation.
BasicMomentumTransportModel::rhoField rhoField
A class for managing temporary objects.
Definition: tmp.H:55
Abstract base class for all fluid physical properties.
Definition: viscosity.H:50
A class for handling words, derived from string.
Definition: word.H:62
const scalar nuEff
Calculate the divergence of the given field.
Calculate the gradient of the given field.
Calculate the matrix for the laplacian of the field.
U
Definition: pEqn.H:72
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
void correct(const RdeltaTType &rDeltaT, const RhoType &rho, volScalarField &psi, const surfaceScalarField &phiCorr, const SpType &Sp, const SuType &Su)
void read(Istream &, label &, const dictionary &)
In-place read with dictionary lookup.
tmp< VolField< typename outerProduct< vector, Type >::type > > grad(const SurfaceField< Type > &ssf)
Definition: fvcGrad.C:46
tmp< fvMatrix< Type > > laplacian(const VolField< Type > &vf, const word &name)
Definition: fvmLaplacian.C:47
dimensionedSymmTensor dev2(const dimensionedSymmTensor &dt)
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)