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-2021 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 "fvc.H"
28 #include "fvm.H"
29 
30 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
31 
32 template<class BasicMomentumTransportModel>
34 (
35  const word& modelName,
36  const alphaField& alpha,
37  const rhoField& rho,
38  const volVectorField& U,
39  const surfaceScalarField& alphaRhoPhi,
40  const surfaceScalarField& phi,
41  const viscosity& viscosity
42 )
43 :
44  BasicMomentumTransportModel
45  (
46  modelName,
47  alpha,
48  rho,
49  U,
50  alphaRhoPhi,
51  phi,
52  viscosity
53  )
54 {}
55 
56 
57 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
58 
59 template<class BasicMomentumTransportModel>
61 {
63 }
64 
65 
66 template<class BasicMomentumTransportModel>
69 {
71  (
72  IOobject::groupName("devTau", this->alphaRhoPhi_.group()),
73  (-(this->alpha_*this->rho_*this->nuEff()))
74  *dev(twoSymm(fvc::grad(this->U_)))
75  );
76 }
77 
78 
79 template<class BasicMomentumTransportModel>
82 (
84 ) const
85 {
86  return
87  (
88  - fvc::div((this->alpha_*this->rho_*this->nuEff())*dev2(T(fvc::grad(U))))
89  - fvm::laplacian(this->alpha_*this->rho_*this->nuEff(), U)
90  );
91 }
92 
93 
94 template<class BasicMomentumTransportModel>
97 (
98  const volScalarField& rho,
100 ) const
101 {
102  return
103  (
104  - fvc::div((this->alpha_*rho*this->nuEff())*dev2(T(fvc::grad(U))))
105  - fvm::laplacian(this->alpha_*rho*this->nuEff(), U)
106  );
107 }
108 
109 
110 template<class BasicMomentumTransportModel>
112 {
114 }
115 
116 
117 // ************************************************************************* //
virtual void correct()=0
Solve the turbulence equations and correct the turbulence viscosity.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(rho0+psi *p, rhoMin);# 1 "/home/ubuntu/OpenFOAM-10/applications/solvers/multiphase/cavitatingFoam/alphavPsi.H" 1{ alphav=max(min((rho - rholSat)/(rhovSat - rholSat), scalar(1)), scalar(0));alphal=1.0 - alphav;Info<< "max-min alphav: "<< max(alphav).value()<< " "<< min(alphav).value()<< endl;psiModel-> correct()
Definition: pEqn.H:68
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.
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
void read(Istream &, label &, const dictionary &)
In-place read with dictionary lookup.
tmp< GeometricField< Type, fvPatchField, volMesh > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcLaplacian.C:45
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
A class for handling words, derived from string.
Definition: word.H:59
BasicMomentumTransportModel::rhoField rhoField
Definition: LESModel.H:107
virtual bool read()=0
Re-read model coefficients if they have changed.
const volScalarField & T
dimensionedSymmTensor dev2(const dimensionedSymmTensor &dt)
virtual tmp< volSymmTensorField > devTau() const
Return the effective stress tensor.
BasicMomentumTransportModel::alphaField alphaField
Definition: LESModel.H:106
A class for managing temporary objects.
Definition: PtrList.H:53
virtual tmp< fvVectorMatrix > divDevTau(volVectorField &U) const
Return the source term for the momentum equation.