123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172 |
- /* multimin/steepest_descent.c
- *
- * Copyright (C) 1996, 1997, 1998, 1999, 2000 Fabrice Rossi
- *
- * This program is free software; you can redistribute it and/or modify
- * it under the terms of the GNU General Public License as published by
- * the Free Software Foundation; either version 3 of the License, or (at
- * your option) any later version.
- *
- * This program is distributed in the hope that it will be useful, but
- * WITHOUT ANY WARRANTY; without even the implied warranty of
- * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
- * General Public License for more details.
- *
- * You should have received a copy of the GNU General Public License
- * along with this program; if not, write to the Free Software
- * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
- */
- /* steepest_descent.c -- the steepest descent algorithm */
- /* Modified by Brian Gough to use single iteration structure */
- #include "gsl__config.h"
- #include "gsl_multimin.h"
- #include "gsl_blas_types.h"
- #include "gsl_blas.h"
- typedef struct
- {
- double step;
- double max_step;
- double tol;
- gsl_vector *x1;
- gsl_vector *g1;
- }
- steepest_descent_state_t;
- static int
- steepest_descent_alloc (void *vstate, size_t n)
- {
- steepest_descent_state_t *state = (steepest_descent_state_t *) vstate;
- state->x1 = gsl_vector_alloc (n);
- if (state->x1 == NULL)
- {
- GSL_ERROR ("failed to allocate space for x1", GSL_ENOMEM);
- }
- state->g1 = gsl_vector_alloc (n);
- if (state->g1 == NULL)
- {
- gsl_vector_free (state->x1);
- GSL_ERROR ("failed to allocate space for g1", GSL_ENOMEM);
- }
- return GSL_SUCCESS;
- }
- static int
- steepest_descent_set (void *vstate, gsl_multimin_function_fdf * fdf,
- const gsl_vector * x, double *f,
- gsl_vector * gradient, double step_size, double tol)
- {
- steepest_descent_state_t *state = (steepest_descent_state_t *) vstate;
- GSL_MULTIMIN_FN_EVAL_F_DF (fdf, x, f, gradient);
- state->step = step_size;
- state->max_step = step_size;
- state->tol = tol;
- return GSL_SUCCESS;
- }
- static void
- steepest_descent_free (void *vstate)
- {
- steepest_descent_state_t *state = (steepest_descent_state_t *) vstate;
- gsl_vector_free (state->x1);
- gsl_vector_free (state->g1);
- }
- static int
- steepest_descent_restart (void *vstate)
- {
- steepest_descent_state_t *state = (steepest_descent_state_t *) vstate;
- state->step = state->max_step;
- return GSL_SUCCESS;
- }
- static int
- steepest_descent_iterate (void *vstate, gsl_multimin_function_fdf * fdf,
- gsl_vector * x, double *f,
- gsl_vector * gradient, gsl_vector * dx)
- {
- steepest_descent_state_t *state = (steepest_descent_state_t *) vstate;
- gsl_vector *x1 = state->x1;
- gsl_vector *g1 = state->g1;
- double f0 = *f;
- double f1;
- double step = state->step, tol = state->tol;
- int failed = 0;
- /* compute new trial point at x1= x - step * dir, where dir is the
- normalized gradient */
- double gnorm = gsl_blas_dnrm2 (gradient);
- if (gnorm == 0.0)
- {
- gsl_vector_set_zero (dx);
- return GSL_ENOPROG;
- }
- trial:
- gsl_vector_set_zero (dx);
- gsl_blas_daxpy (-step / gnorm, gradient, dx);
- gsl_vector_memcpy (x1, x);
- gsl_blas_daxpy (1.0, dx, x1);
- /* evaluate function and gradient at new point x1 */
- GSL_MULTIMIN_FN_EVAL_F_DF (fdf, x1, &f1, g1);
- if (f1 > f0)
- {
- /* downhill step failed, reduce step-size and try again */
- failed = 1;
- step *= tol;
- goto trial;
- }
- if (failed)
- step *= tol;
- else
- step *= 2.0;
- state->step = step;
- gsl_vector_memcpy (x, x1);
- gsl_vector_memcpy (gradient, g1);
- *f = f1;
- return GSL_SUCCESS;
- }
- static const gsl_multimin_fdfminimizer_type steepest_descent_type =
- { "steepest_descent", /* name */
- sizeof (steepest_descent_state_t),
- &steepest_descent_alloc,
- &steepest_descent_set,
- &steepest_descent_iterate,
- &steepest_descent_restart,
- &steepest_descent_free
- };
- const gsl_multimin_fdfminimizer_type
- * gsl_multimin_fdfminimizer_steepest_descent = &steepest_descent_type;
|