
Go to the first, previous, next, last section, table of contents.
One dimensional MinimizationThis chapter describes routines for finding minima of arbitrary onedimensional 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, with 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 header file `gsl_min.h' contains prototypes for the minimization functions and related declarations. To use the minimization algorithms to find the maximum of a function simply invert its sign. OverviewThe minimization algorithms begin with a bounded region known to contain a minimum. The region is described by an lower bound a and an upper bound b, with an estimate of the location of the minimum x. The value of the function at x must be less than the value of the function at the ends of the interval, f(a) > f(x) < f(b) This condition guarantees that a minimum is contained somewhere within the interval. On each iteration a new point x' is selected using one of the available algorithms. If the new point is a better estimate of the minimum, f(x') < f(x), then the current estimate of the minimum x is updated. The new point also allows the size of the bounded interval to be reduced, by choosing the most compact set of points which satisfies the constraint f(a) > f(x) < f(b). The interval is reduced until it encloses the true minimum to a desired tolerance. This provides a best estimate of the location of the minimum and a rigorous error estimate. Several bracketing algorithms are available within a single framework. The user provides a highlevel driver for the algorithm, and the library provides the individual functions necessary for each of the steps. There are three main phases of the iteration. The steps are,
The state for the minimizers is held in a CaveatsNote that minimization functions can only search for one minimum at a time. When there are several 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 minimum in an area where there is more than one. With all minimization algorithms it can be difficult to determine the location of the minimum to full numerical precision. The behavior of the function in the region of the minimum x^* can be approximated by a Taylor expansion, y = f(x^*) + (1/2) f"(x^*) (x  x^*)^2 and the second term of this expansion can be lost when added to the first term at finite precision. This magnifies the error in locating x^*, making it proportional to \sqrt \epsilon (where \epsilon is the relative accuracy of the floating point numbers). For functions with higher order minima, such as x^4, the magnification of the error is correspondingly worse. The best that can be achieved is to converge to the limit of numerical accuracy in the function values, rather than the location of the minimum itself. Initializing the Minimizer
Providing the function to minimize
You must provide a continuous function of one variable for the
minimizers to operate on. In order to allow for general parameters the
functions are defined by a IterationThe following functions drive the iteration of each algorithm. Each function performs one iteration to update the state of any minimizer of the corresponding type. The same functions work for all minimizers so that different methods can be substituted at runtime without modifications to the code.
The minimizer maintains a current best estimate of the position of the minimum at all times, and the current interval bounding the minimum. This information can be accessed with the following auxiliary functions,
Stopping ParametersA minimization procedure should stop when one of the following conditions is true:
The handling of these conditions is under user control. The function below allows the user to test the precision of the current result.
Minimization AlgorithmsThe minimization algorithms described in this section require an initial interval which is guaranteed to contain a minimum  if a and b are the endpoints of the interval and x is an estimate of the minimum then f(a) > f(x) < f(b). This ensures that the function has at least one minimum somewhere in the interval. If a valid initial interval is used then these algorithm cannot fail, provided the function is wellbehaved.
ExamplesThe following program uses the Brent algorithm to find the minimum of the function f(x) = \cos(x) + 1, which occurs at x = \pi. The starting interval is (0,6), with an initial guess for the minimum of 2. #include <stdio.h> #include <gsl/gsl_errno.h> #include <gsl/gsl_math.h> #include <gsl/gsl_min.h> double fn1 (double x, void * params) { return cos(x) + 1.0; } int main (void) { int status; int iter = 0, max_iter = 100; const gsl_min_fminimizer_type *T; gsl_min_fminimizer *s; double m = 2.0, m_expected = M_PI; double a = 0.0, b = 6.0; gsl_function F; F.function = &fn1; F.params = 0; T = gsl_min_fminimizer_brent; s = gsl_min_fminimizer_alloc (T); gsl_min_fminimizer_set (s, &F, m, a, b); printf ("using %s method\n", gsl_min_fminimizer_name (s)); printf ("%5s [%9s, %9s] %9s %10s %9s\n", "iter", "lower", "upper", "min", "err", "err(est)"); printf ("%5d [%.7f, %.7f] %.7f %+.7f %.7f\n", iter, a, b, m, m  m_expected, b  a); do { iter++; status = gsl_min_fminimizer_iterate (s); m = gsl_min_fminimizer_x_minimum (s); a = gsl_min_fminimizer_x_lower (s); b = gsl_min_fminimizer_x_upper (s); status = gsl_min_test_interval (a, b, 0.001, 0.0); if (status == GSL_SUCCESS) printf ("Converged:\n"); printf ("%5d [%.7f, %.7f] " "%.7f %.7f %+.7f %.7f\n", iter, a, b, m, m_expected, m  m_expected, b  a); } while (status == GSL_CONTINUE && iter < max_iter); return status; } Here are the results of the minimization procedure. bash$ ./a.out 0 [0.0000000, 6.0000000] 2.0000000 1.1415927 6.0000000 1 [2.0000000, 6.0000000] 3.2758640 +0.1342713 4.0000000 2 [2.0000000, 3.2831929] 3.2758640 +0.1342713 1.2831929 3 [2.8689068, 3.2831929] 3.2758640 +0.1342713 0.4142862 4 [2.8689068, 3.2831929] 3.2758640 +0.1342713 0.4142862 5 [2.8689068, 3.2758640] 3.1460585 +0.0044658 0.4069572 6 [3.1346075, 3.2758640] 3.1460585 +0.0044658 0.1412565 7 [3.1346075, 3.1874620] 3.1460585 +0.0044658 0.0528545 8 [3.1346075, 3.1460585] 3.1460585 +0.0044658 0.0114510 9 [3.1346075, 3.1460585] 3.1424060 +0.0008133 0.0114510 10 [3.1346075, 3.1424060] 3.1415885 0.0000041 0.0077985 Converged: 11 [3.1415885, 3.1424060] 3.1415927 0.0000000 0.0008175 References and Further ReadingFurther information on Brent's algorithm is available in the following book,
Go to the first, previous, next, last section, table of contents. 