go home Home | Main Page | Modules | Namespace List | Class Hierarchy | Alphabetical List | Data Structures | File List | Namespace Members | Data Fields | Globals | Related Pages
itkGenericConjugateGradientOptimizer.h
Go to the documentation of this file.
1 /*=========================================================================
2  *
3  * Copyright UMC Utrecht and contributors
4  *
5  * Licensed under the Apache License, Version 2.0 (the "License");
6  * you may not use this file except in compliance with the License.
7  * You may obtain a copy of the License at
8  *
9  * http://www.apache.org/licenses/LICENSE-2.0.txt
10  *
11  * Unless required by applicable law or agreed to in writing, software
12  * distributed under the License is distributed on an "AS IS" BASIS,
13  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14  * See the License for the specific language governing permissions and
15  * limitations under the License.
16  *
17  *=========================================================================*/
18 
19 #ifndef __itkGenericConjugateGradientOptimizer_h
20 #define __itkGenericConjugateGradientOptimizer_h
21 
23 #include "itkLineSearchOptimizer.h"
24 #include <vector>
25 #include <map>
26 
27 namespace itk
28 {
42 {
43 public:
44 
48  typedef SmartPointer< const Self > ConstPointer;
49 
50  itkNewMacro( Self );
53 
60 
63 
66  typedef double (Self::* ComputeBetaFunctionType
67  )( const DerivativeType &,
68  const DerivativeType &,
69  const ParametersType & );
70  typedef std::string BetaDefinitionType;
71  typedef std::map< BetaDefinitionType,
73 
74  typedef enum {
83 
84  virtual void StartOptimization( void );
85 
86  virtual void ResumeOptimization( void );
87 
88  virtual void StopOptimization( void );
89 
91  itkGetConstMacro( CurrentIteration, unsigned long );
92  itkGetConstMacro( CurrentValue, MeasureType );
93  itkGetConstReferenceMacro( CurrentGradient, DerivativeType );
94  itkGetConstMacro( InLineSearch, bool );
95  itkGetConstReferenceMacro( StopCondition, StopConditionType );
96  itkGetConstMacro( CurrentStepLength, double );
97 
99  itkSetObjectMacro( LineSearchOptimizer, LineSearchOptimizerType );
100  itkGetObjectMacro( LineSearchOptimizer, LineSearchOptimizerType );
101 
103  itkGetConstMacro( MaximumNumberOfIterations, unsigned long );
104  itkSetClampMacro( MaximumNumberOfIterations, unsigned long,
106 
113  itkGetConstMacro( GradientMagnitudeTolerance, double );
114  itkSetMacro( GradientMagnitudeTolerance, double )
115 
116 
123  itkGetConstMacro( ValueTolerance, double );
124  itkSetMacro( ValueTolerance, double );
125 
129  virtual void SetMaxNrOfItWithoutImprovement( unsigned long arg );
130 
131  itkGetConstMacro( MaxNrOfItWithoutImprovement, unsigned long );
132 
134  virtual void SetBetaDefinition( const BetaDefinitionType & arg );
135 
136  itkGetConstReferenceMacro( BetaDefinition, BetaDefinitionType );
137 
138 protected:
139 
141  virtual ~GenericConjugateGradientOptimizer(){}
142 
143  void PrintSelf( std::ostream & os, Indent indent ) const;
144 
145  DerivativeType m_CurrentGradient;
146  MeasureType m_CurrentValue;
147  unsigned long m_CurrentIteration;
149  bool m_Stop;
151 
155 
158  itkSetMacro( InLineSearch, bool );
159 
165 
167  BetaDefinitionType m_BetaDefinition;
168 
171  BetaDefinitionMapType m_BetaDefinitionMap;
172 
178  virtual void AddBetaDefinition(
179  const BetaDefinitionType & name,
180  ComputeBetaFunctionType function );
181 
192  virtual void ComputeSearchDirection(
193  const DerivativeType & previousGradient,
194  const DerivativeType & gradient,
195  ParametersType & searchDir );
196 
201  virtual void LineSearch(
202  const ParametersType searchDir,
203  double & step,
204  ParametersType & x,
205  MeasureType & f,
206  DerivativeType & g );
207 
212  virtual bool TestConvergence( bool firstLineSearchDone );
213 
215  virtual double ComputeBeta(
216  const DerivativeType & previousGradient,
217  const DerivativeType & gradient,
218  const ParametersType & previousSearchDir );
219 
223  double ComputeBetaSD(
224  const DerivativeType & previousGradient,
225  const DerivativeType & gradient,
226  const ParametersType & previousSearchDir );
227 
229  double ComputeBetaFR(
230  const DerivativeType & previousGradient,
231  const DerivativeType & gradient,
232  const ParametersType & previousSearchDir );
233 
235  double ComputeBetaPR(
236  const DerivativeType & previousGradient,
237  const DerivativeType & gradient,
238  const ParametersType & previousSearchDir );
239 
241  double ComputeBetaDY(
242  const DerivativeType & previousGradient,
243  const DerivativeType & gradient,
244  const ParametersType & previousSearchDir );
245 
247  double ComputeBetaHS(
248  const DerivativeType & previousGradient,
249  const DerivativeType & gradient,
250  const ParametersType & previousSearchDir );
251 
253  double ComputeBetaDYHS(
254  const DerivativeType & previousGradient,
255  const DerivativeType & gradient,
256  const ParametersType & previousSearchDir );
257 
258 private:
259 
260  GenericConjugateGradientOptimizer( const Self & ); // purposely not implemented
261  void operator=( const Self & ); // purposely not implemented
262 
267 
268  LineSearchOptimizerPointer m_LineSearchOptimizer;
269 
270 };
271 
272 } // end namespace itk
273 
274 #endif //#ifndef __itkGenericConjugateGradientOptimizer_h
double ComputeBetaSD(const DerivativeType &previousGradient, const DerivativeType &gradient, const ParametersType &previousSearchDir)
virtual void SetBetaDefinition(const BetaDefinitionType &arg)
Superclass::ScaledCostFunctionType ScaledCostFunctionType
double ComputeBetaDY(const DerivativeType &previousGradient, const DerivativeType &gradient, const ParametersType &previousSearchDir)
virtual void AddBetaDefinition(const BetaDefinitionType &name, ComputeBetaFunctionType function)
double ComputeBetaDYHS(const DerivativeType &previousGradient, const DerivativeType &gradient, const ParametersType &previousSearchDir)
double ComputeBetaPR(const DerivativeType &previousGradient, const DerivativeType &gradient, const ParametersType &previousSearchDir)
double ComputeBetaHS(const DerivativeType &previousGradient, const DerivativeType &gradient, const ParametersType &previousSearchDir)
int max(int a, int b)
virtual void SetMaxNrOfItWithoutImprovement(unsigned long arg)
virtual void ComputeSearchDirection(const DerivativeType &previousGradient, const DerivativeType &gradient, ParametersType &searchDir)
double(Self::* ComputeBetaFunctionType)(const DerivativeType &, const DerivativeType &, const ParametersType &)
A base class for LineSearch optimizers.
void PrintSelf(std::ostream &os, Indent indent) const
virtual bool TestConvergence(bool firstLineSearchDone)
double ComputeBetaFR(const DerivativeType &previousGradient, const DerivativeType &gradient, const ParametersType &previousSearchDir)
virtual void LineSearch(const ParametersType searchDir, double &step, ParametersType &x, MeasureType &f, DerivativeType &g)
virtual double ComputeBeta(const DerivativeType &previousGradient, const DerivativeType &gradient, const ParametersType &previousSearchDir)
std::map< BetaDefinitionType, ComputeBetaFunctionType > BetaDefinitionMapType


Generated on 04-09-2015 for elastix by doxygen 1.8.9.1 elastix logo