Re: QR factorization


In article <f78df7e8-ead9-40ed-a846-78fb0972cf6d@xxxxxxxxxxxxxxxxxxxxxxxxxxxx>,
patricia.andrusca@xxxxxxxxx writes:
Hello!

If I have a matrix A, which I factorize as QR and then I want I want
to extract the matrix Q, how do I do this? What is the algorithm?
Thank you!
depends on your software:

in matlab you get it already explicitly
[Q,R]=qr(A)

if you have a fortran or c-code
(e.g. from lapack or linpack or NR ): look how they store the information
about Q. The most usual way is :
don't store Q explicitly but the information used for it in the original
matrix: often the strict upper triangle of R is stored above the diagonal
in the original A, the diagonal of R being stored in a separate vector.
Q is obtained as a product of n (or n-1 if row number=column number for A)
Householder matrices: these have the form

H(i)= I - (2/(u'*u))*u*u' , u=u(i)

u(i) is stored on and below the diagonal of A since its first i-1 entries are
zero. The process is
H(n)*...*H(1)*A = R

H(i) annihilates all elements in the current column i below the diagonal.
(A being overwritten by R and the u(i) ). Hence your Q is

H(1)*...*H(n)

hence you can obtain the explicit Q if desired (mostly it is not needed
at all) by multiplying the H(i) from the right in the identity:

Q=H(1)
for i=2 to n
Q=Q*H(i)

In doing so it is necessary to use the special form of the H(i).
This gives

Q=Q - (2/u(i)'*u(i)) * (Q*u(i)) *u(i)'

and in computing Q*u(i) you also make use of the fact that the first i-1
components of u(i) are zero. Hence this "matrixmultiply" is indeed a
simple rank-one change. finally you must look whether the u(i) are
already normalized (such that the scalar becomes simply 2) or whether even
the scalar 2/(u(i)'*u(i)) has already been totally absorbed in the stored u(i)
in order to save work.

hth
peter

.

Quantcast