[SciPy-dev] BUG gradient in leastsq (http://projects.scipy.org/scipy/scipy/ticket/416)
Andy Jennings
scipy2mdjhs78c@jenningsstory....
Mon Jul 30 00:58:01 CDT 2007
On 7/26/07, Alan G Isaac <aisaac@american.edu> wrote:
> On Thu, 26 Jul 2007, dmitrey apparently wrote:
> > 1. What's the difference between these 2 funcs from __minpack.h:
> > int jac_multipack_calling_function(int *n, double *x, double *fvec,
> > double *fjac, int *ldfjac, int *iflag)
> > int jac_multipack_lm_function(int *m, int *n, double *x, double *fvec,
> > double *fjac, int *ldfjac, int *iflag)
>
> > They have same description.
>
> > /* This is the function called from the Fortran code it should
> > -- use call_python_function to get a multiarrayobject result
> > -- check for errors and return -1 if any
> > -- otherwise place result of calculation in *fvec or *fjac.
>
> > If iflag = 1 this should compute the function.
> > If iflag = 2 this should compute the jacobian (derivative matrix)
> > */
>
> > 2. So patch assigned to the ticket proposes to rewrite the line 152
> > MATRIXC2F(fjac, result_array->data, *n, *ldfjac)
> > as
> > MATRIXC2F(fjac, result_array->data, *ldfjac, *n)
>
> > however, line 92 is the same. so maybe it needs same patch.
>
> > 3. The MATRIXC2F is defined in minpack.h w/o any description:
> > #define MATRIXC2F(jac,data,n,m) {double *p1=(double *)(jac), *p2,
> > *p3=(double *)(data);\
> > int i,j;\
> > for (j=0;j<(m);p3++,j++) \
> > for (p2=p3,i=0;i<(n);p2+=(m),i++,p1++) \
> > *p1 = *p2; }
>
> > I have no idea what does it do.
> > So I replaced (jac,data,n,m) by (jac,data,m,n), and user's example works
> > correctly for all cases 1-3:
> > 1) w/o gradient info
> > 2) with gradient info, col_deriv=0
> > 3) with gradient info, col_deriv=1 (in the case I modified the user's
> > gradient func so that it returns transposed gradient)
>
> > scipy.test(1) also didn't yield any bugs related to leastsq.
>
>
>
> Andy, and others, can you comment on this?
> Also, I would hope for Travis to comment
> before committing these changes, since he
> contributed the code.
>
> Here is the situation:
> http://projects.scipy.org/scipy/scipy/ticket/416
> Andy identified a bug and a possible fix.
> (But see details above.)
> Neither Dmitrey nor I are familiar with this code.
> Dmitrey is willing to apply the patch if there is support
> for his doing so, but it will be done "mechancially".
> In any case, he will add a unit test exposing the problem.
>
> Cheers,
> Alan Isaac
>
>
>
>
>
Thanks for looking at this bug. The fix looks fine to me.
I don't really have an opinion on the jac_multipack_calling_function
question. I think it's likely that it has the same issue, but on the
chance that it doesn't, you hate to risk breaking something that's
working correctly.
More information about the Scipy-dev
mailing list