Re: [eigen] notes on Householder and Givens

[Benoit Jacob <jacob.benoit.1@xxxxxxxxx>]

Something just occured to me.

If x_0 is "zero", that is muchSmallerThan(norm(x)), then there is
nothing wrong with this formula:

v = x + norm(x) e_0

Indeed, no cancellation can occur since x0 is much smaller than
norm(x), and after that v0 is guaranteed to be approximately equal to
norm(x).

So, i'm sorry to not have thought of that earlier, but it seems that
no pivoting is needed: i think that the following
ei_compute_householder will work for all nonzero vectors x even if x0
= 0:

void ei_compute_householder(const Vector& x of size n, Householder *h)
{
 // compute natural householder vector v
 Vector v = x;
 if(ei_isMuchSmallerThan(x0, norm(x))
 {
   v0 += norm(x);
 }
 else
 {
   v0 += ( x_0 / abs(x_0) ) norm(x);
 }

 Vector v_starting_with_1 = v/v0;

 h->essential = v_starting_with_1.end(n-1);
 h->beta = 2 / v_starting_with_1.squaredNorm();
};

Notice that this is still very different from Golub's code since our
code takes good care of the case x0=0 while they don't.

Cheers,
Benoit