1 /*
2  * Medical Image Registration ToolKit (MIRTK)
3  *
4  * Copyright 2013-2015 Imperial College London
5  * Copyright 2013-2015 Andreas Schuh
6  *
7  * Licensed under the Apache License, Version 2.0 (the "License");
8  * you may not use this file except in compliance with the License.
9  * You may obtain a copy of the License at
10  *
11  *     http://www.apache.org/licenses/LICENSE-2.0
12  *
13  * Unless required by applicable law or agreed to in writing, software
14  * distributed under the License is distributed on an "AS IS" BASIS,
15  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
16  * See the License for the specific language governing permissions and
17  * limitations under the License.
18  */
19 
20 #include "mirtk/LineSearch.h"
21 
22 #include "mirtk/String.h"
23 
24 #include "mirtk/MaxStepLineSearch.h"
25 #include "mirtk/AdaptiveLineSearch.h"
26 #include "mirtk/BrentLineSearch.h"
27 
28 #include <limits>
29 
30 
31 namespace mirtk {
32 
33 
34 // =============================================================================
35 // Factory method
36 // =============================================================================
37 
38 // -----------------------------------------------------------------------------
New(LineSearchStrategy & strategy,ObjectiveFunction * f)39 LineSearch *LineSearch::New(LineSearchStrategy &strategy, ObjectiveFunction *f)
40 {
41   switch (strategy) {
42     case LS_None:     return new MaxStepLineSearch (f);
43     case LS_Adaptive: return new AdaptiveLineSearch(f);
44     case LS_Brent:    return new BrentLineSearch(f);
45     default:
46       cerr << "LineSearch::New: Unknown line search strategy: " << strategy << endl;
47       exit(1);
48   }
49   return NULL;
50 }
51 
52 // =============================================================================
53 // Construction/Destruction
54 // =============================================================================
55 
56 // -----------------------------------------------------------------------------
LineSearch(ObjectiveFunction * f)57 LineSearch::LineSearch(ObjectiveFunction *f)
58 :
59   LocalOptimizer(f),
60   _Direction         (NULL),
61   _Revert            (false),
62   _CurrentValue      (numeric_limits<double>::quiet_NaN()),
63   _NumberOfIterations(20),
64   _MinStepLength     ( 0.1),
65   _MaxStepLength     (10.0),
66   _StepLengthUnit    ( 1.0),
67   _StepLength        ( 0.0)
68 {
69 }
70 
71 // -----------------------------------------------------------------------------
CopyAttributes(const LineSearch & other)72 void LineSearch::CopyAttributes(const LineSearch &other)
73 {
74   _Direction          = other._Direction;
75   _Revert             = other._Revert;
76   _CurrentValue       = other._CurrentValue;
77   _NumberOfIterations = other._NumberOfIterations;
78   _MinStepLength      = other._MinStepLength;
79   _MaxStepLength      = other._MaxStepLength;
80   _StepLengthUnit     = other._StepLengthUnit;
81   _StepLength         = other._StepLength;
82 }
83 
84 // -----------------------------------------------------------------------------
LineSearch(const LineSearch & other)85 LineSearch::LineSearch(const LineSearch &other)
86 :
87   LocalOptimizer(other)
88 {
89   CopyAttributes(other);
90 }
91 
92 // -----------------------------------------------------------------------------
operator =(const LineSearch & other)93 LineSearch &LineSearch::operator =(const LineSearch &other)
94 {
95   if (this != &other) {
96     LocalOptimizer::operator =(other);
97     CopyAttributes(other);
98   }
99   return *this;
100 }
101 
102 // -----------------------------------------------------------------------------
~LineSearch()103 LineSearch::~LineSearch()
104 {
105 }
106 
107 // =============================================================================
108 // Parameters
109 // =============================================================================
110 
111 // -----------------------------------------------------------------------------
Set(const char * name,const char * value)112 bool LineSearch::Set(const char *name, const char *value)
113 {
114   // Number of line search iterations
115   if (strcmp(name, "Maximum no. of line search iterations")    == 0 ||
116       strcmp(name, "Maximum number of line search iterations") == 0 ||
117       strcmp(name, "Maximum no. of line iterations")           == 0 ||
118       strcmp(name, "Maximum number of line iterations")        == 0 ||
119       strcmp(name,    "No. of line search iterations")         == 0 ||
120       strcmp(name, "Number of line search iterations")         == 0 ||
121       strcmp(name,    "No. of line iterations")                == 0 ||
122       strcmp(name, "Number of line iterations")                == 0) {
123     return FromString(value, _NumberOfIterations) && _NumberOfIterations > 0;
124   // Minimum length of step
125   } else if (strcmp(name, "Minimum length of steps") == 0) {
126     return FromString(value, _MinStepLength) && _MinStepLength > .0;
127   // Maximum length of step
128   } else if (strcmp(name, "Maximum length of steps") == 0) {
129     return FromString(value, _MaxStepLength) && _MaxStepLength > .0;
130   // Length of steps
131   } else if (strcmp(name, "Length of steps") == 0) {
132     if (FromString(value, _MaxStepLength) && _MaxStepLength > .0) {
133       _MinStepLength = _MaxStepLength;
134     }
135   // Range of step lengths (e.g., "[0.1 1]" or "0.1 1")
136   } else if (strcmp(name, "Range of steps")                  == 0 ||
137              strcmp(name, "Range of step lengths")           == 0 ||
138              strcmp(name, "Min-/Maximum length of steps")    == 0 ||
139              strcmp(name, "Min-/maximum length of steps")    == 0 ||
140              strcmp(name, "Minimum/Maximum length of steps") == 0 ||
141              strcmp(name, "Minimum/maximum length of steps") == 0) {
142     istringstream ss(value);
143     char c1, c2;
144     string str;
145     ss >> c1;
146     if (c1 != '[') ss.putback(c1);
147     ss >> str;
148     if (!FromString(str.c_str(), _MinStepLength) || _MinStepLength < .0) return false;
149     ss >> str;
150     if (!FromString(str.c_str(), _MaxStepLength) || _MaxStepLength < _MinStepLength) return false;
151     ss >> c2;
152     if (c1 == '[' && c2 != ']') return false;
153     return true;
154   }
155   return false;
156 }
157 
158 // -----------------------------------------------------------------------------
Parameter() const159 ParameterList LineSearch::Parameter() const
160 {
161   ParameterList params;
162   Insert(params, "Maximum no. of line iterations", ToString(_NumberOfIterations));
163   if (_MinStepLength == _MaxStepLength) {
164     Insert(params, "Length of steps", ToString(_MaxStepLength));
165   } else {
166     Insert(params, "Minimum length of steps", ToString(_MinStepLength));
167     Insert(params, "Maximum length of steps", ToString(_MaxStepLength));
168   }
169   return params;
170 }
171 
172 // =============================================================================
173 // Optimization
174 // =============================================================================
175 
176 // -----------------------------------------------------------------------------
Initialize()177 void LineSearch::Initialize()
178 {
179   // Initialize base class
180   LocalOptimizer::Initialize();
181 
182   // Check parameters
183   if (_MinStepLength <= 0.) {
184     Throw(ERR_InvalidArgument, __FUNCTION__, "Minimum length of steps must be positive, but is ", _MinStepLength);
185   }
186   if (_MaxStepLength < _MinStepLength) {
187     Throw(ERR_InvalidArgument, __FUNCTION__, "Maximum length of steps must be greater or equal minimum length, but is ", _MaxStepLength);
188   }
189 
190   // Set initial step length
191   _StepLength = _MaxStepLength;
192 }
193 
194 
195 } // namespace mirtk
196