The WY representation for products of Householder matrices.

*(English)*Zbl 0628.65033The authors present a block Householder algorithm for obtaining the QR factorization of an \(m\times n\) real matrix A. Their procedure is suitable for computers with an architecture that favors matrix-matrix multiplications such as FPS-164/MAX, (FPS stands for Floating Point Systems), which was used by the authors. The standard Householder orthogonalization to obtain the QR factorization is given by the algorithm

For \(k=1\) to n

For \(j=1\) to k-1

\(a_ k:=a_ k+(u^ T_ ja_ k)v_ j\)

end for j

Compute \(v_ k\) such that \(\| v_ k\| =1\) and \(P_ k=I+(-2v_ k)v^ T_ k=I+u_ kv^ T_ k\)

is the Householder matrix that will zero components \(k+1\) to m of the

column \(a_ k\)

end for k.

This algorithm is rich in inner products \((u^ T_ ja_ k)\) and the operation \(y:=\alpha x+y\) but mut not matrix multiplications. The new way devised by the authors hinges on the fact that the Householder product \(Q_ k=P_ 1...P_ k\) can be written as \[ (1)\quad Q_ k=I+W_ kY^ T_ k=I+[W_{k-1},Q_{k-1}u_ k][Y_{k-1},v_ k] \] where \(W_ 1=u_ 1\), \(Y_ 1=v_ 1.\)

Let \(A(r:s,t:u)\) be the block of A defined by rows r to s and columns t to u. Let A be partitioned in column blocks \(A=[A_ 1,...,A_ n]\) where \(A_ i\) has p columns with the exception of \(A_ n\) which may have less than p columns and p is an integer dictated by the architecture of the computer. The authors use the following algorithm to overwrite A with R

for \(k=1\) to n

\(s:=(k-1)p+1\)

Compute [using the procedure outlined below] W and Y such that

\(I+WY^ T\) is orthogonal and \((I+WY^ T)A(i:m\), \(s:s+p-1)\) is upper triangular

\(A(s:m,s:n):=(I+WY^ T)A(s:m,s:n)\)

end for k.

At the beginning of step k, A is overwritten by \(\begin{pmatrix} R_{11} & R_{12} & R_{13} \\ 0 & A_ k & B \end{pmatrix}\). The block \([A_ kB]\) has \(m-(k-1)p\) rows, \(A_ k\) has p columns and B has \(1+(k-1)p\) columns. The procedure for modifying \([A_ kB]\) has the following steps:

1) calculate the QR factorization of \(A_ k\) using the WY procedure given in (1),

2) calculate \(Z^ T=W^ TB\), 3) overwrite B with \(B+YZ^ T\).

For \(k=1\) to n

For \(j=1\) to k-1

\(a_ k:=a_ k+(u^ T_ ja_ k)v_ j\)

end for j

Compute \(v_ k\) such that \(\| v_ k\| =1\) and \(P_ k=I+(-2v_ k)v^ T_ k=I+u_ kv^ T_ k\)

is the Householder matrix that will zero components \(k+1\) to m of the

column \(a_ k\)

end for k.

This algorithm is rich in inner products \((u^ T_ ja_ k)\) and the operation \(y:=\alpha x+y\) but mut not matrix multiplications. The new way devised by the authors hinges on the fact that the Householder product \(Q_ k=P_ 1...P_ k\) can be written as \[ (1)\quad Q_ k=I+W_ kY^ T_ k=I+[W_{k-1},Q_{k-1}u_ k][Y_{k-1},v_ k] \] where \(W_ 1=u_ 1\), \(Y_ 1=v_ 1.\)

Let \(A(r:s,t:u)\) be the block of A defined by rows r to s and columns t to u. Let A be partitioned in column blocks \(A=[A_ 1,...,A_ n]\) where \(A_ i\) has p columns with the exception of \(A_ n\) which may have less than p columns and p is an integer dictated by the architecture of the computer. The authors use the following algorithm to overwrite A with R

for \(k=1\) to n

\(s:=(k-1)p+1\)

Compute [using the procedure outlined below] W and Y such that

\(I+WY^ T\) is orthogonal and \((I+WY^ T)A(i:m\), \(s:s+p-1)\) is upper triangular

\(A(s:m,s:n):=(I+WY^ T)A(s:m,s:n)\)

end for k.

At the beginning of step k, A is overwritten by \(\begin{pmatrix} R_{11} & R_{12} & R_{13} \\ 0 & A_ k & B \end{pmatrix}\). The block \([A_ kB]\) has \(m-(k-1)p\) rows, \(A_ k\) has p columns and B has \(1+(k-1)p\) columns. The procedure for modifying \([A_ kB]\) has the following steps:

1) calculate the QR factorization of \(A_ k\) using the WY procedure given in (1),

2) calculate \(Z^ T=W^ TB\), 3) overwrite B with \(B+YZ^ T\).

Reviewer: I.Imam

##### MSC:

65F25 | Orthogonalization in numerical linear algebra |

65F05 | Direct numerical methods for linear systems and matrix inversion |

15A23 | Factorization of matrices |

65Y05 | Parallel numerical computation |