0

consider the following example from the GNU website on minimzation. Suppose instead of the given objective fn1, I need to interpolate the objective from a set of points (xa,ya) to obtain my objective function double fn2(double x){gsl_interp_eval( myinterp, xa[],ya[],x). The question is how to avoid having to setup the interpolation myinterp object inside the objective function fn2 each time this is called and rather pass it to the objective as an argument. I've tried with a struct but no luck.

I guess a more general way to ask is to say "how to pass an object of type gsl_interp to a function of type gsl_function?"

#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\n",
              iter, a, b,
              m, m - m_expected, b - a);
    }
  while (status == GSL_CONTINUE && iter < max_iter);

  gsl_min_fminimizer_free (s);

  return status;
}
Florian Oswald
  • 5,054
  • 5
  • 30
  • 38

1 Answers1

0

sorry i made a mistake, hopefully not wasted anybody's time. i declared in the struct something like

struct myint = { double a; gsl_interp int;};

i.e. my name of the gsl_interp object was int. that's not very clever since int is of course a type declaration.

Florian Oswald
  • 5,054
  • 5
  • 30
  • 38