27 #include <queso/Defines.h>
28 #include <queso/GslVector.h>
29 #include <queso/VectorSpace.h>
30 #include <queso/ScalarFunction.h>
31 #include <queso/GslOptimizer.h>
32 #include <queso/OptimizerMonitor.h>
34 #include <gsl/gsl_multimin.h>
35 #include <gsl/gsl_blas.h>
43 double c_evaluate(
const gsl_vector * x,
void * context) {
51 for (
unsigned int i = 0; i < state.sizeLocal(); i++) {
52 state[i] = gsl_vector_get(x, i);
70 gsl_vector * derivative) {
79 for (
unsigned int i = 0; i < state.sizeLocal(); i++) {
80 state[i] = gsl_vector_get(x, i);
90 for (
unsigned int i = 0; i < deriv.sizeLocal(); i++) {
91 gsl_vector_set(derivative, i, GSL_NAN);
106 bool userComputedDerivative =
true;
107 for (
unsigned int i = 0; i < deriv.sizeLocal(); i++) {
110 if (gsl_isnan(deriv[i])) {
111 userComputedDerivative =
false;
116 if (userComputedDerivative) {
117 for (
unsigned int i = 0; i < deriv.sizeLocal(); i++) {
118 gsl_vector_set(derivative, i, -deriv[i]);
126 for (
unsigned int i = 0; i < deriv.sizeLocal(); i++) {
127 double tempState = state[i];
136 state[i] = tempState;
139 if (!gsl_isnan(fx) && !gsl_isnan(fxph)) {
140 gsl_vector_set(derivative, i, (fxph - fx) / h);
143 gsl_vector_set(derivative, i, GSL_NAN);
152 double * f, gsl_vector * derivative) {
162 m_objectiveFunction(objectiveFunction),
163 m_initialPoint(new
GslVector(objectiveFunction.domainSet().
164 vectorSpace().zeroVector())),
165 m_minimizer(new
GslVector(this->m_objectiveFunction.domainSet().
166 vectorSpace().zeroVector())),
167 m_solver_type(BFGS2),
168 m_fstep_size(this->m_objectiveFunction.domainSet().vectorSpace().zeroVector()),
169 m_fdfstep_size(getFdfstepSize()),
170 m_line_tol(getLineTolerance())
186 m_objectiveFunction(objectiveFunction),
187 m_initialPoint(new
GslVector(objectiveFunction.domainSet().
188 vectorSpace().zeroVector())),
189 m_minimizer(new
GslVector(this->m_objectiveFunction.domainSet().
190 vectorSpace().zeroVector())),
191 m_solver_type(BFGS2),
192 m_fstep_size(this->m_objectiveFunction.domainSet().vectorSpace().zeroVector()),
193 m_fdfstep_size(getFdfstepSize()),
194 m_line_tol(getLineTolerance())
218 std::cerr <<
"Minimization was given initial point outside of domain"
225 zeroVector().sizeLocal();
249 for (
unsigned int i = 0; i < initialPoint.
sizeLocal(); i++) {
268 bool gradient_needed =
false;
278 gradient_needed =
true;
294 return gradient_needed;
300 gsl_vector * x = gsl_vector_alloc(dim);
301 for (
unsigned int i = 0; i <
dim; i++) {
306 const gsl_multimin_fdfminimizer_type* type = NULL;
312 type = gsl_multimin_fdfminimizer_conjugate_fr;
315 type = gsl_multimin_fdfminimizer_conjugate_pr;
318 type = gsl_multimin_fdfminimizer_vector_bfgs;
321 type = gsl_multimin_fdfminimizer_vector_bfgs2;
324 type = gsl_multimin_fdfminimizer_steepest_descent;
335 gsl_multimin_fdfminimizer * solver =
336 gsl_multimin_fdfminimizer_alloc(type, dim);
339 gsl_multimin_function_fdf minusLogPosterior;
340 minusLogPosterior.n =
dim;
344 minusLogPosterior.params = (
void *)(
this);
353 status = gsl_multimin_fdfminimizer_iterate(solver);
358 std::cerr <<
"Error while GSL does optimisation. "
359 <<
"See below for GSL error type." << std::endl;
360 std::cerr <<
"Gsl error: " << gsl_strerror(status) << std::endl;
365 status = gsl_multimin_test_gradient(solver->gradient, this->getTolerance());
369 gsl_vector* x = gsl_multimin_fdfminimizer_x(solver);
370 std::vector<double> x_min(dim);
371 for(
unsigned int i = 0; i <
dim; i++)
372 x_min[i] = gsl_vector_get(x,i);
374 double f = gsl_multimin_fdfminimizer_minimum(solver);
376 gsl_vector* grad = gsl_multimin_fdfminimizer_gradient(solver);
377 double grad_norm = gsl_blas_dnrm2(grad);
379 monitor->
append( x_min, f, grad_norm );
384 for (
unsigned int i = 0; i <
dim; i++) {
385 (*m_minimizer)[i] = gsl_vector_get(solver->x, i);
389 gsl_multimin_fdfminimizer_free(solver);
398 gsl_vector* x = gsl_vector_alloc(dim);
399 for (
unsigned int i = 0; i <
dim; i++) {
404 const gsl_multimin_fminimizer_type* type = NULL;
409 type = gsl_multimin_fminimizer_nmsimplex;
412 type = gsl_multimin_fminimizer_nmsimplex2;
415 type = gsl_multimin_fminimizer_nmsimplex2rand;
428 gsl_multimin_fminimizer* solver =
429 gsl_multimin_fminimizer_alloc(type, dim);
432 gsl_multimin_function minusLogPosterior;
433 minusLogPosterior.n =
dim;
435 minusLogPosterior.params = (
void *)(
this);
438 gsl_vector* step_size = gsl_vector_alloc(dim);
440 for(
unsigned int i = 0; i <
dim; i++) {
444 gsl_multimin_fminimizer_set(solver, &minusLogPosterior, x, step_size);
453 status = gsl_multimin_fminimizer_iterate(solver);
458 std::cerr <<
"Error while GSL does optimisation. "
459 <<
"See below for GSL error type." << std::endl;
460 std::cerr <<
"Gsl error: " << gsl_strerror(status) << std::endl;
465 size = gsl_multimin_fminimizer_size(solver);
467 status = gsl_multimin_test_size (size, this->
getTolerance());
471 gsl_vector* x = gsl_multimin_fminimizer_x(solver);
472 std::vector<double> x_min(dim);
473 for(
unsigned int i = 0; i <
dim; i++)
474 x_min[i] = gsl_vector_get(x,i);
476 double f = gsl_multimin_fminimizer_minimum(solver);
478 monitor->
append( x_min, f, size );
485 for (
unsigned int i = 0; i <
dim; i++) {
486 (*m_minimizer)[i] = gsl_vector_get(solver->x, i);
490 gsl_vector_free(step_size);
491 gsl_multimin_fminimizer_free(solver);
513 if( solver == std::string(
"fletcher_reeves_cg") )
515 else if( solver == std::string(
"polak_ribiere_cg") )
517 else if( solver == std::string(
"bfgs") )
519 else if( solver == std::string(
"bfgs2") )
521 else if( solver == std::string(
"steepest_decent") ) {
525 else if( solver == std::string(
"steepest_descent") )
527 else if( solver == std::string(
"nelder_mead") )
529 else if( solver == std::string(
"nelder_mead2") )
531 else if( solver == std::string(
"nelder_mead2_rand") )
537 std::cerr <<
"Error: Invalid GslOptimizer solver name: " << solver << std::endl
538 <<
" Valids choices are: fletcher_reeves_cg" << std::endl
539 <<
" polak_ribiere_cg" << std::endl
540 <<
" bfgs" << std::endl
541 <<
" bfgs2" << std::endl
542 <<
" steepest_descent" << std::endl
543 <<
" nelder_mead" << std::endl
544 <<
" nelder_mead2" << std::endl
545 <<
" nelder_mead2_rand" << std::endl;
573 fstepSizeVector.
cwSet(fstepSize);
virtual void setFstepSize(double fstepSize)
Sets the step size to use in gradient-free solvers.
void set_solver_type(SolverType solver)
void append(std::vector< double > &x_min, double objective, double norm)
Add current value of minimizer, objective, and error norm.
SolverType string_to_enum(std::string &solver)
A templated (base) class for handling scalar functions.
void setInitialPoint(const GslVector &intialPoint)
Set the point at which the optimization starts.
ScopedPtr< OptimizerOptions >::Type m_optionsObj
A base class for handling optimisation of scalar functions.
GslOptimizer(const BaseScalarFunction< GslVector, GslMatrix > &objectiveFunction)
Constructs an object that will maximize a scalar function.
#define queso_deprecated()
virtual void minimize(OptimizerMonitor *monitor=NULL)
Minimize the objective function, starting at m_initialPoint.
virtual void setLineTolerance(double lineTolerance)
Sets the tolerance to use for line minimisation.
Class for vector operations using GSL library.
void minimize_with_gradient(unsigned int dim, OptimizerMonitor *monitor)
This class provides options for a Optimizer.
void cwSet(double value)
Component-wise sets all values to this with value.
bool solver_needs_gradient(SolverType solver)
Helper function.
const GslVector & minimizer() const
Return the point that minimizes the objective function.
virtual ~GslOptimizer()
Destructor.
virtual double getLineTolerance() const
Gets the tolerance to use for line minimisation.
double getTolerance() const
Returns the tolerance used to test for an extremum in the optimizer.
const BaseScalarFunction< GslVector, GslMatrix > & m_objectiveFunction
A base class for handling optimisation of scalar functions.
GslVector m_fstep_size
For use in gradient-free algorithms.
void c_evaluate_derivative(const gsl_vector *x, void *context, gsl_vector *derivative)
virtual void setFdfstepSize(double fdfstepSize)
Sets the step to use in gradient-based solvers.
double c_evaluate(const gsl_vector *x, void *context)
double m_fdfstep_size
For use in gradient-based algorithms.
unsigned int sizeLocal() const
Returns the length of this vector.
void c_evaluate_with_derivative(const gsl_vector *x, void *context, double *f, gsl_vector *derivative)
void set_step_size(const GslVector &step_size)
Sets step size used in gradient-free solvers.
void minimize_no_gradient(unsigned int dim, OptimizerMonitor *monitor)
Object to monitor convergence of optimizers.
virtual double getFdfstepSize() const
Gets the step to use in gradient-based solvers.
virtual std::string getSolverType() const
Gets the algorithm to use for minimisation.
virtual void setSolverType(std::string solverType)
Sets the algorithm to use for minimisation.
const BaseScalarFunction< GslVector, GslMatrix > & objectiveFunction() const
Returns the objective function.
double getFiniteDifferenceStepSize() const
Returns the step size used in the finite difference formula.
virtual double getFstepSize() const
Gets the step size to use in gradient-free solvers.
unsigned int getMaxIterations() const
Returns the maximum number of iterations the optimizer will do.
GslVector * m_initialPoint