Re: [eigen] How to obtain the matrices Q and R of a QR decomposition

[Gael Guennebaud <gael.guennebaud@xxxxxxxxx>]

Hi,

here it works well:

#include <iostream>
#include <Eigen/QR>

using namespace Eigen;

int main (int argc, char* argv[]){
  int m = 4;
  int n = 3;

  MatrixXcd A (m,n), Q, R;
  A.setRandom();


  std::cout << "A = \n" << A << "\n\n";
  Eigen::HouseholderQR<MatrixXcd> qr(A);
  Q = qr.householderQ();
  R = qr.matrixQR().triangularView<Upper>();
  std::cout << Q << "\n\n" << R << "\n\n" << Q * R - A << "\n";
}

give me:

A =
 (-0.211234,0.680375)   (0.10794,-0.444451)   (0.434594,0.271423)
   (0.59688,0.566198) (0.257742,-0.0452059)  (0.213938,-0.716795)
 (-0.604897,0.823295) (0.0268018,-0.270431) (-0.514226,-0.967399)
 (0.536459,-0.329554)    (0.83239,0.904459)  (0.608354,-0.725537)

     (-0.130389,0.419978)      (-0.219276,0.243524)
(0.467691,0.399364)       (0.498566,0.269141)
      (0.368438,0.349499)      (-0.381022,0.120666)
(0.226617,-0.199732)      (-0.636019,0.294273)
     (-0.373387,0.508198)    (-0.153642,0.00610515)
(-0.422318,-0.58557)      (0.226197,0.0791447)
     (0.331142,-0.203425)     (-0.626771,-0.562486)
(-0.000392371,-0.0963378)     (0.358459,-0.0442264)

            (1.62003,0)    (-0.177362,0.462072) (-0.0649518,-0.0507344)
                  (0,0)             (-1.2718,0)    (-0.0972986,1.03066)
                  (0,0)                   (0,0)             (1.35661,0)
                  (0,0)                   (0,0)                   (0,0)

           (5.55112e-17,0)           (-2.77556e-17,0)            (0,1.66533e-16)
          (0,-1.11022e-16)                      (0,0)            (8.32667e-17,0)
 (2.22045e-16,1.11022e-16)           (-5.55112e-17,0)                      (0,0)
           (1.11022e-16,0) (-2.22045e-16,1.11022e-16) (-1.11022e-16,1.11022e-16)


gael

On Thu, Nov 11, 2010 at 10:40 AM, Sebastian Birk
<birk@xxxxxxxxxxxxxxxxxxxxx> wrote:
> Hi,
>
> using Eigen 3 I would expect the following code to give me Q and R:
>
>    MatrixXcd A = MatrixXcd::Random(m,m), Q, R;
>    HouseholderQR<MatrixXcd> hh(A);
>
>    // MATLAB economy size QR decomp.: [Q,R] = qr( A, 0 );
>    Q = hh.householderQ();
>    R = hh.matrixQR().triangularView<Upper>();
>
> But when I take a look (via cout) at Q and R they both are exactly the
> same as A.
> Where is the flaw in this approach?
>
> Thanks,
> Sebastian
>
>