
Go to the first, previous, next, last section, table of contents.
Multidimensional MinimizationThis chapter describes routines for finding minima of arbitrary multidimensional functions. The library provides low level components for a variety of iterative minimizers and convergence tests. These can be combined by the user to achieve the desired solution, while providing full access to the intermediate steps of the algorithms. Each class of methods uses the same framework, so that you can switch between minimizers at runtime without needing to recompile your program. Each instance of a minimizer keeps track of its own state, allowing the minimizers to be used in multithreaded programs. The minimization algorithms can be used to maximize a function by inverting its sign. The header file `gsl_multimin.h' contains prototypes for the minimization functions and related declarations. OverviewThe problem of multidimensional minimization requires finding a point x such that the scalar function, f(x_1, ..., x_n) takes a value which is lower than at any neighboring point. For smooth functions the gradient g = \nabla f vanishes at the minimum. In general there are no bracketing methods available for the minimization of ndimensional functions. All algorithms proceed from an initial guess using a search algorithm which attempts to move in a downhill direction. A onedimensional line minimisation is performed along this direction until the lowest point is found to a suitable tolerance. The search direction is then updated with local information from the function and its derivatives, and the whole process repeated until the true ndimensional minimum is found. Several minimization algorithms are available within a single framework. The user provides a highlevel driver for the algorithms, and the library provides the individual functions necessary for each of the steps. There are three main phases of the iteration. The steps are,
Each iteration step consists either of an improvement to the
linemimisation in the current direction or an update to the search
direction itself. The state for the minimizers is held in a
CaveatsNote that the minimization algorithms can only search for one local minimum at a time. When there are several local minima in the search area, the first minimum to be found will be returned; however it is difficult to predict which of the minima this will be. In most cases, no error will be reported if you try to find a local minimum in an area where there is more than one. It is also important to note that the minimization algorithms find local minima; there is no way to determine whether a minimum is a global minimum of the function in question. Initializing the Multidimensional MinimizerThe following function initializes a multidimensional minimizer. The minimizer itself depends only on the dimension of the problem and the algorithm and can be reused for different problems.
Providing a function to minimizeYou must provide a parametric function of n variables for the minimizers to operate on. You also need to provide a routine which calculates the gradient of the function and a third routine which calculates both the function value and the gradient together. In order to allow for general parameters the functions are defined by the following data type:
The following example function defines a simple paraboloid with two parameters, /* Paraboloid centered on (dp[0],dp[1]) */ double my_f (const gsl_vector *v, void *params) { double x, y; double *dp = (double *)params; x = gsl_vector_get(v, 0); y = gsl_vector_get(v, 1); return 10.0 * (x  dp[0]) * (x  dp[0]) + 20.0 * (y  dp[1]) * (y  dp[1]) + 30.0; } /* The gradient of f, df = (df/dx, df/dy). */ void my_df (const gsl_vector *v, void *params, gsl_vector *df) { double x, y; double *dp = (double *)params; x = gsl_vector_get(v, 0); y = gsl_vector_get(v, 1); gsl_vector_set(df, 0, 20.0 * (x  dp[0])); gsl_vector_set(df, 1, 40.0 * (y  dp[1])); } /* Compute both f and df together. */ void my_fdf (const gsl_vector *x, void *params, double *f, gsl_vector *df) { *f = my_f(x, params); my_df(x, params, df); } The function can be initialized using the following code, gsl_multimin_function_fdf my_func; double p[2] = { 1.0, 2.0 }; /* center at (1,2) */ my_func.f = &my_f; my_func.df = &my_df; my_func.fdf = &my_fdf; my_func.n = 2; my_func.params = (void *)p; IterationThe following function drives the iteration of each algorithm. The function performs one iteration to update the state of the minimizer. The same function works for all minimizers so that different methods can be substituted at runtime without modifications to the code.
Stopping CriteriaA minimization procedure should stop when one of the following conditions is true:
The handling of these conditions is under user control. The functions below allow the user to test the precision of the current result.
AlgorithmsThere are several minimization methods available. The best choice of algorithm depends on the problem. Each of the algorithms uses the value of the function and its gradient at each evaluation point.
ExamplesThis example program finds the minimum of the paraboloid function defined earlier. The location of the minimum is offset from the origin in x and y, and the function value at the minimum is nonzero. The main program is given below, it requires the example function given earlier in this chapter. int main (void) { size_t iter = 0; int status; const gsl_multimin_fdfminimizer_type *T; gsl_multimin_fdfminimizer *s; /* Position of the minimum (1,2). */ double par[2] = { 1.0, 2.0 }; gsl_vector *x; gsl_multimin_function_fdf my_func; my_func.f = &my_f; my_func.df = &my_df; my_func.fdf = &my_fdf; my_func.n = 2; my_func.params = ∥ /* Starting point, x = (5,7) */ x = gsl_vector_alloc (2); gsl_vector_set (x, 0, 5.0); gsl_vector_set (x, 1, 7.0); T = gsl_multimin_fdfminimizer_conjugate_fr; s = gsl_multimin_fdfminimizer_alloc (T, 2); gsl_multimin_fdfminimizer_set (s, &my_func, x, 0.01, 1e4); do { iter++; status = gsl_multimin_fdfminimizer_iterate (s); if (status) break; status = gsl_multimin_test_gradient (s>gradient, 1e3); if (status == GSL_SUCCESS) printf ("Minimum found at:\n"); printf ("%5d %.5f %.5f %10.5f\n", iter, gsl_vector_get (s>x, 0), gsl_vector_get (s>x, 1), s>f); } while (status == GSL_CONTINUE && iter < 100); gsl_multimin_fdfminimizer_free (s); gsl_vector_free (x); return 0; } The initial stepsize is chosen as 0.01, a conservative estimate in this case, and the line minimization parameter is set at 0.0001. The program terminates when the norm of the gradient has been reduced below 0.001. The output of the program is shown below, x y f 1 4.99629 6.99072 687.84780 2 4.98886 6.97215 683.55456 3 4.97400 6.93501 675.01278 4 4.94429 6.86073 658.10798 5 4.88487 6.71217 625.01340 6 4.76602 6.41506 561.68440 7 4.52833 5.82083 446.46694 8 4.05295 4.63238 261.79422 9 3.10219 2.25548 75.49762 10 2.85185 1.62963 67.03704 11 2.19088 1.76182 45.31640 12 0.86892 2.02622 30.18555 Minimum found at: 13 1.00000 2.00000 30.00000 Note that the algorithm gradually increases the step size as it successfully moves downhill, as can be seen by plotting the successive points. The conjugate gradient algorithm finds the minimum on its second direction because the function is purely quadratic. Additional iterations would be needed for a more complicated function. References and Further ReadingA brief description of multidimensional minimization algorithms and further references can be found in the following book,
Go to the first, previous, next, last section, table of contents. 