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  This file is part of the elastix software.
4 
5  Copyright (c) University Medical Center Utrecht. All rights reserved.
6  See src/CopyrightElastix.txt or http://elastix.isi.uu.nl/legal.php for
7  details.
8 
9  This software is distributed WITHOUT ANY WARRANTY; without even
10  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
11  PURPOSE. See the above copyright notices for more information.
12 
13 ======================================================================*/
14 
15 #ifndef __itkGenericConjugateGradientOptimizer_h
16 #define __itkGenericConjugateGradientOptimizer_h
17 
19 #include "itkLineSearchOptimizer.h"
20 #include <vector>
21 #include <map>
22 
23 namespace itk
24 {
38 {
39 public:
40 
43  typedef SmartPointer< Self > Pointer;
44  typedef SmartPointer< const Self > ConstPointer;
45 
46  itkNewMacro( Self );
49 
56 
59 
63  )( const DerivativeType &,
64  const DerivativeType &,
65  const ParametersType & );
66  typedef std::string BetaDefinitionType;
67  typedef std::map< BetaDefinitionType,
69 
70  typedef enum {
79 
80  virtual void StartOptimization( void );
81 
82  virtual void ResumeOptimization( void );
83 
84  virtual void StopOptimization( void );
85 
87  itkGetConstMacro( CurrentIteration, unsigned long );
88  itkGetConstMacro( CurrentValue, MeasureType );
89  itkGetConstReferenceMacro( CurrentGradient, DerivativeType );
90  itkGetConstMacro( InLineSearch, bool );
91  itkGetConstReferenceMacro( StopCondition, StopConditionType );
92  itkGetConstMacro( CurrentStepLength, double );
93 
95  itkSetObjectMacro( LineSearchOptimizer, LineSearchOptimizerType );
96  itkGetObjectMacro( LineSearchOptimizer, LineSearchOptimizerType );
97 
99  itkGetConstMacro( MaximumNumberOfIterations, unsigned long );
100  itkSetClampMacro( MaximumNumberOfIterations, unsigned long,
102 
109  itkGetConstMacro( GradientMagnitudeTolerance, double );
110  itkSetMacro( GradientMagnitudeTolerance, double )
111 
112 
119  itkGetConstMacro( ValueTolerance, double );
120  itkSetMacro( ValueTolerance, double );
121 
125  virtual void SetMaxNrOfItWithoutImprovement( unsigned long arg );
126 
127  itkGetConstMacro( MaxNrOfItWithoutImprovement, unsigned long );
128 
130  virtual void SetBetaDefinition( const BetaDefinitionType & arg );
131 
132  itkGetConstReferenceMacro( BetaDefinition, BetaDefinitionType );
133 
134 protected:
135 
137  virtual ~GenericConjugateGradientOptimizer(){}
138 
139  void PrintSelf( std::ostream & os, Indent indent ) const;
140 
143  unsigned long m_CurrentIteration;
145  bool m_Stop;
147 
151 
154  itkSetMacro( InLineSearch, bool );
155 
161 
164 
168 
174  virtual void AddBetaDefinition(
175  const BetaDefinitionType & name,
176  ComputeBetaFunctionType function );
177 
188  virtual void ComputeSearchDirection(
189  const DerivativeType & previousGradient,
190  const DerivativeType & gradient,
191  ParametersType & searchDir );
192 
197  virtual void LineSearch(
198  const ParametersType searchDir,
199  double & step,
200  ParametersType & x,
201  MeasureType & f,
202  DerivativeType & g );
203 
208  virtual bool TestConvergence( bool firstLineSearchDone );
209 
211  virtual double ComputeBeta(
212  const DerivativeType & previousGradient,
213  const DerivativeType & gradient,
214  const ParametersType & previousSearchDir );
215 
219  double ComputeBetaSD(
220  const DerivativeType & previousGradient,
221  const DerivativeType & gradient,
222  const ParametersType & previousSearchDir );
223 
225  double ComputeBetaFR(
226  const DerivativeType & previousGradient,
227  const DerivativeType & gradient,
228  const ParametersType & previousSearchDir );
229 
231  double ComputeBetaPR(
232  const DerivativeType & previousGradient,
233  const DerivativeType & gradient,
234  const ParametersType & previousSearchDir );
235 
237  double ComputeBetaDY(
238  const DerivativeType & previousGradient,
239  const DerivativeType & gradient,
240  const ParametersType & previousSearchDir );
241 
243  double ComputeBetaHS(
244  const DerivativeType & previousGradient,
245  const DerivativeType & gradient,
246  const ParametersType & previousSearchDir );
247 
249  double ComputeBetaDYHS(
250  const DerivativeType & previousGradient,
251  const DerivativeType & gradient,
252  const ParametersType & previousSearchDir );
253 
254 private:
255 
256  GenericConjugateGradientOptimizer( const Self & ); // purposely not implemented
257  void operator=( const Self & ); // purposely not implemented
258 
263 
265 
266 };
267 
268 } // end namespace itk
269 
270 #endif //#ifndef __itkGenericConjugateGradientOptimizer_h
double ComputeBetaSD(const DerivativeType &previousGradient, const DerivativeType &gradient, const ParametersType &previousSearchDir)
SmartPointer< Self > Pointer
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)
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)
double(Self::* ComputeBetaFunctionType)(const DerivativeType &, const DerivativeType &, const ParametersType &)
std::map< BetaDefinitionType, ComputeBetaFunctionType > BetaDefinitionMapType


Generated on 11-03-2014 for elastix by doxygen 1.8.6 elastix logo