  The QZ algorithm for computing eigenvalues and eigenvectors of a matrix
pencil $A - \lambda B$ requires that the matrices first be reduced to
Hessenberg-triangular (HT) form. The current method of choice for HT reduction
relies entirely on Givens rotations regrouped and accumulated into small dense
matrices which are subsequently applied using matrix multiplication routines. A
non-vanishing fraction of the total flop count must nevertheless still be
performed as sequences of overlapping Givens rotations alternately applied from
the left and from the right. The many data dependencies associated with this
computational pattern leads to inefficient use of the processor and poor
scalability.
  In this paper, we therefore introduce a fundamentally different approach that
relies entirely on (large) Householder reflectors partially accumulated into
block reflectors, by using (compact) WY representations. Even though the new
algorithm requires more floating point operations than the state of the art
algorithm, extensive experiments on both real and synthetic data indicate that
it is still competitive, even in a sequential setting. The new algorithm is
conjectured to have better parallel scalability, an idea which is partially
supported by early small-scale experiments using multi-threaded BLAS. The
design and evaluation of a parallel formulation is future work.
