DASimpleFoam.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 DASimpleFoam
8 
9  This class is modified from OpenFOAM's source code
10  applications/solvers/incompressible/simpleFoam
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 DASimpleFoam_H
34 #define DASimpleFoam_H
35 
36 #include "DASolver.H"
37 #include "addToRunTimeSelectionTable.H"
38 #include "singlePhaseTransportModel.H"
39 #include "turbulentTransportModel.H"
40 #include "simpleControl.H"
41 #include "DAFvSource.H"
42 
43 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
44 
45 namespace Foam
46 {
47 
48 /*---------------------------------------------------------------------------*\
49  Class DASimpleFoam Declaration
50 \*---------------------------------------------------------------------------*/
51 
53  : public DASolver
54 {
55 
56 protected:
58  autoPtr<simpleControl> simplePtr_;
59 
61  autoPtr<volScalarField> pPtr_;
62 
64  autoPtr<volVectorField> UPtr_;
65 
67  autoPtr<surfaceScalarField> phiPtr_;
68 
70  autoPtr<volScalarField> alphaPorosityPtr_;
71 
73  autoPtr<singlePhaseTransportModel> laminarTransportPtr_;
74 
76  autoPtr<incompressible::turbulenceModel> turbulencePtr_;
77 
79  autoPtr<DATurbulenceModel> daTurbulenceModelPtr_;
80 
82  autoPtr<DAFvSource> daFvSourcePtr_;
83 
85  autoPtr<volVectorField> fvSourcePtr_;
86 
88  autoPtr<IOMRFZoneListDF> MRFPtr_;
89 
91  label hasFvSource_ = 0;
92 
94  scalar cumulativeContErr_ = 0.0;
95 
97  label pRefCell_ = 0;
98 
100  scalar pRefValue_ = 0.0;
101 
103  autoPtr<dimensionedScalar> PrPtr_;
104 
106  autoPtr<dimensionedScalar> PrtPtr_;
107 
109  autoPtr<volScalarField> TPtr_;
110 
112  autoPtr<volScalarField> alphatPtr_;
113 
115  label hasTField_;
116 
117 #ifdef CODI_ADR
118  using Position = typename codi::RealReverse::Tape::Position;
120  Position adjResStart_;
121  Position adjResEnd_;
122  Position gradPStart_;
123  Position gradPEnd_;
124  Position pEqnFluxStart_;
125  Position pEqnFluxEnd_;
126  Position divPhiStart_;
127  Position divPhiEnd_;
128  Position fluxUStart_;
129  Position fluxUEnd_;
130 #endif
131 
133  scalar relaxUEqn_ = 1.0;
134  dictionary solverDictU_;
135  dictionary solverDictP_;
136 
137 public:
138  TypeName("DASimpleFoam");
139  // Constructors
140 
141  //- Construct from components
142  DASimpleFoam(
143  char* argsAll,
144  PyObject* pyOptions);
145 
146  //- Destructor
147  virtual ~DASimpleFoam()
148  {
149  }
150 
152  virtual void initSolver();
153 
155  virtual label solvePrimal();
156 
158  virtual label runFPAdj(
159  Vec dFdW,
160  Vec psi);
161 
163  void invTranProdUEqn(
164  const volVectorField& mySource,
165  volVectorField& pseudoU);
166 
168  void invTranProdPEqn(
169  const volScalarField& mySource,
170  volScalarField& pseudoP);
171 
173  void calcLduResiduals(
174  volVectorField& URes,
175  volScalarField& pRes,
176  surfaceScalarField& phiRes);
177 
179  void vec2Fields(
180  const word mode,
181  Vec cVec,
182  volVectorField& UField,
183  volScalarField& pField,
184  surfaceScalarField& phiField,
185  volScalarField& nuTildaField);
186 
187  void calcAdjointResidual(
188  volVectorField& URes,
189  volScalarField& pRes,
190  surfaceScalarField& phiRes,
191  volScalarField& nuTildaRes,
192  Vec dFdW,
193  volVectorField& UPsi,
194  volScalarField& pPsi,
195  surfaceScalarField& phiPsi,
196  volScalarField& nuTildaPsi,
197  volVectorField& adjURes,
198  volScalarField& adjPRes,
199  surfaceScalarField& adjPhiRes,
200  volScalarField& adjNuTildaRes,
201  label& cnt);
202 
203  scalar L2norm(const volScalarField& v);
204 
205  vector L2norm(const volVectorField& U);
206 
207  scalar L2norm(const surfaceScalarField& Phi);
208 
209  void swap(List<scalar>& a, List<scalar>& b);
210 };
211 
212 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
213 
214 } // End namespace Foam
215 
216 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
217 
218 #endif
219 
220 // ************************************************************************* //
Foam::DASimpleFoam::phiPtr_
autoPtr< surfaceScalarField > phiPtr_
surface flux field pointer
Definition: DASimpleFoam.H:67
Foam::DASimpleFoam::DASimpleFoam
DASimpleFoam(char *argsAll, PyObject *pyOptions)
Definition: DASimpleFoam.C:41
Foam::DASimpleFoam::relaxUEqn_
scalar relaxUEqn_
fvSolution parameters for fixed-point adjoint
Definition: DASimpleFoam.H:133
Foam::DASimpleFoam::swap
void swap(List< scalar > &a, List< scalar > &b)
Definition: DASimpleFoam.C:1481
U
volVectorField & U
Definition: createRefsPimpleDyM.H:7
Foam::DASimpleFoam::pRefCell_
label pRefCell_
pressure referefence cell id
Definition: DASimpleFoam.H:97
pseudoP
volScalarField pseudoP("pseudoP", p)
Foam::DASimpleFoam::invTranProdUEqn
void invTranProdUEqn(const volVectorField &mySource, volVectorField &pseudoU)
Inverse transpose product, MU^(-T)
Definition: DASimpleFoam.C:1017
Foam::DASimpleFoam::solvePrimal
virtual label solvePrimal()
solve the primal equations
Definition: DASimpleFoam.C:123
Foam::DASimpleFoam::fvSourcePtr_
autoPtr< volVectorField > fvSourcePtr_
fvSource term
Definition: DASimpleFoam.H:85
Foam::DASimpleFoam::solverDictU_
dictionary solverDictU_
Definition: DASimpleFoam.H:134
Foam::DASimpleFoam::daTurbulenceModelPtr_
autoPtr< DATurbulenceModel > daTurbulenceModelPtr_
DATurbulenceModel pointer.
Definition: DASimpleFoam.H:79
Foam::DASolver
Definition: DASolver.H:54
pseudoU
volVectorField pseudoU("pseudoU", U)
Foam::DASimpleFoam::hasTField_
label hasTField_
whether to have the temperature field
Definition: DASimpleFoam.H:115
Foam::DASimpleFoam::PrPtr_
autoPtr< dimensionedScalar > PrPtr_
Prandtl number pointer.
Definition: DASimpleFoam.H:103
Foam::DASimpleFoam::alphaPorosityPtr_
autoPtr< volScalarField > alphaPorosityPtr_
alpha porosity
Definition: DASimpleFoam.H:70
DAFvSource.H
Foam::DASimpleFoam::laminarTransportPtr_
autoPtr< singlePhaseTransportModel > laminarTransportPtr_
laminar transport properties pointer
Definition: DASimpleFoam.H:73
Foam::DASimpleFoam::hasFvSource_
label hasFvSource_
whether to have fvSource term
Definition: DASimpleFoam.H:91
Foam::DASimpleFoam::initSolver
virtual void initSolver()
initialize fields and variables
Definition: DASimpleFoam.C:81
Foam::DASimpleFoam::TypeName
TypeName("DASimpleFoam")
Foam::DASimpleFoam::solverDictP_
dictionary solverDictP_
Definition: DASimpleFoam.H:135
Foam::DASimpleFoam::TPtr_
autoPtr< volScalarField > TPtr_
temperature field pointer
Definition: DASimpleFoam.H:109
Foam::DASimpleFoam::daFvSourcePtr_
autoPtr< DAFvSource > daFvSourcePtr_
DASource pointer.
Definition: DASimpleFoam.H:82
Foam::DASimpleFoam::cumulativeContErr_
scalar cumulativeContErr_
continuity error
Definition: DASimpleFoam.H:94
Foam::DASimpleFoam::PrtPtr_
autoPtr< dimensionedScalar > PrtPtr_
Turbulence Prandtl pointer.
Definition: DASimpleFoam.H:106
Foam::DASimpleFoam::runFPAdj
virtual label runFPAdj(Vec dFdW, Vec psi)
solve the adjoint equation using the fixed-point iteration method
Definition: DASimpleFoam.C:189
Foam::DASimpleFoam::~DASimpleFoam
virtual ~DASimpleFoam()
Definition: DASimpleFoam.H:147
Foam::DASimpleFoam::turbulencePtr_
autoPtr< incompressible::turbulenceModel > turbulencePtr_
turbulence pointer
Definition: DASimpleFoam.H:76
Foam::DASimpleFoam::L2norm
scalar L2norm(const volScalarField &v)
Definition: DASimpleFoam.C:1428
Foam
Definition: checkGeometry.C:32
Foam::DASimpleFoam::calcLduResiduals
void calcLduResiduals(volVectorField &URes, volScalarField &pRes, surfaceScalarField &phiRes)
calculate all the flow residuals except for the turbulence one
Definition: DASimpleFoam.C:1140
Foam::DASimpleFoam::pPtr_
autoPtr< volScalarField > pPtr_
pressure field pointer
Definition: DASimpleFoam.H:61
Foam::DASimpleFoam::invTranProdPEqn
void invTranProdPEqn(const volScalarField &mySource, volScalarField &pseudoP)
Inverse transpose product, Mp^(-T)
Definition: DASimpleFoam.C:1073
Foam::DASimpleFoam
Definition: DASimpleFoam.H:52
Foam::DASimpleFoam::simplePtr_
autoPtr< simpleControl > simplePtr_
simple pointer
Definition: DASimpleFoam.H:58
Foam::DASimpleFoam::vec2Fields
void vec2Fields(const word mode, Vec cVec, volVectorField &UField, volScalarField &pField, surfaceScalarField &phiField, volScalarField &nuTildaField)
assign coupled petsc vec to fields or vice versa
Definition: DASimpleFoam.C:910
Foam::DASimpleFoam::calcAdjointResidual
void calcAdjointResidual(volVectorField &URes, volScalarField &pRes, surfaceScalarField &phiRes, volScalarField &nuTildaRes, Vec dFdW, volVectorField &UPsi, volScalarField &pPsi, surfaceScalarField &phiPsi, volScalarField &nuTildaPsi, volVectorField &adjURes, volScalarField &adjPRes, surfaceScalarField &adjPhiRes, volScalarField &adjNuTildaRes, label &cnt)
Definition: DASimpleFoam.C:1239
DASolver.H
Foam::DASimpleFoam::UPtr_
autoPtr< volVectorField > UPtr_
velocity field pointer
Definition: DASimpleFoam.H:64
Foam::DASimpleFoam::pRefValue_
scalar pRefValue_
pressure reference value
Definition: DASimpleFoam.H:100
Foam::DASimpleFoam::alphatPtr_
autoPtr< volScalarField > alphatPtr_
thermal diffusivity field pointer
Definition: DASimpleFoam.H:112
Foam::DASimpleFoam::MRFPtr_
autoPtr< IOMRFZoneListDF > MRFPtr_
MRF pointer.
Definition: DASimpleFoam.H:88
psi
const volScalarField & psi
Definition: createRefsRhoPimple.H:14
dafoam_plot3dtransform.mode
mode
Definition: dafoam_plot3dtransform.py:21