objective.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-2020 OpenCFD Ltd.
11 -------------------------------------------------------------------------------
12 License
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 Class
29  Foam::objective
30 
31 Description
32  Abstract base class for objective functions. No point in making this
33  runTime selectable since its children will have different constructors.
34 
35 SourceFiles
36  objective.C
37 
38 \*---------------------------------------------------------------------------*/
39 
40 #ifndef objective_H
41 #define objective_H
42 
43 #include "localIOdictionary.H"
44 #include "autoPtr.H"
45 #include "runTimeSelectionTables.H"
46 #include "OFstream.H"
47 #include "boundaryFieldsFwd.H"
48 #include "solverControl.H"
49 #include "objectiveFwd.H"
50 
51 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
52 
53 namespace Foam
54 {
55 
56 /*---------------------------------------------------------------------------*\
57  Class objective Declaration
58 \*---------------------------------------------------------------------------*/
59 
60 class objective
61 :
62  public localIOdictionary
63 {
64 protected:
65 
66  // Protected data
67 
68  const fvMesh& mesh_;
74  bool nullified_;
75  bool normalize_;
76 
77  //- Objective function value and weight
78  scalar J_;
79 
80  //- Average objective value
81  scalar JMean_;
82 
83  //- Objective weight
84  scalar weight_;
85 
86  //- Normalization factor
88 
89  //- Target value, in case the objective is used as a constraint
90  // Should be used in caution and with updateMethods than get affected
91  // by the target value, without requiring a sqr (e.g. SQP, MMA)
92  scalar target_;
93 
94  //- Objective integration start and end times (for unsteady flows)
97 
98  //- Contribution to field sensitivity derivatives
99  // Topology optimisation or other variants with
100  // as many design variables as the mesh cells
101  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
102 
103  autoPtr<volScalarField> dJdbPtr_;
104 
105  // Contribution to surface sensitivity derivatives
106  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
107 
108  //- Term from material derivative
110 
111  //- Term multiplying delta(n dS)/delta b
112  autoPtr<boundaryVectorField> bdSdbMultPtr_;
113 
114  //- Term multiplying delta(n)/delta b
115  autoPtr<boundaryVectorField> bdndbMultPtr_;
116 
117  //- Term multiplying delta(x)/delta b at the boundary
118  autoPtr<boundaryVectorField> bdxdbMultPtr_;
119 
120  //- Term multiplying delta(x)/delta b at the boundary
121  //- for objectives that directly depend on x, e.g. moment
122  //- Needed in both FI and SI computations
123  autoPtr<boundaryVectorField> bdxdbDirectMultPtr_;
124 
125  //- Contribution located in specific parts of a patch.
126  //- Mainly intended for patch boundary edges contributions, e.g.
127  //- NURBS surfaces G1 continuity
128  autoPtr<vectorField3> bEdgeContribution_;
129 
130  //- For use with discrete-like sensitivities
131  autoPtr<boundaryTensorField> bdJdStressPtr_;
132 
133  // Contribution to volume-based sensitivities from volume-based
134  // objective functions
135  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
136 
137  //- Multiplier of d(Volume)/db
138  autoPtr<volScalarField> divDxDbMultPtr_;
139 
140  //- Emerging from volume objectives that include spatial derivatives
141  autoPtr<volTensorField> gradDxDbMultPtr_;
142 
143  //- Output file variables
144  fileName objFunctionFolder_;
145 
146  //- File to keep the objective values after the end of the primal solver
147  mutable autoPtr<OFstream> objFunctionFilePtr_;
148 
149  //- File to keep the objective values at each iteration of the primal
150  //- solver
151  mutable autoPtr<OFstream> instantValueFilePtr_;
152 
153  //- File to keep the average objective values after the end of the
154  //- primal solver
155  mutable autoPtr<OFstream> meanValueFilePtr_;
156 
157 
158  // Protected Member Functions
159 
160  //- Return objective dictionary
161  const dictionary& dict() const;
162 
163  //- Set the output file ptr
164  void setObjectiveFilePtr() const;
165 
166  //- Set the output file ptr for the instantaneous value
167  void setInstantValueFilePtr() const;
168 
169  //- Set the output file ptr for the mean value
170  void setMeanValueFilePtr() const;
171 
172 
173 private:
174 
175  // Private Member Functions
176 
177  //- No copy construct
178  objective(const objective&) = delete;
179 
180  //- No copy assignment
181  void operator=(const objective&) = delete;
182 
183  //- Make objective Function Folder
184  void makeFolder();
185 
186 
187 public:
188 
189  //- Runtime type information
190  TypeName("objective");
191 
192 
193  // Declare run-time constructor selection table
194 
196  (
197  autoPtr,
198  objective,
199  objective,
200  (
201  const fvMesh& mesh,
202  const dictionary& dict,
203  const word& adjointSolverName,
204  const word& primalSolverName
205  ),
206  (mesh, dict, adjointSolverName, primalSolverName)
207  );
208 
209 
210  // Constructors
211 
212  //- Construct from components
213  objective
214  (
215  const fvMesh& mesh,
216  const dictionary& dict,
217  const word& adjointSolverName,
218  const word& primalSolverName
219  );
220 
221 
222  // Selectors
223 
224  //- Return a reference to the selected turbulence model
225  static autoPtr<objective> New
226  (
227  const fvMesh& mesh,
228  const dictionary& dict,
229  const word& objectiveType,
230  const word& adjointSolverName,
231  const word& primalSolverName
232  );
233 
234 
235  //- Destructor
236  virtual ~objective() = default;
237 
238 
239  // Member Functions
240 
241  virtual bool readDict(const dictionary& dict);
242 
243  //- Return the instantaneous objective function value
244  virtual scalar J() = 0;
245 
246  //- Return the mean objective function value, if it exists,
247  //- otherwise the mean one
248  scalar JCycle() const;
249 
250  //- Accumulate contribution for the mean objective value
251  // For steady-state runs
252  void accumulateJMean(solverControl& solverControl);
253 
254  //- Accumulate contribution for the mean objective value
255  // For unsteady runs
256  void accumulateJMean();
257 
258  //- Return the objective function weight
259  scalar weight() const;
260 
261  //- Is the objective normalized
262  bool normalize() const;
263 
264  //- Normalize all fields allocated by the objective
265  virtual void doNormalization();
266 
267  //- Check whether this is an objective integration time
268  bool isWithinIntegrationTime() const;
269 
270  //- Increment integration times
271  void incrementIntegrationTimes(const scalar timeSpan);
272 
273  //- Contribution to field sensitivities
274  const volScalarField& dJdb();
275 
276  //- Contribution to surface sensitivities for a specific patch
277  const fvPatchVectorField& boundarydJdb(const label);
278 
279  //- Multiplier of delta(n dS)/delta b
280  const fvPatchVectorField& dSdbMultiplier(const label);
281 
282  //- Multiplier of delta(n dS)/delta b
283  const fvPatchVectorField& dndbMultiplier(const label);
284 
285  //- Multiplier of delta(x)/delta b
286  const fvPatchVectorField& dxdbMultiplier(const label);
287 
288  //- Multiplier of delta(x)/delta b
289  const fvPatchVectorField& dxdbDirectMultiplier(const label);
290 
291  //- Multiplier located at patch boundary edges
292  const vectorField& boundaryEdgeMultiplier
293  (
294  const label patchI,
295  const label edgeI
296  );
297 
298  //- Objective partial deriv wrt stress tensor
299  const fvPatchTensorField& boundarydJdStress(const label);
300 
301  //- Contribution to surface sensitivities for all patches
302  const boundaryVectorField& boundarydJdb();
303 
304  //- Multiplier of delta(n dS)/delta b for all patches
305  const boundaryVectorField& dSdbMultiplier();
306 
307  //- Multiplier of delta(n dS)/delta b for all patches
308  const boundaryVectorField& dndbMultiplier();
309 
310  //- Multiplier of delta(x)/delta b for all patches
311  const boundaryVectorField& dxdbMultiplier();
312 
313  //- Multiplier of delta(x)/delta b for all patches
314  const boundaryVectorField& dxdbDirectMultiplier();
315 
316  //- Multiplier located at patch boundary edges
317  const vectorField3& boundaryEdgeMultiplier();
318 
319  //- Objective partial deriv wrt stress tensor
320  const boundaryTensorField& boundarydJdStress();
321 
322  //- Multiplier of grad( delta(x)/delta b) for volume-based sensitivities
323  const volScalarField& divDxDbMultiplier();
324 
325  //- Multiplier of grad( delta(x)/delta b) for volume-based sensitivities
326  const volTensorField& gradDxDbMultiplier();
327 
328  //- Update objective function derivatives
329  virtual void update() = 0;
330 
331  //- Nullify adjoint contributions
332  virtual void nullify();
333 
334  //- Update normalization factors, for objectives in
335  //- which the factor is not known a priori
336  virtual void updateNormalizationFactor();
337 
338  //- Update objective function derivative term
339  virtual void update_boundarydJdb()
340  {}
341 
342  //- Update d (normal dS) / db multiplier. Surface-based sensitivity term
343  virtual void update_dSdbMultiplier()
344  {}
345 
346  //- Update d (normal) / db multiplier. Surface-based sensitivity term
347  virtual void update_dndbMultiplier()
348  {}
349 
350  //- Update d (x) / db multiplier. Surface-based sensitivity term
351  virtual void update_dxdbMultiplier()
352  {}
353 
354  //- Update d (x) / db multiplier. Surface and volume-based sensitivity
355  //- term
356  virtual void update_dxdbDirectMultiplier()
357  {}
358 
359  //- Update boundary edge contributions
360  virtual void update_boundaryEdgeContribution()
361  {}
362 
363  //- Update dJ/dStress field
364  virtual void update_dJdStressMultiplier()
365  {}
366 
367  //- Update div( dx/db multiplier). Volume-based sensitivity term
368  virtual void update_divDxDbMultiplier()
369  {}
370 
371  //- Update grad( dx/db multiplier). Volume-based sensitivity term
372  virtual void update_gradDxDbMultiplier()
373  {}
374 
375  //- Write objective function history
376  virtual bool write(const bool valid = true) const;
377 
378  //- Write objective function history at each primal solver iteration
379  virtual void writeInstantaneousValue() const;
380 
381  //- Write objective function history at each primal solver iteration
382  virtual void writeInstantaneousSeparator() const;
383 
384  //- Write mean objective function history
385  virtual void writeMeanValue() const;
386 
387  //- Write averaged objective for continuation
388  virtual bool writeData(Ostream& os) const;
389 
390  // Inline functions for checking whether pointers are set or not
391  inline const word& objectiveName() const;
392  inline bool hasdJdb() const;
393  inline bool hasBoundarydJdb() const;
394  inline bool hasdSdbMult() const;
395  inline bool hasdndbMult() const;
396  inline bool hasdxdbMult() const;
397  inline bool hasdxdbDirectMult() const;
398  inline bool hasBoundaryEdgeContribution() const;
399  inline bool hasBoundarydJdStress() const;
400  inline bool hasDivDxDbMult() const;
401  inline bool hasGradDxDbMult() const;
402 
403  // Inline functions for checking whether integration times are set
404  inline bool hasIntegrationStartTime() const;
405  inline bool hasIntegrationEndTime() const;
406 };
407 
408 
409 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
410 
411 } // End namespace Foam
412 
413 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
414 
415 #include "objectiveI.H"
416 
417 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
418 
419 #endif
420 
421 // ************************************************************************* //
Foam::fvPatchField< vector >
Foam::objective::normFactor_
autoPtr< scalar > normFactor_
Normalization factor.
Definition: objective.H:86
Foam::objective::integrationStartTimePtr_
autoPtr< scalar > integrationStartTimePtr_
Objective integration start and end times (for unsteady flows)
Definition: objective.H:94
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:62
objectiveI.H
Foam::objective::normalize_
bool normalize_
Definition: objective.H:74
Foam::fileName
A class for handling file names.
Definition: fileName.H:69
Foam::baseIOdictionary::writeData
virtual bool writeData(Ostream &) const
The writeData function required by regIOobject write operation.
Definition: baseIOdictionaryIO.C:51
Foam::objective::objectiveName_
const word objectiveName_
Definition: objective.H:71
Foam::dictionary::New
static autoPtr< dictionary > New(Istream &is)
Construct top-level dictionary on freestore from Istream.
Definition: dictionaryIO.C:69
update
mesh update()
Foam::objective::target_
scalar target_
Target value, in case the objective is used as a constraint.
Definition: objective.H:91
Foam::solverControl
Base class for solver control classes.
Definition: solverControl.H:51
objectiveFwd.H
Foam::objective::adjointSolverName_
const word adjointSolverName_
Definition: objective.H:69
OFstream.H
Foam::Field< vector >
Foam::objective::J_
scalar J_
Objective function value and weight.
Definition: objective.H:77
Foam::boundaryTensorField
volTensorField::Boundary boundaryTensorField
Definition: boundaryFieldsFwd.H:56
Foam::baseIOdictionary::TypeName
TypeName("dictionary")
Declare type-name, virtual type (with debug switch)
Foam::objective::primalSolverName_
const word primalSolverName_
Definition: objective.H:70
Foam::objective::weight_
scalar weight_
Objective weight.
Definition: objective.H:83
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:121
Foam::objective::integrationEndTimePtr_
autoPtr< scalar > integrationEndTimePtr_
Definition: objective.H:95
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
Foam::fvMesh
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:84
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::autoPtr< scalar >
Foam::objective::nullified_
bool nullified_
Definition: objective.H:73
Foam::objective::computeMeanFields_
bool computeMeanFields_
Definition: objective.H:72
Foam::GeometricField::Boundary
Definition: GeometricField.H:114
boundaryFieldsFwd.H
Useful typenames for fields defined only at the boundaries.
runTimeSelectionTables.H
Macros to ease declaration of run-time selection tables.
Foam::dictionary::write
void write(Ostream &os, const bool subDict=true) const
Write dictionary, normally with sub-dictionary formatting.
Definition: dictionaryIO.C:206
Foam::objective::mesh_
const fvMesh & mesh_
Definition: objective.H:67
localIOdictionary.H
Foam::objective::dict_
dictionary dict_
Definition: objective.H:68
Foam::objective::JMean_
scalar JMean_
Average objective value.
Definition: objective.H:80
Foam::baseIOdictionary::operator=
void operator=(const baseIOdictionary &)
Copy assignment of dictionary entries (leave regIOobject untouched)
Definition: baseIOdictionary.C:90
Foam::GeometricField< scalar, fvPatchField, volMesh >
Foam::objective
Abstract base class for objective functions. No point in making this runTime selectable since its chi...
Definition: objective.H:59
Foam::localIOdictionary
localIOdictionary is derived from IOdictionary but excludes parallel master reading.
Definition: localIOdictionary.H:51
declareRunTimeNewSelectionTable
#define declareRunTimeNewSelectionTable( autoPtr, baseType, argNames, argList, parList)
Declare a run-time selection for derived classes.
Definition: runTimeSelectionTables.H:148
Foam::normalize
quaternion normalize(const quaternion &q)
Return the normalized (unit) quaternion of the given quaternion.
Definition: quaternionI.H:661
solverControl.H
autoPtr.H