Re: problem with gsl_eigen_hermv(), possibly due to bug in gsl_linalg_complex_householder_transform()

[Brian Gough <bjg at network-theory dot co dot uk>]

 > On Tue, 15 Oct 2002 17:26:32 -0500 (CDT)
 > Steve Martin <spm@zippy.physics.niu.edu> wrote:
 > > I think there may be a bug in gsl 1.2 which affects the computation of
 > > eigenvectors of Hermitian matrices using gsl_eigen_hermv(), but quite
 > > possibly other things too.

Here is a patch to gsl/linalg/householdercomplex.c which should fix
the problem.

Index: householdercomplex.c
===================================================================
RCS file: /cvs/gsl/gsl/linalg/householdercomplex.c,v
retrieving revision 1.5
diff -c -r1.5 householdercomplex.c
*** householdercomplex.c	19 Nov 2001 22:32:06 -0000	1.5
--- householdercomplex.c	17 Oct 2002 16:25:02 -0000
***************
*** 42,54 ****
        double beta_r = - (GSL_REAL(alpha) >= 0 ? +1 : -1) * absa ;
  
        gsl_complex tau;
-       GSL_REAL(tau) = (beta_r - GSL_REAL(alpha)) / beta_r ;
-       GSL_IMAG(tau) = - GSL_IMAG(alpha) / beta_r ;
  
!       {
!         gsl_complex beta = gsl_complex_rect (beta_r, 0.0);
!         gsl_vector_complex_set (v, 0, beta) ;
!       }
        
        return tau;
      }
--- 42,63 ----
        double beta_r = - (GSL_REAL(alpha) >= 0 ? +1 : -1) * absa ;
  
        gsl_complex tau;
  
!       if (beta_r == 0.0)
!         {
!           GSL_REAL(tau) = 0.0;
!           GSL_IMAG(tau) = 0.0;
!         }
!       else 
!         {
!           GSL_REAL(tau) = (beta_r - GSL_REAL(alpha)) / beta_r ;
!           GSL_IMAG(tau) = - GSL_IMAG(alpha) / beta_r ;
! 
!           {
!             gsl_complex beta = gsl_complex_rect (beta_r, 0.0);
!             gsl_vector_complex_set (v, 0, beta) ;
!           }
!         }
        
        return tau;
      }