lineSearch.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 -------------------------------------------------------------------------------
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 
29 Class
30  Foam::lineSearch
31 
32 Description
33  Abstract base class for line search methods
34 
35 SourceFiles
36  lineSearch.C
37 
38 \*---------------------------------------------------------------------------*/
39 
40 #ifndef lineSearch_H
41 #define lineSearch_H
42 
43 #include "runTimeSelectionTables.H"
44 #include "IOdictionary.H"
45 #include "scalarField.H"
46 #include "stepUpdate.H"
47 
48 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
49 
50 namespace Foam
51 {
52 
53 /*---------------------------------------------------------------------------*\
54  Class lineSearch Declaration
55 \*---------------------------------------------------------------------------*/
56 
57 class lineSearch
58 {
59 protected:
60 
61  // Protected data
62 
63  //- Subdict within updateMethod
64  const dictionary dict_;
65 
66  //- IOdictionary under time/uniform for continuation
68 
69  //- Directional derivative of the merit function
70  scalar directionalDeriv_;
71 
72  //- Update direction
74 
75  //- Old merit value from this opt cycle
76  scalar oldMeritValue_;
77 
78  //- New merit value from this opt cycle
79  scalar newMeritValue_;
80 
81  //- Merit directional deriv from the previous opt cycle
82  scalar prevMeritDeriv_;
83 
84  //- Correction multiplier at the first step of line search
85  scalar initialStep_;
86 
87  //- Minimum allowed correction multiplier
88  scalar minStep_;
89 
90  //- Correction multiplier
91  scalar step_;
92 
93  //- Inner line search iteration
94  label iter_;
95 
96  //- Maximum line search iterations
97  label maxIters_;
98 
99  //- Whether to extrapolate the correction multiplier for
100  //- this optimisation cycle based on the previous ones.
101  // Useful for non-quasi Newton methods
103 
104  //- Mechanism to update method if line search conditions are not set
106 
107 
108  // Protected Member Functions
109 
110  //- Optional coeffs dict
111  const dictionary& coeffsDict();
112 
113 
114 private:
115 
116  // Private Member Functions
117 
118  //- No copy construct
119  lineSearch(const lineSearch&) = delete;
120 
121  //- No copy assignment
122  void operator=(const lineSearch&) = delete;
123 
124 
125 public:
126 
127  //- Runtime type information
128  TypeName("lineSearch");
129 
130 
131  // Declare run-time constructor selection table
132 
134  (
135  autoPtr,
136  lineSearch,
137  dictionary,
138  (
139  const dictionary& dict,
140  const Time& time
141  ),
142  (dict, time)
143  );
144 
145 
146  // Constructors
147 
148  //- Construct from components
149  lineSearch(const dictionary& dict, const Time& time);
150 
151  // Selectors
152 
153  //- Return a reference to the selected turbulence model
154  static autoPtr<lineSearch> New
155  (
156  const dictionary& dict,
157  const Time& time
158  );
159 
160 
161  //- Destructor
162  virtual ~lineSearch() = default;
163 
164 
165  // Member Functions
166 
167  //- Set objective derivative
168  virtual void setDeriv(const scalar deriv);
169 
170  //- Set direction
171  void setDirection(const scalarField& direction);
172 
173  //- Set new objective value
174  void setNewMeritValue(const scalar value);
175 
176  //- Set old objective value
177  void setOldMeritValue(const scalar value);
178 
179  //- Reset step to initial value
180  virtual void reset();
181 
182  //- Get max number of iterations
183  label maxIters() const;
184 
185  //- Get current step
186  scalar step() const;
187 
188  //- Return the correction of the design variables
189  virtual bool converged() = 0;
190 
191  //- Update the line search step based on the specific line search
192  //- strategy, e.g. bisection, quadratic fit, etc.
193  virtual void updateStep() = 0;
194 
195  //- Update the step using the supplied value
196  virtual void updateStep(const scalar newStep);
197 
198 
199  // Member operators
200 
201  //- Increment iteration number and store merit value corresponding to
202  //- the previous optimisation cycle
203  virtual lineSearch& operator++();
204 
205  //- Postfix increment. Necessary for compilation
206  virtual lineSearch& operator++(int);
207 };
208 
209 
210 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
211 
212 } // End namespace Foam
213 
214 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
215 
216 #endif
217 
218 // ************************************************************************* //
Foam::lineSearch::declareRunTimeSelectionTable
declareRunTimeSelectionTable(autoPtr, lineSearch, dictionary,(const dictionary &dict, const Time &time),(dict, time))
Foam::IOdictionary
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:54
Foam::lineSearch::reset
virtual void reset()
Reset step to initial value.
Definition: lineSearch.C:163
Foam::lineSearch::setOldMeritValue
void setOldMeritValue(const scalar value)
Set old objective value.
Definition: lineSearch.C:156
Foam::lineSearch::prevMeritDeriv_
scalar prevMeritDeriv_
Merit directional deriv from the previous opt cycle.
Definition: lineSearch.H:81
Foam::lineSearch::setNewMeritValue
void setNewMeritValue(const scalar value)
Set new objective value.
Definition: lineSearch.C:149
Foam::Time
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:73
Foam::lineSearch::~lineSearch
virtual ~lineSearch()=default
Destructor.
scalarField.H
Foam::lineSearch::step
scalar step() const
Get current step.
Definition: lineSearch.C:189
Foam::lineSearch::initialStep_
scalar initialStep_
Correction multiplier at the first step of line search.
Definition: lineSearch.H:84
Foam::lineSearch::setDeriv
virtual void setDeriv(const scalar deriv)
Set objective derivative.
Definition: lineSearch.C:136
Foam::lineSearch::directionalDeriv_
scalar directionalDeriv_
Directional derivative of the merit function.
Definition: lineSearch.H:69
Foam::lineSearch::operator++
virtual lineSearch & operator++()
Definition: lineSearch.C:201
Foam::lineSearch::iter_
label iter_
Inner line search iteration.
Definition: lineSearch.H:93
Foam::lineSearch::dict_
const dictionary dict_
Subdict within updateMethod.
Definition: lineSearch.H:63
Foam::lineSearch::maxIters
label maxIters() const
Get max number of iterations.
Definition: lineSearch.C:183
Foam::lineSearch::step_
scalar step_
Correction multiplier.
Definition: lineSearch.H:90
Foam::lineSearch::setDirection
void setDirection(const scalarField &direction)
Set direction.
Definition: lineSearch.C:143
Foam::lineSearch::TypeName
TypeName("lineSearch")
Runtime type information.
Foam::Field< scalar >
Foam::lineSearch::New
static autoPtr< lineSearch > New(const dictionary &dict, const Time &time)
Return a reference to the selected turbulence model.
Definition: lineSearch.C:96
stepUpdate.H
Foam::lineSearch::lineSearchDict_
IOdictionary lineSearchDict_
IOdictionary under time/uniform for continuation.
Definition: lineSearch.H:66
Foam::lineSearch::converged
virtual bool converged()=0
Return the correction of the design variables.
Foam::lineSearch::direction_
scalarField direction_
Update direction.
Definition: lineSearch.H:72
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:123
Foam::lineSearch::newMeritValue_
scalar newMeritValue_
New merit value from this opt cycle.
Definition: lineSearch.H:78
Foam::lineSearch::coeffsDict
const dictionary & coeffsDict()
Optional coeffs dict.
Definition: lineSearch.C:44
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
IOdictionary.H
Foam::autoPtr
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: HashPtrTable.H:53
Foam::lineSearch::maxIters_
label maxIters_
Maximum line search iterations.
Definition: lineSearch.H:96
Foam::lineSearch::oldMeritValue_
scalar oldMeritValue_
Old merit value from this opt cycle.
Definition: lineSearch.H:75
runTimeSelectionTables.H
Macros to ease declaration of run-time selection tables.
Foam::lineSearch::minStep_
scalar minStep_
Minimum allowed correction multiplier.
Definition: lineSearch.H:87
Foam::direction
uint8_t direction
Definition: direction.H:52
Foam::lineSearch::updateStep
virtual void updateStep()=0
Foam::lineSearch::extrapolateInitialStep_
bool extrapolateInitialStep_
Definition: lineSearch.H:101
Foam::lineSearch
Abstract base class for line search methods.
Definition: lineSearch.H:56
Foam::lineSearch::stepUpdate_
autoPtr< stepUpdate > stepUpdate_
Mechanism to update method if line search conditions are not set.
Definition: lineSearch.H:104