Full Q out of QR

[jsievers]

Hi - does anyone know if it's possible to get the full Q matrix out of QR for a rectangular input matrix? The Lapack user's guide makes it seem possible: "...associated routines are provided to work with Q: xORGQR can generate all or part of Q," (node40) where earlier in the page it's made clear that Q refers to a square matrix. If A is rectangular, then the User's Guide says that Q is made up of Q1 and Q2, where Q1 is the orthogonal basis that spans A and Q2 spans the rest of the space. I need Q2, but if I read the documentation on dorgqr correctly, it seems to only want to make Q1. I suspect it's possible, given the wording in the User's Guide (and also the fact that MATLAB defaults to producing the full square Q for a rectangular A), but can't figure out how to do it. If anyone knows how to get the full Q, I'd really appreciate any info. Or, if not, than any tips on efficient ways to get an orthogonal basis spanning the complement of the space spanned by a rectangular matrix would also be appreciated (one could always fill out A with zeros to make it square and call, e.g., SVD, but that seems rather wastefull).

Thanks very much,

Jon

[Julien Langou]

Hello Jon,

Does anyone know if it's possible to get the full Q matrix out of QR for a
rectangular input matrix?

Instead of using xORGQR/xUNGQR, you can use xORMQR/xUNMQR where C is the
identity matrix (M-by-M), SIDE='L' and TRANS='N', this will compute C <- Q * I. So this
will compute Q.

If you just want Q2, then just provide C, the (M-N) last columns of the M-by-M
identity matrix.

You can also compute Q1 by this mean: you provide xORMQR/xUNMQR with C the N
first columns of the M-by-M identity matrix. However in this particular case,
you are better off using xORGQR/xUNGQR (which justifies the existence of the
xORGQR/xUNGQR). Two good reasons: first if you use xORGQR, Q1 overwrites A so
you can save memory space. Second, xORGQR utilizes the fact that "C" is the
identity matrix to reduce the number fo flops. To give a rough idea, if A is
square (highest gain for ORGQR versus ORMQR), xORGQR uses 4/3*N^3 flops while
xORMQR uses 4*N^3 flops.

The Lapack user's guide makes it seem possible: "...associated routines are
provided to work with Q: xORGQR can generate all or part of Q," (node40) where
earlier in the page it's made clear that Q refers to a square matrix. If A is
rectangular, then the User's Guide says that Q is made up of Q1 and Q2, where
Q1 is the orthogonal basis that spans A and Q2 spans the rest of the space. I
need Q2, but if I read the documentation on dorgqr correctly, it seems to only
want to make Q1.

Well, you are right, we'll need to fix the documentation. The LUG's explanation of xORGQR/xUNGQR let me think as well that the xORGQR/xUNGQR provides a full QR
decomposition. (See: http://www.netlib.org/lapack/lug/node40.html.)

Julien

[jsievers]

I'll give it a go - thanks very much!

Jon

[sven]

Dear Jon and Julien,

Actually you can obtain the full Q from xORGRQ.

If A is an m by p matrix with m >= p then The whole of Q may be computed by the call:

CALL xORGQR (M,M,P,A,LDA,TAU,WORK,LWORK,INFO)

You need to make sure that A has been declared with at least m columns as it will be overwritten by the m by m matrix Q.

You might find the NAG documentation for F08AFF at:

http://www.nag.co.uk/numeric/FL/manual/ ... conts.html

helpful. Julien is correct about the best way to obtain Q2.

Best wishes,

Sven Hammarling.

[vdehdari]

Hi Julien,

I used your method (using DORMQR) because I'm using fortran.
My code is like this:

Code: Select all

! Matrix A is Mxp

SIDE='L'
TRANS='N'

CALL DGEQRF (M, p, A, max(1,M), TAU, WORK, max(1,p), INFO)

Call DORMQR( SIDE, TRANS, M, p, min(M,p), A, max(1,M), TAU, C, max(1,M), WORK, max(1,p), INFO )

do IROW = 1, M
  write(*,*) (C(IROW,ICOL),ICOL = 1, M)
enddo

but C matrix is a MxM matrix with zero in all of its array.
Do you know what is my problem?

Thanks