adjointEikonalSolverIncompressible.H
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2007-2020 PCOpt/NTUA
9 Copyright (C) 2013-2020 FOSS GP
10 Copyright (C) 2019 OpenCFD Ltd.
11-------------------------------------------------------------------------------
12License
13 This file is part of OpenFOAM.
14
15 OpenFOAM is free software: you can redistribute it and/or modify it
16 under the terms of the GNU General Public License as published by
17 the Free Software Foundation, either version 3 of the License, or
18 (at your option) any later version.
19
20 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
21 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
22 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
23 for more details.
24
25 You should have received a copy of the GNU General Public License
26 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
27
28
29Class
30 Foam::incompressible::adjointEikonalSolver
31
32Description
33 Solver of the adjoint to the eikonal PDE
34
35 Reference:
36 \verbatim
37 For the development of the adjoint eikonal PDE and its boundary
38 conditions
39
40 Papoutsis-Kiachagias, E. M., & Giannakoglou, K. C. (2014).
41 Continuous Adjoint Methods for Turbulent Flows, Applied to Shape
42 and Topology Optimization: Industrial Applications.
43 Archives of Computational Methods in Engineering, 23(2), 255-299.
44 http://doi.org/10.1007/s11831-014-9141-9
45 \endverbatim
46
47 To be as consistent as possible, it is recommended to use the
48 advectionDiffusion wallDist method in fvSchemes, instead of the more widely
49 used meshWave
50
51 Example of the adjointEikonalSolver specification in optimisationDict:
52 \verbatim
53 optimisation
54 {
55 sensitivities
56 {
57 includeDistance true;
58 adjointEikonalSolver
59 {
60 // epsilon should be the same as the one used
61 // in fvSchemes/wallDist/advectionDiffusionCoeffs
62 epsilon 0.1;
63 iters 1000;
64 tolerance 1e-6;
65 }
66 }
67 }
68 \endverbatim
69
70 Example of the entries in fvSchemes:
71 \verbatim
72 divSchemes
73 {
74 .
75 .
76 // avoid bounded schemes since yPhi is not conservative
77 div(-yPhi,da) Gauss linearUpwind grad(da);
78 .
79 .
80 }
81 laplacianSchemes
82 {
83 .
84 .
85 laplacian(yWall,da) Gauss linear corrected;
86 .
87 .
88 }
89 \endverbatim
90
91 Also, the solver specification and a relaxation factor for da are required
92 in fvSolution
93
94 \verbatim
95 da
96 {
97 solver PBiCGStab;
98 preconditioner DILU;
99 tolerance 1e-9;
100 relTol 0.1;
101 }
102
103 relaxationFactors
104 {
105 equations
106 {
107 .
108 .
109 da 0.5;
110 .
111 .
112 }
113 }
114 \endverbatim
115
116See also
117 Foam::patchDistMethod::advectionDiffusion
118 Foam::wallDist
119
120
121SourceFiles
122 adjointEikonalSolver.C
123
124\*---------------------------------------------------------------------------*/
125
126#ifndef adjointEikonalSolverIncompressible_H
127#define adjointEikonalSolverIncompressible_H
128
129#include "IOdictionary.H"
131#include "createZeroField.H"
132#include "boundaryFieldsFwd.H"
133#include "RASModelVariables.H"
134
135// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
136
137namespace Foam
138{
139
140namespace incompressible
141{
142
143/*---------------------------------------------------------------------------*\
144 Class adjointEikonalSolver Declaration
145\*---------------------------------------------------------------------------*/
148{
149private:
150
151 // Private Member Functions
152
153 //- No copy construct
155
156 //- No copy assignment
157 void operator=(const adjointEikonalSolver&) = delete;
158
159
160protected:
161
162 // Protected data
164 const fvMesh& mesh_;
169
175 label nEikonalIters_;
177 scalar tolerance_;
179 scalar epsilon_;
186
187 //- Wall face sens w.r.t. (x,y.z)
189
190
191 // Protected functions
192
193 //- Return the boundary condition types for da
194 wordList patchTypes() const;
195
196 //- Compute convecting velocity
198
199 //- Read options each time a new solution is found
200 void read();
201
202
203
204public:
205
206 //- Runtime type information
207 TypeName("adjointEikonalSolver");
208
209
210 // Constructors
211
212 //- Construct from components
214 (
215 const fvMesh& mesh,
216 const dictionary& dict,
218 incompressibleAdjointVars& adjointVars,
219 const labelHashSet& sensitivityPatchIDs
220 );
221
222 // Destructor
224 virtual ~adjointEikonalSolver() = default;
225
226
227 // Member Functions
228
229 //- Read dict if changed
230 virtual bool readDict(const dictionary& dict);
231
232 //- Accumulate source term
233 void accumulateIntegrand(const scalar dt);
234
235 //- Calculate the adjoint distance field
236 void solve();
237
238 //- Reset source term
239 void reset();
240
241 //- Return the sensitivity term depending on da
243
244 //- Return the volume-based sensitivity term depending on da
246
247 //- Return the adjoint distance field
248 const volScalarField& da();
249
250 //- Return the distance field
251 const volScalarField& d();
252
253 //- Return the gradient of the eikonal equation
255};
256
257
258// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
259
260} // End namespace incompressible
261} // End namespace Foam
262
263// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
264
265#endif
266
267// ************************************************************************* //
Useful typenames for fields defined only at the boundaries.
Generic GeometricBoundaryField class.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: autoPtr.H:66
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
Class including all adjoint fields for incompressible flows.
tmp< volVectorField > gradEikonal()
Return the gradient of the eikonal equation.
boundaryVectorField & distanceSensitivities()
Return the sensitivity term depending on da.
const volScalarField & da()
Return the adjoint distance field.
autoPtr< boundaryVectorField > distanceSensPtr_
Wall face sens w.r.t. (x,y.z)
tmp< surfaceScalarField > computeYPhi()
Compute convecting velocity.
autoPtr< Foam::incompressibleAdjoint::adjointRASModel > & adjointTurbulence_
const autoPtr< incompressible::RASModelVariables > & RASModelVars_
const volScalarField & d()
Return the distance field.
virtual bool readDict(const dictionary &dict)
Read dict if changed.
tmp< volTensorField > getFISensitivityTerm() const
Return the volume-based sensitivity term depending on da.
void read()
Read options each time a new solution is found.
void accumulateIntegrand(const scalar dt)
Accumulate source term.
TypeName("adjointEikonalSolver")
Runtime type information.
wordList patchTypes() const
Return the boundary condition types for da.
void solve()
Calculate the adjoint distance field.
A class for managing temporary objects.
Definition: tmp.H:65
dynamicFvMesh & mesh
Namespace for OpenFOAM.
dictionary dict
#define TypeName(TypeNameString)
Declare a ClassName() with extra virtual type info.
Definition: typeInfo.H:73