forum.alglib.net
http://forum.alglib.net/

QR decomposition, how to unpack Q in full
http://forum.alglib.net/viewtopic.php?f=2&t=4
Page 1 of 2

Author:  tam [ Mon May 10, 2010 3:38 pm ]
Post subject:  QR decomposition, how to unpack Q in full

Hi,

First of all, thanks a lot for letting me use your ALGLIB. It has helped me a lot for my math project.

However, your note on this website said that you can unpack the Q from the QR decomposition in full by using RMatrixQRUnpackQ but when I use RMatrixQRUnpackQ, I still do not get the full Q but only partial of the Q. Example, if I have a matrix size 5000(myRow) by 16(myCol). When I unpack the Q, I get the matrix size 5000 by 5000 but the value is correct for the first 16 columns of each row. I want to see the full Q, which mean, each row of Q will have all 5000 values. Could you show me how do I unpack Q in full using your ALGLIB? Thanks a lot.

Here is partial of my code in C++. Note: minRowCol is the minimum number between myRow and myCol:


ap::real_2d_array a;
ap::real_1d_array tau;
ap::real_2d_array q;
ap::real_2d_array r;

a.setlength(myRow, myCol);
tau.setlength(myRow);
q.setlength(myRow, myRow);
r.setlength(myRow, myCol);

//some code to assign values to a

rmatrixqr(q,myRow,myCol,tau);
rmatrixqrunpackq(a,myRow,myCol,tau,myRow,q);

Author:  Sergey.Bochkanov [ Tue May 11, 2010 4:14 am ]
Post subject:  Re: QR decomposition, how to unpack Q in full

Hello!

tam wrote:
Hi,
Example, if I have a matrix size 5000(myRow) by 16(myCol). When I unpack the Q, I get the matrix size 5000 by 5000 but the value is correct for the first 16 columns of each row.

What do you mean when you say that other 4984 columns are "incorrect"? The only requirement they have to satisfy is orthogonality condition. Actually, they never influence result of Q*R multiplication.

Author:  tam [ Tue May 11, 2010 1:09 pm ]
Post subject:  Re: QR decomposition, how to unpack Q in full

The other values are all zeros. When I mean 'first 16 of each row', I mean the 5000 rows of the first 16 columns values are correct but the other 5000 rows of the rest other 4984 columns are just zero. Your function, rmatrixQRUnpackQ, said it will just 'partial unpack Q' and didn't mention that it will fully unpack Q but then I read your website, you have mentioned that you can either partial or full unpack Q so I was just wondering how can I do the "full" unpack.

The reason I know the other values are wrong is because I unpacked this Q in matlab and I got all nonzero values. I also do this unpack using GNU libraries (gsl_linalg_QR_unpack) and it is also the same results as matlab, which is nonzero values. I just thought I will just give a try on your library since I your library code run fast but I didn't know how.

Author:  Sergey.Bochkanov [ Tue May 11, 2010 1:54 pm ]
Post subject:  Re: QR decomposition, how to unpack Q in full

tam wrote:
//some code to assign values to a
rmatrixqr(q,myRow,myCol,tau);
rmatrixqrunpackq(a,myRow,myCol,tau,myRow,q);

I've made QR decomposition of a 5000x16 matrix and unpacked "full" Q. All columns (before 16th and after) were non-zero. It works as exprected.

Then I've noticed small error in your example: you initialize elements of a, then you use q (whose elements are uninitialized) as first parameter of rmatrixqr(), and then you use a (whose elements are uninitialized too) as first parameter of rmatrixqrunpackq(). If your actual code looks like this, then probably your result is just a consequence of your operations with uninitialized matrices.

Correct code will look as follows:
Quote:
//some code to assign values to a
rmatrixqr(a,myRow,myCol,tau);
rmatrixqrunpackq(a,myRow,myCol,tau,myRow,q);


Am I right? If no, i.e. if your actual code is correct and doesn't contain this error, please, send me actual code and your data (this 5000x16 matrix).

Author:  tam [ Tue May 11, 2010 2:23 pm ]
Post subject:  Re: QR decomposition, how to unpack Q in full

Thanks for pointing the error out but I doubled check my code and my code was right. I must mistype when I posted here. Anyway, my code and sample data are included here for your check to see if I am doing it right. My A matrix in the code is the matrix sample I am included.

Hey, I cannot attach text files to my reply. Why? Do you have an email address?

Author:  Sergey.Bochkanov [ Tue May 11, 2010 3:51 pm ]
Post subject:  Re: QR decomposition, how to unpack Q in full

My e-mail is sergey.bochkanov at alglib.net

Please, send me your matrix, I'll check what's wrong with it.

P.S. That's very strange that you can't attach files, you should have all permissions necessary to do it. I'll try to solve this issue too :)

Author:  tam [ Tue May 11, 2010 4:01 pm ]
Post subject:  Re: QR decomposition, how to unpack Q in full

Just sent you email from my gmail. Hopefully it won't go straight to your junk mail :-)

Author:  admin [ Tue May 11, 2010 4:08 pm ]
Post subject:  Re: QR decomposition, how to unpack Q in full

OK, *.txt issue is solved. Now proceeding to the linear algebra... :)

Author:  tam [ Tue May 11, 2010 4:12 pm ]
Post subject:  Re: QR decomposition, how to unpack Q in full

Test attachment. Wow, it worked. Awesome Sergey!!

Attachments:
testAttachment.txt [10 Bytes]
Downloaded 2478 times

Author:  Sergey.Bochkanov [ Tue May 11, 2010 6:37 pm ]
Post subject:  Re: QR decomposition, how to unpack Q in full

That's very strange - I've checked rmatrixqr() on your problem, and it worked.

I was unable to completely reproduce your code, because there were external variables A/Q/R. So I've made a constant array from the matrix you've provided. My code looks as follows:
Code:
#include "ortfac.h"
#include "_data.h"
int main(int argc, char* argv[])
{
    int myRow = 5849;
    int myCol = 16;
    int m = myRow;
    int n = myCol;
    int i, j;
    ap::real_2d_array a;
    ap::real_1d_array tau;
    ap::real_2d_array q;
    ap::real_2d_array r;
    a.setcontent(0, myRow-1, 0, myCol-1, _data);
    rmatrixqr(a, m, n, tau);
    rmatrixqrunpackq(a, m, n, tau, m, q);
   
    // testing: print several elements of q past the 16th column
    for(i=0; i<5; i++)
    {
        for(j=16; j<20; j++)
            printf("q(%ld,%ld) = %7.4lf\n", (long)i, (long)j, (double)q(i,j));
    }
    printf("q(%ld,%ld) = %7.4lf\n", (long)(myRow-1), (long)(myRow-1), (double)q(myRow-1,myRow-1));
    getch();
    return 0;
}

As long as I can see, it is essentially the same code that you've executed. But it works OK.

I'd be glad to investigate it further, if you can answer several questions:
1. Could you execute my code (see attachment) and tell me its output?
2. How did you decide that Q were not "fully" unpacked? Using debugger to examine q (lower case)? Or through examining external variable Q (upper case)?
3. What compiler/OS/optimization settings you used for your project?
4. Did you modify ALGLIB files, or you running unmodified library? What version?

Attachments:
my-results.zip [861.59 KiB]
Downloaded 2649 times

Page 1 of 2 All times are UTC
Powered by phpBB © 2000, 2002, 2005, 2007 phpBB Group
http://www.phpbb.com/