[Bug-gsl] multimin example (BUG#9)

[Fabio Brugnara]

Dear GSL maintainers,

I have to solve a nonlinear optimization problem and, not wanting to
reinvent the wheel, I looked around for a suitable numerical computation
library. In the end I chose to use the GSL, as it looked quite reasonably
structured and featured for my needs. Unfortunately, I was immediately
presented with a bug, at the first trial. The problem is that the very same
example program given in the manual does not produce the expected results.
I've then seen that this was a "known bug", BUG#9. Actually, the behavior
is incorrect both for the Fletcher-Reeves and the Polak-Ribiere minimizers.
Looking at the code and comparing it to other descriptions of the methods,
I think the error is in the following code fragment of
conjugate_[fr,pr]_iterate():

  if (fc < fa)
    {
      /* Success, reduced the function value */
      state->step = stepc * 2.0;
      *f = fc;
      gsl_vector_memcpy (x, x1);
      GSL_MULTIMIN_FN_EVAL_DF (fdf, x1, gradient);
      state->g0norm = gsl_blas_dnrm2 (gradient);       <<<<<<<< PROBLEM HERE
      return GSL_SUCCESS;
    }

g0norm will be used to find the new search direction but, to my
understanding, g0norm should keep the value given at the last
change of direction. Actually, simply removing the marked statement
make both methods work as expected, at least on the example.
The same statement appears in vector_bfgs_iterate(), but in that method
g0norm is not actually used.

I'm not an expert, so please verify if this is really the problem, and
that this fix does not introduce other problems.

Best regards, and thank you for offering this software,
Fabio Brugnara

PS: I'm referring to GSL 1.3 on a Red Hat Linux system, but the problem is
    not related to the operating system or compiler.