[Bug-gsl] gsl_multifit_linear() covariances incorrect due to sorted singular values?

[Doug Rogers]

I've noticed that the singular values returned by SVD are sorted,
destroying the relationship that they had with their column vector. This
means, I believe, that the covariances calculated in
gsl_multifit_wlinear(), et al., are not in the proper order - they
correspond to different solution coefficients.

I've noticed that the sorting code at the end of gsl_linalg_SV_decomp()
is quite independent of the main body and could easily be moved or
performed as a separate step when enabled.

I removed the sorting (in a separate routine), and sure enough the
covariance matrix changed, swapping rows and columns along with the
singular values.

My modification to svd.c was something like this, to maintain
compatibility with the old name:

gsl_linalg_SV_decomp_unsorted(A, S, V, work)  /* New! */
{
  /*
   * All but the final sorting code.
   */
}

gsl_linalg_SV_decomp(A, S, V, work)
{
  int result = gsl_linalg_SV_decomp_unsorted(A, S, V, work);

  if (result != 0)
  {
    return result;
  }

  /*
   * Sorting code here as before, but all the other code
   * copied into gsl_linalg_SV_decomp_unsorted().
   */
}

I have a recursive patch with this change, which passes the unit test. I
added a prototype for gsl_linalg_SV_decomp_unsorted() in gsl_linalg.h. I
did the same for gsl_linalg_SV_decomp_mod() since that would allow it to
be used by gsl_multifit_linear(). I noticed that gsl_multifit_linear()
has no unit test.

Some questions:

Is this something that could be added for 1.9?

Is there a development list or forum for GSL? I didn't see anything
mentioned for developers at the GNU GSL page.

Would it be better to create a work structure with attributes for the
SVD API?

Thanks,

Doug Rogers

-- 
Innovative Concepts, Inc. www.innocon.com 703-893-2007 x220