0
$\begingroup$

I'm using the subroutine sgeqr2 from Lapack. This subroutine solves the QR-factorization

$$A = QR$$

It's easy to find the $R$ matrix, because the in-out argument $A$ of subroutine sgeqr2 will return a matrix, where the upper values from the diagonal is the $R$-values. Eeasy.

But how can I find the $Q$ matrix? According to the "manual".

The matrix Q is represented as a product of elementary reflectors

     Q = H(1) H(2) . . . H(k), where k = min(m,n).

  Each H(i) has the form

     H(i) = I - tau * v * v**T

  where tau is a real scalar, and v is a real vector with
  v(1:i-1) = 0 and v(i) = 1; v(i+1:m) is stored on exit in A(i+1:m,i),
  and tau in TAU(i).

The conclusion here is to take the vector $tau$ and pick the first value of $tau$. Then multiply it with a matrix $vv^T$.

Have I interprent this text correct?

For a matrix $A$

0.674878,   0.151285,   0.875139,   0.150518,
0.828102,   0.150747,   0.934674,   0.474325,
0.476510,   0.914686,   0.740681,   0.060455,
0.792594,   0.471488,   0.529343,   0.743405,
0.084739,   0.475160,   0.419307,   0.628999,
0.674878,   0.151285,   0.875139,   0.150518

I get the $R$ matrix

-1.568159 -0.751743 -1.762103 -0.808132 
0.000000 0.887756 0.233565 0.241302 
0.000000 0.000000 0.500422 -0.142971 
0.000000 0.000000 0.000000 -0.700355 
0.000000 0.000000 0.000000 0.000000 
0.000000 0.000000 0.000000 0.000000

And $tau$ vector

1.430363 1.205732 1.007224 1.655577 

How should I do, to find $Q$ matrix if I know the dimension of $A$ and the values from $tau$?

I wrote som MALAB/Octave code to find the $Q$-matrix from $H$, but it won't work.

tau = [1.430363 1.205732 1.007224 1.655577];

H = eye(4); % Initial 

for i = 1:4
  v(1:i-1) = 0;
  v(i) = 1;
  H = H*(eye(4) - tau(i)*v'*v);
  H
end
$\endgroup$
  • $\begingroup$ You may use *ORMQR function from Lapack. You should also use *geqrf instead of *geqr2. $\endgroup$ – Pawel Kowal Feb 17 at 15:17
0
$\begingroup$

Now I got the answer!

Assume that we have a matrix $A$

0.674878,   0.151285,   0.875139,   0.150518,
0.828102,   0.150747,   0.934674,   0.474325,
0.476510,   0.914686,   0.740681,   0.060455,
0.792594,   0.471488,   0.529343,   0.743405,
0.084739,   0.475160,   0.419307,   0.628999,
0.674878,   0.151285,   0.875139,   0.150518

And we implement it into the subroutine.

Then we get the matrix $A$ back, but it has change into this

-1.568159 -0.751743 -1.762103 -0.808132 
0.369188 0.887756 0.233565 0.241302 
0.212440 -0.675308 0.500422 -0.142971 
0.353358 -0.142374 0.875625 -0.700355 
0.037779 -0.412039 -0.411443 0.439228 
0.300877 0.112496 -0.222825 -0.122951 

The $R$ matrix is the upper triangular above the diagonal.

-1.568159 -0.751743 -1.762103 -0.808132 
0.000000 0.887756 0.233565 0.241302 
0.000000 0.000000 0.500422 -0.142971 
0.000000 0.000000 0.000000 -0.700355 
0.000000 0.000000 0.000000 0.000000 
0.000000 0.000000 0.000000 0.000000 

And then the lower traingular under the diagonal, we want to use those in this loop.

Q = eye(6); % Initial 
v = [ 0 0 0 0 0 0];
m = 6;
for i = 1:4
  v(1:i-1) = 0;
  v(i) = 1;
  v(i+1:m) = A(i+1:m,i);
  A(i+1:m,i)
  Q = Q*(eye(6) - tau(i)*v'*v);
end

Where tau is this vector

1.430363 1.205732 1.007224 1.655577

Then our $Q$ matrix will be.

  -0.430363  -0.194015   0.323945   0.148698  -0.562826  -0.577340
  -0.528073  -0.277360   0.137758  -0.191610   0.739060  -0.205691
  -0.303866   0.773024   0.049331   0.520577   0.184471  -0.051341
  -0.505430   0.103108  -0.770067  -0.285532  -0.234751   0.065334
  -0.054038   0.489479   0.419171  -0.752685  -0.118974   0.033112
  -0.430363  -0.194015   0.323945   0.148698  -0.183642   0.785092

I post code instead of matrix math, because it will be easier to read.

$\endgroup$
0
$\begingroup$

Your implementation is very inefficient.

Consider

Q = Q*(eye(6) - tau(i)*v'*v);

This would create the matrix (eye(6) - tau(i)*v'*v) and then multiply it by Q.

If your matrices are n x n, then this would require $ 2 n^3 $ computations for that multiply.

If you instead update Q with

Q = Q - tau(i) * ( Q * v' ) * v then the computation is $ O( n^2 ) $ per iteration.

There are other inefficiencies in what you do: You have a lot of zeroes in your vector v that you could take advantage of.

But, of course, in your example you are only working with a $ 6 x 6 $ matrix, so it doesn't really matter so much.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.