adjointSolver.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-2019 PCOpt/NTUA
9 Copyright (C) 2013-2019 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::adjointSolver
31
32Description
33 Base class for adjoint solvers
34
35\*---------------------------------------------------------------------------*/
36
37#ifndef adjointSolver_H
38#define adjointSolver_H
39
40#include "fvMesh.H"
41#include "Time.H"
42#include "IOdictionary.H"
43#include "solver.H"
44#include "objectiveManager.H"
45#include "sensitivity.H"
46#include "primalSolver.H"
48
49// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
50
51namespace Foam
52{
53
54/*---------------------------------------------------------------------------*\
55 Class adjointSolver Declaration
56\*---------------------------------------------------------------------------*/
58class adjointSolver
59:
60 public solver
61{
62private:
63
64 // Private Member Functions
65
66 //- No copy construct
67 adjointSolver(const adjointSolver&) = delete;
68
69 //- No copy assignment
70 void operator=(const adjointSolver&) = delete;
71
72
73protected:
74
75 // Protected data
76
77 //- Name of primal solver
79
80 //- Object to manage objective functions
82
83 //- Sensitivities field
85
86 //- Are sensitivities computed
88
89 //- Is the adjoint solver used to tackle a constraint
90 bool isConstraint_;
91
92
93public:
94
95 // Static Data Members
96
97 //- Run-time type information
98 TypeName("adjointSolver");
99
100
101 // Declare run-time constructor selection table
104 (
105 autoPtr,
108 (
109 fvMesh& mesh,
110 const word& managerType,
111 const dictionary& dict,
113 ),
114 (mesh, managerType, dict, primalSolverName)
115 );
116
117
118 // Constructors
119
120 //- Construct from mesh, dictionary, and primal solver name
122 (
123 fvMesh& mesh,
124 const word& managerType,
125 const dictionary& dict,
127 );
128
129
130 // Selectors
131
132 //- Return a reference to the selected turbulence model
134 (
135 fvMesh& mesh,
136 const word& managerType,
137 const dictionary& dict,
139 );
140
141
142 //- Destructor
143 virtual ~adjointSolver() = default;
144
145
146 // Member Functions
147
148 // Access
149
150 virtual bool readDict(const dictionary& dict);
151
152 //- Return the primal solver name
153 const word& primalSolverName() const
154 {
155 return primalSolverName_;
156 }
157
158 //- Return a const-reference to the primal solver
159 //- corresponding to this adjoint solver
160 const primalSolver& getPrimalSolver() const;
161
162 //- Return a non const-reference to the primal solver
163 //- corresponding to this adjoint solver
165
166 //- Return a const reference to the objective manager
168
169 //- Return a reference to the objective manager
171
172
173 // Evolution
174
175 //- Functions to be called after loop
176 // Computes adjoint sensitivities
177 virtual void postLoop();
178
179 //- Compute sensitivities of the underlaying objectives
180 virtual void computeObjectiveSensitivities() = 0;
181
182 //- Is the solving referring to a constraint
183 virtual bool isConstraint();
184
185 //- Grab a reference to the computed sensitivities
186 virtual const scalarField& getObjectiveSensitivities() = 0;
187
188 //- Clears the sensitivity field known by the adjoint solver
189 virtual void clearSensitivities();
190
191 //- Return the base sensitivity object
192 virtual sensitivity& getSensitivityBase() = 0;
193
194 //- Update primal based quantities, e.g. the primal fields
195 //- in adjoint turbulence models
196 // Does nothing in the base
197 virtual void updatePrimalBasedQuantities();
198
199 //- Write the sensitivity derivatives
200 virtual bool writeData(Ostream& os) const;
201};
202
203
204// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
205
206} // End namespace Foam
207
208// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
209
210#endif
211
212// ************************************************************************* //
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
Base class for adjoint solvers.
Definition: adjointSolver.H:60
virtual void clearSensitivities()
Clears the sensitivity field known by the adjoint solver.
virtual const scalarField & getObjectiveSensitivities()=0
Grab a reference to the computed sensitivities.
tmp< scalarField > sensitivities_
Sensitivities field.
Definition: adjointSolver.H:83
autoPtr< objectiveManager > objectiveManagerPtr_
Object to manage objective functions.
Definition: adjointSolver.H:80
static autoPtr< adjointSolver > New(fvMesh &mesh, const word &managerType, const dictionary &dict, const word &primalSolverName)
Return a reference to the selected turbulence model.
Definition: adjointSolver.C:79
virtual bool writeData(Ostream &os) const
Write the sensitivity derivatives.
virtual ~adjointSolver()=default
Destructor.
virtual bool readDict(const dictionary &dict)
const primalSolver & getPrimalSolver() const
bool computeSensitivities_
Are sensitivities computed.
Definition: adjointSolver.H:86
virtual void computeObjectiveSensitivities()=0
Compute sensitivities of the underlaying objectives.
virtual void updatePrimalBasedQuantities()
declareRunTimeNewSelectionTable(autoPtr, adjointSolver, adjointSolver,(fvMesh &mesh, const word &managerType, const dictionary &dict, const word &primalSolverName),(mesh, managerType, dict, primalSolverName))
TypeName("adjointSolver")
Run-time type information.
virtual bool isConstraint()
Is the solving referring to a constraint.
const word primalSolverName_
Name of primal solver.
Definition: adjointSolver.H:77
virtual void postLoop()
Functions to be called after loop.
bool isConstraint_
Is the adjoint solver used to tackle a constraint.
Definition: adjointSolver.H:89
virtual sensitivity & getSensitivityBase()=0
Return the base sensitivity object.
const objectiveManager & getObjectiveManager() const
Return a const reference to the objective manager.
const word & primalSolverName() const
Return the primal solver name.
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 for managing incompressible objective functions.
Base class for primal solvers.
Definition: primalSolver.H:55
Abstract base class for adjoint sensitivities.
Definition: sensitivity.H:64
Base class for solution control classes.
Definition: solver.H:54
virtual const dictionary & dict() const
Return the solver dictionary.
Definition: solver.C:113
const fvMesh & mesh() const
Return the solver mesh.
Definition: solver.C:96
A class for managing temporary objects.
Definition: tmp.H:65
A class for handling words, derived from Foam::string.
Definition: word.H:68
OBJstream os(runTime.globalPath()/outputName)
Namespace for OpenFOAM.
Macros to ease declaration of run-time selection tables.
#define declareRunTimeNewSelectionTable(ptrWrapper, baseType, argNames, argList, parList)
Declare a run-time selection for derived classes.
#define TypeName(TypeNameString)
Declare a ClassName() with extra virtual type info.
Definition: typeInfo.H:73