topic Hi, in Intel® oneAPI Math Kernel Library

QR decomposition for tall matrices

hemantp_p_ — Thu, 22 May 2014 20:32:57 GMT

Sir,

i use the routines *qwgrf and the *orgqr to generate the Q and R .

However if the matrix is tall ie more rows than coloumns these routines just dont work (my mistake ) also while using *orgqr

it requires the mtrix to passed as an arguement. given that for tall matrices Q will be square which matrix should be passed as an arguement.

The routine works just fine for square and wide but for tall i am not able to fiqure it out.

in the code attached the matrix size is 36 elemnts, 6*6 square, 9*4 tall, 4*9 wide.

arguents passed to the executable can be either of the three eg. ./a.out 9 4

There is no Lapack-E routine

mecej4 — Thu, 22 May 2014 22:21:47 GMT

There is no Lapack-E routine containing "qwgrf" in its name, and I am unable to decipher your second sentence.

Here is an example program that performs the QR decomposition of the 6 X 4 matrix in your example:

#include <stdio.h>
#include <stdlib.h>
#include <mkl.h>

int main(){
float A[]={-0.57, -1.28, -0.39,  0.25,
           -1.93,  1.08, -0.31, -2.14,
            2.30,  0.24,  0.40, -0.35,
           -1.93,  0.64, -0.66,  0.08,
            0.15,  0.30,  0.15, -2.13,
           -0.02,  1.03, -1.43,  0.50};
float *tau;
int iret,m,n,stride,i,j;

m=6; n=4; stride=n; tau=(float *)malloc(min(m,n)*sizeof(float));

iret=LAPACKE_sgeqrf(LAPACK_ROW_MAJOR,m,n,A,stride,tau);
printf("SGEQRF return code = %d\n",iret);
iret=LAPACKE_sorgqr(LAPACK_ROW_MAJOR,m,n,n,A,stride,tau);
printf("SORGQR return code = %d\n",iret);
for(i=0; i<m; i++){
   printf("\n"); for(j=0; j<n; j++)printf(" %8.4f",A[i*stride+j]);
   }
}

The results are in agreement with those given by Matlab.

Hi,

Ying_H_Intel — Fri, 23 May 2014 03:26:01 GMT

Hi,

As mecej4 point out, the stride in your code is wrong. It is stride of row (row-major). So it is N in your cases.

And the sgeqrf and sorgqr can work square, tall and wide.

for example, for M=9, N=4, the leading 4 column of Q are same as Matlab.

if (M>N)

K=N;

else
K=M;

tau = (float *)malloc(K*sizeof(float));

printf("yo2\n");

LAPACKE_sgeqrf( LAPACK_ROW_MAJOR, M, N, q_mat, N, tau );

LAPACKE_sorgqr( LAPACK_ROW_MAJOR, M, K, K, q_mat, N, tau );

printf("\nQ\n");
for(i=0; i<M; i++){
printf("\n"); for(j=0; j<N; j++)printf(" %8.4f",q_mat[i*N+j]);
}

If you'd like to get the whole of Q when M=9 and N=4.

Then you may need to a temp space A (M, M).

A=(float *)malloc(M*M*sizeof(float));

LAPACKE_sgeqrf( LAPACK_ROW_MAJOR, M, N, q_mat, N, tau );

for(i=0; i<M; i++){ for(j=0; j<N; j++) A[i*M+j]=q_mat[i*N+j]; }

LAPACKE_sorgqr( LAPACK_ROW_MAJOR, M, M, K, A, M, tau );

for(i=0; i<M; i++){ printf("\n"); for(j=0; j<M; j++) printf(" %6.3f",A[i*M+j]);
}

Best Regards,
Ying

yo2

6.8000 -6.0500 -0.4500 8.3200
-9.6700 -2.1100 -3.3000 2.5800
2.7100 -5.1400 5.6600 5.3600
-2.7000 4.3500 -7.2600 5.9700
-4.4400 0.2700 -7.1700 6.0800
8.2300 1.0800 9.0400 2.1400
-6.8700 5.6600 5.3600 -2.7000
4.3500 -7.2600 5.9700 -4.4400
0.2700 -7.1700 6.0800 4.3500

-0.3875 -0.2654 0.3130 -0.3702
0.5510 -0.4131 -0.0310 -0.1207
-0.1544 -0.3071 -0.2144 -0.3834
0.1538 0.2491 0.3303 -0.4066
0.2530 -0.0983 0.3349 -0.3051
-0.4689 0.2987 -0.3680 -0.3813
0.3914 0.2347 -0.6283 -0.1845
-0.2479 -0.4199 -0.1563 0.4326
-0.0154 -0.5218 -0.2815 -0.2765

Case 2:

-0.387 -0.265 0.313 -0.370 -0.139 -0.105 0.679 -0.154 -0.171
0.551 -0.413 -0.031 -0.121 -0.308 0.509 -0.001 -0.095 -0.384
-0.154 -0.307 -0.214 -0.383 0.115 -0.470 -0.488 -0.017 -0.464
0.154 0.249 0.330 -0.407 -0.426 -0.114 -0.125 0.650 0.085
0.253 -0.098 0.335 -0.305 0.810 0.168 0.034 0.175 0.068
-0.469 0.299 -0.368 -0.381 0.055 0.621 -0.084 0.085 -0.086
0.391 0.235 -0.628 -0.185 0.106 -0.275 0.499 0.148 -0.069
-0.248 -0.420 -0.156 0.433 0.094 0.082 0.150 0.696 -0.167
-0.015 -0.522 -0.281 -0.277 -0.106 0.007 -0.067 -0.005 0.746Press any key to
continue . . .

Sir,

hemantp_p_ — Fri, 23 May 2014 18:30:13 GMT

Sir,

Thank You all for your prompt help.