gsl_linalg_complex_householder_transform
Mario Pernici
Mario.Pernici@mi.infn.it
Thu May 29 20:40:00 GMT 2003
Hi Brian,
< I think I used the LAPACK method described in Lapack Working Note 72
thank you for the reference!
Now my examples for the functions in the section
Tridiagonal Decomposition of Hermitian Matrices
work without modifications of the code of gsl.
The key point is that
gsl_complex
gsl_linalg_complex_householder_transform (gsl_vector_complex * v)
returns a complex number tau such that
(I - tau^* v v^H) x = -/+ ||x|| e_1
(see Working Note 72)
while I thought the relation was (I - tau v v^H) x = -/+ ||x|| e_1
Maybe this should be added in the comments to the source code of
this function.
And you could add to householdercomplex.c the function
gsl_linalg_complex_householder_hv in analogy with householder.c .
======================
int
gsl_linalg_complex_householder_hv(gsl_complex tau, const
gsl_vector_complex * v, gsl_vector_complex * w)
{
size_t i;
gsl_complex z;
if (GSL_REAL(tau) == 0.0 && GSL_IMAG(tau) == 0.0)
return GSL_SUCCESS;
/* z = v^H w */
gsl_blas_zdotc(v, w, &z);
/* w = w - tau * (v^H w) * v */
z = gsl_complex_mul(z, tau);
for (i = 0; i < v->size; i++)
{
gsl_complex vi, wi;
vi = gsl_vector_complex_get(v,i);
wi = gsl_vector_complex_get(w,i);
vi = gsl_complex_mul(vi, z);
wi = gsl_complex_sub(wi, vi);
gsl_vector_complex_set(w, i, wi);
}
}
=============================
Bye
Mario
More information about the Gsl-discuss
mailing list