29 #include <vnl/vnl_math.h> 31 #include <itkTimeProbe.h> 45 NLOPTOptimizers::NLOPTOptimizers()
49 m_ErrorCode = NLOPT_FAILURE;
50 m_Algorithm = NLOPT_LN_BOBYQA;
51 m_LocalOptimizer = NLOPT_LN_BOBYQA;
61 m_PopulationSize = -1;
62 m_VectorStorageSize = -1;
74 NLOPTOptimizers::~NLOPTOptimizers()
95 double NLOPTOptimizers::NloptFunctionWrapper(
unsigned n,
const double *x,
double *grad,
void *data)
104 for (
unsigned int i=0; i<n ; i++ )
105 itkCurrentPosition[i] = x[i]/optimizer->GetScales()[i];
107 double f = (optimizer->GetCostFunction()->
GetValue(itkCurrentPosition));
115 DerivativeType derivative;
116 optimizer->GetCostFunction()->GetDerivative (itkCurrentPosition, derivative);
117 for (
unsigned int i=0; i<n ; i++ )
118 grad[i] = derivative[i];
134 void NLOPTOptimizers::StartOptimization()
136 if( m_CostFunction.IsNull() )
137 throw itk::ExceptionObject(__FILE__, __LINE__,
"Error. The cost function has not been set",
"NLOPTOptimizers");
142 unsigned int n = m_CostFunction->GetNumberOfParameters();
148 if ( this->GetInitialPosition().GetSize() != n )
149 throw itk::ExceptionObject(__FILE__, __LINE__,
"Invalid initial position parameter. Its size should be equal to the number of parameters",
"NLOPTOptimizers");
150 m_CurrentPosition.set_size(n);
152 double *x =
new double[n];
153 const double *in_x = GetInitialPosition().data_block();
154 for (
unsigned int i=0; i<n ; i++ )
156 x[i] = in_x[i] * GetScales()[i];
157 m_CurrentPosition[i] = in_x[i];
166 if ( m_LowerBoundParameters.GetSize()!=0 )
168 if ( m_LowerBoundParameters.GetSize()!=n )
169 throw itk::ExceptionObject(__FILE__, __LINE__,
"Invalid lower bound parameter. Its size should be equal to the number of parameters",
"NLOPTOptimizers");
172 for (
unsigned int i=0; i<n; i++ )
173 lb[i] = m_LowerBoundParameters[i]*GetScales()[i];
175 if ( m_UpperBoundParameters.GetSize()!=0 )
177 if ( m_UpperBoundParameters.GetSize()!=n )
178 throw itk::ExceptionObject(__FILE__, __LINE__,
"Invalid upper bound parameter. Its size should be equal to the number of parameters",
"NLOPTOptimizers");
181 for (
unsigned int i=0; i<n; i++ )
182 ub[i] = m_UpperBoundParameters[i]*GetScales()[i];
188 m_NloptOptions = nlopt_create((::nlopt_algorithm)(
int)m_Algorithm, n);
189 nlopt_set_lower_bounds(m_NloptOptions, lb);
190 nlopt_set_upper_bounds(m_NloptOptions, ub);
192 nlopt_set_max_objective(m_NloptOptions, (nlopt_func)this->NloptFunctionWrapper, (
void *)
this);
194 nlopt_set_min_objective(m_NloptOptions, (nlopt_func)this->NloptFunctionWrapper, (
void *)
this);
196 if ( m_StopValSet ) nlopt_set_stopval(m_NloptOptions, m_StopVal);
198 nlopt_set_ftol_rel(m_NloptOptions, m_FTolRel);
199 nlopt_set_ftol_abs(m_NloptOptions, m_FTolAbs);
200 nlopt_set_xtol_rel(m_NloptOptions, m_XTolRel);
201 nlopt_set_xtol_abs1(m_NloptOptions, m_XTolAbs);
202 nlopt_set_maxtime(m_NloptOptions, m_MaxTime);
203 nlopt_set_vector_storage(m_NloptOptions, m_VectorStorageSize);
204 nlopt_set_maxeval(m_NloptOptions, m_MaxEval);
205 nlopt_set_population(m_NloptOptions, m_PopulationSize);
207 for (
unsigned int i = 0;i < m_InequalityConstraints.size();++i)
208 nlopt_add_inequality_constraint(m_NloptOptions, (nlopt_func)ConstraintsFunctionType::GetConstraintValue, m_InequalityConstraints[i]->GetAdditionalData(), m_InequalityConstraints[i]->GetTolerance());
210 for (
unsigned int i = 0;i < m_EqualityConstraints.size();++i)
211 nlopt_add_equality_constraint(m_NloptOptions, (nlopt_func)ConstraintsFunctionType::GetConstraintValue, m_EqualityConstraints[i]->GetAdditionalData(), m_EqualityConstraints[i]->GetTolerance());
217 m_NloptLocalOptions = nlopt_create((::nlopt_algorithm)(
int)m_LocalOptimizer, n);
219 if ( m_StopValSet ) nlopt_set_stopval(m_NloptLocalOptions, m_StopVal);
221 nlopt_set_ftol_rel(m_NloptLocalOptions, m_FTolRel);
222 nlopt_set_ftol_abs(m_NloptLocalOptions, m_FTolAbs);
223 nlopt_set_xtol_rel(m_NloptLocalOptions, m_XTolRel);
224 nlopt_set_xtol_abs1(m_NloptLocalOptions, m_XTolAbs);
225 nlopt_set_maxtime(m_NloptLocalOptions, m_MaxTime);
226 nlopt_set_vector_storage(m_NloptLocalOptions, m_VectorStorageSize);
227 nlopt_set_maxeval(m_NloptLocalOptions, m_MaxEval);
232 nlopt_set_local_optimizer(m_NloptOptions, m_NloptLocalOptions);
234 this->InvokeEvent( itk::StartEvent() );
240 m_ErrorCode =
static_cast<nlopt_result
>((int)nlopt_optimize (m_NloptOptions, x, &valf));
241 SetCurrentCost(valf);
247 for (
unsigned int i=0; i<n ; i++ )
248 p[i] = x[i]/GetScales()[i];
249 this->SetCurrentPosition(p);
255 if ( ub!=NULL )
delete[] ub;
256 if ( lb!=NULL )
delete[] lb;
258 nlopt_destroy(m_NloptOptions);
259 nlopt_destroy(m_NloptLocalOptions);
261 this->InvokeEvent( itk::EndEvent() );
275 bool NLOPTOptimizers::isSuccessful()
const 278 && m_ErrorCode != NLOPT_ROUNDOFF_LIMITED
279 && m_ErrorCode != NLOPT_FORCED_STOP)
295 std::string NLOPTOptimizers::GetErrorCodeDescription()
const 297 switch ( m_ErrorCode )
302 case NLOPT_STOPVAL_REACHED:
303 return "Optimization stopped because 'StopVal' was reached.";
304 case NLOPT_FTOL_REACHED:
305 return "Optimization stopped because 'FTolRel' or 'FTolAbs' was reached.";
306 case NLOPT_XTOL_REACHED:
307 return "Optimization stopped because 'XTolRel' or 'XTolAbs' was reached.";
308 case NLOPT_MAXEVAL_REACHED:
309 return "Optimization stopped because 'MaxEval' was reached.";
310 case NLOPT_MAXTIME_REACHED:
311 return "Optimization stopped because 'MaxTime' was reached.";
315 return "Generic failure code.";
316 case NLOPT_INVALID_ARGS:
317 return "Invalid arguments (e.g. lower bounds are bigger than upper bounds, an unknown algorithm was specified, etcetera).";
318 case NLOPT_OUT_OF_MEMORY:
319 return "Out of memory.";
320 case NLOPT_ROUNDOFF_LIMITED:
321 return "Halted because roundoff errors limited progress.";
322 case NLOPT_FORCED_STOP:
323 return "Halted because of a forced termination.";
328 return "Unkown error code";
348 m_InequalityConstraints.push_back(constraint);
359 void NLOPTOptimizers::ClearInequalityConstraints()
361 m_InequalityConstraints.clear();
383 m_EqualityConstraints.push_back(constraint);
394 void NLOPTOptimizers::ClearEqualityConstraints()
396 m_EqualityConstraints.clear();
411 m_LowerBoundParameters.SetSize(p.GetSize());
412 for (
unsigned int i=0; i<p.GetSize(); i++ )
413 m_LowerBoundParameters[i] = p[i];
428 m_UpperBoundParameters.SetSize(p.GetSize());
429 for (
unsigned int i=0; i<p.GetSize(); i++ )
430 m_UpperBoundParameters[i] = p[i];
436 void NLOPTOptimizers::PrintSelf(std::ostream& os, itk::Indent indent)
const 438 Superclass::PrintSelf(os,indent);
440 int major, minor, bugfix;
441 nlopt_version(&major, &minor, &bugfix);
442 os << indent <<
"NLOPT v"<<major<<
"."<<minor<<
"."<<bugfix<<std::endl;
444 os << indent <<
"Algorithm: " << ::nlopt_algorithm_name((::nlopt_algorithm)(
int)m_Algorithm) << std::endl;
445 os << indent <<
"Maximize On/Off " << m_Maximize << std::endl;
447 os << indent <<
"StopVal: "<<m_StopVal<<std::endl;
448 os << indent <<
"FTolRel: "<<m_FTolRel<<std::endl;
449 os << indent <<
"FTolAbs: "<<m_FTolAbs<<std::endl;
450 os << indent <<
"XTolRel: "<<m_XTolRel<<std::endl;
451 os << indent <<
"MaxEval: "<<m_MaxEval<<std::endl;
452 os << indent <<
"MaxTime: "<<m_MaxTime<<std::endl;
Implements an ITK wrapper for the NLOPT library.
MeasureType GetValue() const
SingleValuedNonLinearOptimizer::ParametersType ParametersType