- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content

I have a matrix A=QR and the QR decomposition is computed using LAPACKE_dgeqrf. To be more concrete A and the output QR are included below. At this point I append a new column to A and it becomes A1, A1 is included below as well. At this point, I try to apply an update which requires computing first Q'*A1 and therefore I use the function LAPACKE_dormqr. This below is the documentation and the way I am invoking it for my given matrices:

matrix A of dimensions 4 x 2:

===================

7.500000e+01 4.000000e+00

4.000000e+01 2.800000e+01

6.600000e+01 5.000000e+00

1.800000e+01 1.000000e+01

matrix QR of dimensions 4 x 2:

===================

-1.091100e+02 -1.768857e+01

2.172614e-01 -2.474095e+01

3.584813e-01 -5.777660e-02

9.776762e-02 1.640589e-01

matrix A1 of dimensions 4 x 3:

===================

7.500000e+01 4.000000e+00 3.900000e+01

4.000000e+01 2.800000e+01 7.700000e+01

6.600000e+01 5.000000e+00 8.000000e+01

1.800000e+01 1.000000e+01 1.900000e+01

mklman.pdf MKL 10.3 Update 9:

===================

C:

lapack_int LAPACKE_ormqr( int matrix_order, char side, char trans, lapack_int m,

lapack_int n, lapack_int k, const

tau,

int l_order = LAPACK_COL_MAJOR;

lapack_int m = 4; // C->rows;

lapack_int n = 3; // C->cols;

lapack_int k = 3; // B->QR->rows;

lapack_int lda = m;

lapack_int ldc = m;

lapack_int info = 0;

lapack_int info = LAPACKE_dormqr(l_order, 'L', 'T', m, n, k, B->QR->data, lda, B->tau->data, C->data, ldc);

So far so good, but now the matrix coming out Q'*A1 is:

matrix Q'*A1 of dimensions 4 x 3

===================

-1.091100e+02 -1.768857e+01 -1.065621e+02

0.000000e+00 -2.474095e+01 -4.110871e+01

0.000000e+00 0.000000e+00 3.281543e+01

0.000000e+00 8.881784e-16 -9.419678e+00

Cross checking with matlab turns out that the elements in column Q'*A1(:,3) are incorrect ... what am I doing wrong here?

This is the same but done using matlab:

>> A = [7.500000e+01 4.000000e+00;

4.000000e+01 2.800000e+01;

6.600000e+01 5.000000e+00;

1.800000e+01 1.000000e+01]

A =

75 4

40 28

66 5

18 10

>> A1 = [A [39 77 80 19]']

A1 =

75 4 39

40 28 77

66 5 80

18 10 19

>> [Q R] = qr(A)

Q =

0.6874 -0.3298 -0.5647 0.3161

0.3666 0.8696 -0.2063 -0.2585

0.6049 -0.2304 0.6340 -0.4232

0.1650 0.2862 0.4865 0.8088

R =

109.1100 17.6886

0 24.7410

0 0

0 0

>> Q'*A1

ans =

109.1100 17.6886 106.5621

0.0000 24.7410 41.1087

0.0000 0.0000 22.0579

-0.0000 -0.0000 -26.0582

The final matlab answer and the final MKL LAPACK answer are different in the last column:

>> Matlab=Q'*A1

Matlab =

109.1100 17.6886 106.5621

0.0000 24.7410 41.1087

0.0000 0.0000 22.0579

-0.0000 -0.0000 -26.0582

>> MKL=[

-1.091100e+02 -1.768857e+01 -1.065621e+02;

0.000000e+00 -2.474095e+01 -4.110871e+01;

0.000000e+00 0.000000e+00 3.281543e+01;

0.000000e+00 8.881784e-16 -9.419678e+00]

MKL =

-109.1100 -17.6886 -106.5621

0 -24.7410 -41.1087

0 0 32.8154

0 0.0000 -9.4197

>> abs(Matlab) - abs(MKL)

ans =

0.0000 -0.0000 0.0000

0.0000 0.0000 -0.0000

0.0000 0.0000 -10.7575

0.0000 0.0000 16.6385

Why is it so different, what am I doing wrong when calling LAPACKE_dormqr?

Thanks in advance,

Best regards,

Giovanni

1 Solution

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content

Your parameter setting should be correct. But the different result between MKL and Matlab should be the QR implementation itself.

for example, the below showA=QR. Theelements in Greenaresignificant. But

Red values could take any value with only one condition Q*Q^T=E.

As a result, MKL and Matlab have different last twocolumns, then different Q^T *A1 in last two row.

Best Regards,

Ying

Link Copied

7 Replies

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content

Since Matlab itself uses MKL for internal operations, I would not expect to see any significant discrepancies.

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content

I know what you mean but that's why I compare taking abs (see the end of my OP). It is also possible that the Q come different by permuting or changing signs for some row elements for a given column. But AFAIK Q'*C should give the same result specially like as you say matlab is just invoking MKL under the hood.

I guess this might have something to do with the way LAPACKE_dormqr works? Note that I am always using the compact QR representation upper trapezoidal R with lower min(M,N) householder reflectors and tau coeffients vector.

Best regards,

Giovanni

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content

*> But AFAIK Q'*C should give the same result specially like as you say matlab > is just invoking MKL under the hood.*

Two points:

1) The version of MKL used by Matlab need not be the same as the one used by your Fortran/C compiler.

2) If Matlab gives you Q

_{1}and Fortran/C gives you Q

_{2}= -Q

_{1}, Q

_{1}'C will be equal to -Q

_{2}'C, since you are supplying the C after the decomposition was performed.

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content

Ok point taken :)

Could you confirm this is the right way to call the function?

QR is computed using DGEQRF over matrix B (m x (n - 1)).

C is a copy of matrix B with one extra column C (m x n).

C (m x n)

B->QR (m x (n - 1)) as returned by DGEQRF having MIN(m, n - 1) reflectors in the lower triangular.

int l_order = LAPACK_COL_MAJOR;

lapack_int m = C->rows;

lapack_int n = C->cols;

lapack_int k = MIN(B->QR->rows, B->QR->cols);

lapack_int lda = m;

lapack_int ldc = m;

lapack_int info = 0;

lapack_int info = LAPACKE_dormqr(l_order, 'L', 'T', m, n, k, B->QR->data, lda, B->tau->data, C->data, ldc);

Many thanks in advance,

Best regards,

Giovanni

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content

Your parameter setting should be correct. But the different result between MKL and Matlab should be the QR implementation itself.

for example, the below showA=QR. Theelements in Greenaresignificant. But

Red values could take any value with only one condition Q*Q^T=E.

As a result, MKL and Matlab have different last twocolumns, then different Q^T *A1 in last two row.

Best Regards,

Ying

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content

Thank you for following up on this :) I kinda of forgot this Thread already :) but is good to know. I got this bit working long ago.

Best!

Giovanni

Topic Options

- Subscribe to RSS Feed
- Mark Topic as New
- Mark Topic as Read
- Float this Topic for Current User
- Bookmark
- Subscribe
- Printer Friendly Page