DASpalartAllmaras.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2 
3  DAFoam : Discrete Adjoint with OpenFOAM
4  Version : v4
5 
6  Description:
7  Child class for the SpalartAllmaras model
8 
9  This file is modified from OpenFOAM's source code
10  src/TurbulenceModels/turbulenceModels/RAS/SpalartAllmaras/SpalartAllmaras.H
11 
12  OpenFOAM: The Open Source CFD Toolbox
13 
14  Copyright (C): 2011-2016 OpenFOAM Foundation
15 
16  OpenFOAM License:
17 
18  OpenFOAM is free software: you can redistribute it and/or modify it
19  under the terms of the GNU General Public License as published by
20  the Free Software Foundation, either version 3 of the License, or
21  (at your option) any later version.
22 
23  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
24  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
25  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
26  for more details.
27 
28  You should have received a copy of the GNU General Public License
29  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
30 
31 \*---------------------------------------------------------------------------*/
32 
33 #ifndef DASpalartAllmaras_H
34 #define DASpalartAllmaras_H
35 
36 #include "DATurbulenceModel.H"
37 #include "addToRunTimeSelectionTable.H"
38 
39 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
40 
41 namespace Foam
42 {
43 
44 /*---------------------------------------------------------------------------*\
45  Class DASpalartAllmaras Declaration
46 \*---------------------------------------------------------------------------*/
47 
49  : public DATurbulenceModel
50 {
51 
52 protected:
54 
55  dimensionedScalar sigmaNut_;
56  dimensionedScalar kappa_;
57  dimensionedScalar Cb1_;
58  dimensionedScalar Cb2_;
59  dimensionedScalar Cw1_;
60  dimensionedScalar Cw2_;
61  dimensionedScalar Cw3_;
62  dimensionedScalar Cv1_;
63  dimensionedScalar Cs_;
65 
67 
68  tmp<volScalarField> chi() const;
69 
70  tmp<volScalarField> fv1(const volScalarField& chi) const;
71 
72  tmp<volScalarField> fv2(
73  const volScalarField& chi,
74  const volScalarField& fv1) const;
75 
76  tmp<volScalarField> Stilda(
77  const volScalarField& chi,
78  const volScalarField& fv1) const;
79 
80  tmp<volScalarField> fw(const volScalarField& Stilda) const;
81 
83 
85 
86  volScalarField& nuTilda_;
87  volScalarField nuTildaRes_;
89 
91  const volScalarField& y_;
92 
94  volScalarField betaFINuTilda_;
95 
97  label solveTurbState_ = 0;
98 
99  autoPtr<fvScalarMatrix> psiNuTildaPC_;
100 
101  volScalarField dPsiNuTilda_;
102 
103 public:
104  TypeName("SpalartAllmaras");
105  // Constructors
106 
107  //- Construct from components
109  const word modelType,
110  const fvMesh& mesh,
111  const DAOption& daOption);
112 
113  //- Destructor
115  {
116  }
117 
118  // Member functions
119 
121  tmp<volScalarField> DnuTildaEff() const;
122 
124  virtual void correctModelStates(wordList& modelStates) const;
125 
127  virtual void correctNut();
128 
130  virtual void correctBoundaryConditions();
131 
133  virtual void updateIntermediateVariables();
134 
136  virtual void correctStateResidualModelCon(List<List<word>>& stateCon) const;
137 
139  virtual void addModelResidualCon(HashTable<List<List<word>>>& allCon) const;
140 
142  virtual void calcResiduals(const dictionary& options);
143 
145  virtual void correct(label printToScreen);
146 
148  virtual void getFvMatrixFields(
149  const word varName,
150  scalarField& diag,
151  scalarField& upper,
152  scalarField& lower);
153 
155  virtual void getTurbProdOverDestruct(volScalarField& PoD) const;
156 
158  virtual void getTurbConvOverProd(volScalarField& CoP) const;
159 
161  virtual void solveAdjointFP(
162  const word varName,
163  const scalarList& rhs,
164  scalarList& dPsi);
165 };
166 
167 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
168 
169 } // End namespace Foam
170 
171 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
172 
173 #endif
174 
175 // ************************************************************************* //
Foam::DASpalartAllmaras::correct
virtual void correct(label printToScreen)
solve the residual equations and update the state
Definition: DASpalartAllmaras.C:386
Foam::DASpalartAllmaras::sigmaNut_
dimensionedScalar sigmaNut_
Definition: DASpalartAllmaras.H:55
Foam::DASpalartAllmaras::Cw1_
dimensionedScalar Cw1_
Definition: DASpalartAllmaras.H:59
Foam::DASpalartAllmaras::DASpalartAllmaras
DASpalartAllmaras(const word modelType, const fvMesh &mesh, const DAOption &daOption)
Definition: DASpalartAllmaras.C:41
Foam::DASpalartAllmaras::solveTurbState_
label solveTurbState_
whether to solve for turb states
Definition: DASpalartAllmaras.H:97
Foam::DASpalartAllmaras::Cv1_
dimensionedScalar Cv1_
Definition: DASpalartAllmaras.H:62
Foam::DAOption
Definition: DAOption.H:29
Foam::DASpalartAllmaras::correctStateResidualModelCon
virtual void correctStateResidualModelCon(List< List< word >> &stateCon) const
update the original variable connectivity for the adjoint state residuals in stateCon
Definition: DASpalartAllmaras.C:257
Foam::DASpalartAllmaras::getTurbProdOverDestruct
virtual void getTurbProdOverDestruct(volScalarField &PoD) const
return the value of the destruction term from the turbulence model
Definition: DASpalartAllmaras.C:601
Foam::DASpalartAllmaras::betaFINuTilda_
volScalarField betaFINuTilda_
beta field for field inversion
Definition: DASpalartAllmaras.H:94
Foam::DASpalartAllmaras::Cb1_
dimensionedScalar Cb1_
Definition: DASpalartAllmaras.H:57
Foam::DASpalartAllmaras::TypeName
TypeName("SpalartAllmaras")
Foam::DASpalartAllmaras::Cw2_
dimensionedScalar Cw2_
Definition: DASpalartAllmaras.H:60
Foam::DASpalartAllmaras::kappa_
dimensionedScalar kappa_
Definition: DASpalartAllmaras.H:56
mesh
fvMesh & mesh
Definition: createRefsHeatTransfer.H:4
Foam::DASpalartAllmaras::Cs_
dimensionedScalar Cs_
Definition: DASpalartAllmaras.H:63
Foam::DASpalartAllmaras::psiNuTildaPC_
autoPtr< fvScalarMatrix > psiNuTildaPC_
Definition: DASpalartAllmaras.H:99
Foam
Definition: checkGeometry.C:32
Foam::DASpalartAllmaras::chi
tmp< volScalarField > chi() const
Definition: DASpalartAllmaras.C:124
Foam::DASpalartAllmaras
Definition: DASpalartAllmaras.H:48
Foam::DASpalartAllmaras::getFvMatrixFields
virtual void getFvMatrixFields(const word varName, scalarField &diag, scalarField &upper, scalarField &lower)
return the diag(), upper(), and lower() scalarFields from the turbulence model's fvMatrix
Definition: DASpalartAllmaras.C:490
Foam::DASpalartAllmaras::calcResiduals
virtual void calcResiduals(const dictionary &options)
compute the turbulence residuals
Definition: DASpalartAllmaras.C:407
Foam::DASpalartAllmaras::Stilda
tmp< volScalarField > Stilda(const volScalarField &chi, const volScalarField &fv1) const
Definition: DASpalartAllmaras.C:143
Foam::DASpalartAllmaras::dPsiNuTilda_
volScalarField dPsiNuTilda_
Definition: DASpalartAllmaras.H:101
Foam::DASpalartAllmaras::nuTildaRes_
volScalarField nuTildaRes_
Definition: DASpalartAllmaras.H:87
Foam::DASpalartAllmaras::addModelResidualCon
virtual void addModelResidualCon(HashTable< List< List< word >>> &allCon) const
add the model residual connectivity to stateCon
Definition: DASpalartAllmaras.C:307
Foam::DASpalartAllmaras::nuTilda_
volScalarField & nuTilda_
Definition: DASpalartAllmaras.H:86
Foam::DATurbulenceModel
Definition: DATurbulenceModel.H:52
Foam::DASpalartAllmaras::Cw3_
dimensionedScalar Cw3_
Definition: DASpalartAllmaras.H:61
Foam::DASpalartAllmaras::correctBoundaryConditions
virtual void correctBoundaryConditions()
update turbulence variable boundary values
Definition: DASpalartAllmaras.C:235
Foam::DASpalartAllmaras::correctNut
virtual void correctNut()
update nut based on other turbulence variables and update the BCs
Definition: DASpalartAllmaras.C:215
Foam::DASpalartAllmaras::fv1
tmp< volScalarField > fv1(const volScalarField &chi) const
Definition: DASpalartAllmaras.C:129
Foam::DASpalartAllmaras::~DASpalartAllmaras
virtual ~DASpalartAllmaras()
Definition: DASpalartAllmaras.H:114
Foam::DASpalartAllmaras::DnuTildaEff
tmp< volScalarField > DnuTildaEff() const
Return the effective diffusivity for nuTilda.
Definition: DASpalartAllmaras.C:174
Foam::DASpalartAllmaras::correctModelStates
virtual void correctModelStates(wordList &modelStates) const
update the turbulence state for DAStateInfo::regStates_
Definition: DASpalartAllmaras.C:181
Foam::DASpalartAllmaras::solveAdjointFP
virtual void solveAdjointFP(const word varName, const scalarList &rhs, scalarList &dPsi)
solve the fvMatrixT field with given rhs and solution
Definition: DASpalartAllmaras.C:531
Foam::DASpalartAllmaras::Cb2_
dimensionedScalar Cb2_
Definition: DASpalartAllmaras.H:58
Foam::DASpalartAllmaras::updateIntermediateVariables
virtual void updateIntermediateVariables()
update any intermediate variables that are dependent on state variables and are used in calcResiduals
Definition: DASpalartAllmaras.C:246
Foam::DASpalartAllmaras::fv2
tmp< volScalarField > fv2(const volScalarField &chi, const volScalarField &fv1) const
Definition: DASpalartAllmaras.C:136
Foam::DASpalartAllmaras::fw
tmp< volScalarField > fw(const volScalarField &Stilda) const
Definition: DASpalartAllmaras.C:156
Foam::DASpalartAllmaras::getTurbConvOverProd
virtual void getTurbConvOverProd(volScalarField &CoP) const
return the value of the convective over production term from the turbulence model
Definition: DASpalartAllmaras.C:624
DATurbulenceModel.H
Foam::DASpalartAllmaras::y_
const volScalarField & y_
3D wall distance
Definition: DASpalartAllmaras.H:91