00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034 #include "lapacke.h"
00035 #include "lapacke_utils.h"
00036
00037 lapack_int LAPACKE_dggsvd_work( int matrix_order, char jobu, char jobv,
00038 char jobq, lapack_int m, lapack_int n,
00039 lapack_int p, lapack_int* k, lapack_int* l,
00040 double* a, lapack_int lda, double* b,
00041 lapack_int ldb, double* alpha, double* beta,
00042 double* u, lapack_int ldu, double* v,
00043 lapack_int ldv, double* q, lapack_int ldq,
00044 double* work, lapack_int* iwork )
00045 {
00046 lapack_int info = 0;
00047 if( matrix_order == LAPACK_COL_MAJOR ) {
00048
00049 LAPACK_dggsvd( &jobu, &jobv, &jobq, &m, &n, &p, k, l, a, &lda, b, &ldb,
00050 alpha, beta, u, &ldu, v, &ldv, q, &ldq, work, iwork,
00051 &info );
00052 if( info < 0 ) {
00053 info = info - 1;
00054 }
00055 } else if( matrix_order == LAPACK_ROW_MAJOR ) {
00056 lapack_int lda_t = MAX(1,m);
00057 lapack_int ldb_t = MAX(1,p);
00058 lapack_int ldq_t = MAX(1,n);
00059 lapack_int ldu_t = MAX(1,m);
00060 lapack_int ldv_t = MAX(1,p);
00061 double* a_t = NULL;
00062 double* b_t = NULL;
00063 double* u_t = NULL;
00064 double* v_t = NULL;
00065 double* q_t = NULL;
00066
00067 if( lda < n ) {
00068 info = -11;
00069 LAPACKE_xerbla( "LAPACKE_dggsvd_work", info );
00070 return info;
00071 }
00072 if( ldb < n ) {
00073 info = -13;
00074 LAPACKE_xerbla( "LAPACKE_dggsvd_work", info );
00075 return info;
00076 }
00077 if( ldq < n ) {
00078 info = -21;
00079 LAPACKE_xerbla( "LAPACKE_dggsvd_work", info );
00080 return info;
00081 }
00082 if( ldu < m ) {
00083 info = -17;
00084 LAPACKE_xerbla( "LAPACKE_dggsvd_work", info );
00085 return info;
00086 }
00087 if( ldv < p ) {
00088 info = -19;
00089 LAPACKE_xerbla( "LAPACKE_dggsvd_work", info );
00090 return info;
00091 }
00092
00093 a_t = (double*)LAPACKE_malloc( sizeof(double) * lda_t * MAX(1,n) );
00094 if( a_t == NULL ) {
00095 info = LAPACK_TRANSPOSE_MEMORY_ERROR;
00096 goto exit_level_0;
00097 }
00098 b_t = (double*)LAPACKE_malloc( sizeof(double) * ldb_t * MAX(1,n) );
00099 if( b_t == NULL ) {
00100 info = LAPACK_TRANSPOSE_MEMORY_ERROR;
00101 goto exit_level_1;
00102 }
00103 if( LAPACKE_lsame( jobu, 'u' ) ) {
00104 u_t = (double*)LAPACKE_malloc( sizeof(double) * ldu_t * MAX(1,m) );
00105 if( u_t == NULL ) {
00106 info = LAPACK_TRANSPOSE_MEMORY_ERROR;
00107 goto exit_level_2;
00108 }
00109 }
00110 if( LAPACKE_lsame( jobv, 'v' ) ) {
00111 v_t = (double*)LAPACKE_malloc( sizeof(double) * ldv_t * MAX(1,p) );
00112 if( v_t == NULL ) {
00113 info = LAPACK_TRANSPOSE_MEMORY_ERROR;
00114 goto exit_level_3;
00115 }
00116 }
00117 if( LAPACKE_lsame( jobq, 'q' ) ) {
00118 q_t = (double*)LAPACKE_malloc( sizeof(double) * ldq_t * MAX(1,n) );
00119 if( q_t == NULL ) {
00120 info = LAPACK_TRANSPOSE_MEMORY_ERROR;
00121 goto exit_level_4;
00122 }
00123 }
00124
00125 LAPACKE_dge_trans( matrix_order, m, n, a, lda, a_t, lda_t );
00126 LAPACKE_dge_trans( matrix_order, p, n, b, ldb, b_t, ldb_t );
00127
00128 LAPACK_dggsvd( &jobu, &jobv, &jobq, &m, &n, &p, k, l, a_t, &lda_t, b_t,
00129 &ldb_t, alpha, beta, u_t, &ldu_t, v_t, &ldv_t, q_t,
00130 &ldq_t, work, iwork, &info );
00131 if( info < 0 ) {
00132 info = info - 1;
00133 }
00134
00135 LAPACKE_dge_trans( LAPACK_COL_MAJOR, m, n, a_t, lda_t, a, lda );
00136 LAPACKE_dge_trans( LAPACK_COL_MAJOR, p, n, b_t, ldb_t, b, ldb );
00137 if( LAPACKE_lsame( jobu, 'u' ) ) {
00138 LAPACKE_dge_trans( LAPACK_COL_MAJOR, m, m, u_t, ldu_t, u, ldu );
00139 }
00140 if( LAPACKE_lsame( jobv, 'v' ) ) {
00141 LAPACKE_dge_trans( LAPACK_COL_MAJOR, p, p, v_t, ldv_t, v, ldv );
00142 }
00143 if( LAPACKE_lsame( jobq, 'q' ) ) {
00144 LAPACKE_dge_trans( LAPACK_COL_MAJOR, n, n, q_t, ldq_t, q, ldq );
00145 }
00146
00147 if( LAPACKE_lsame( jobq, 'q' ) ) {
00148 LAPACKE_free( q_t );
00149 }
00150 exit_level_4:
00151 if( LAPACKE_lsame( jobv, 'v' ) ) {
00152 LAPACKE_free( v_t );
00153 }
00154 exit_level_3:
00155 if( LAPACKE_lsame( jobu, 'u' ) ) {
00156 LAPACKE_free( u_t );
00157 }
00158 exit_level_2:
00159 LAPACKE_free( b_t );
00160 exit_level_1:
00161 LAPACKE_free( a_t );
00162 exit_level_0:
00163 if( info == LAPACK_TRANSPOSE_MEMORY_ERROR ) {
00164 LAPACKE_xerbla( "LAPACKE_dggsvd_work", info );
00165 }
00166 } else {
00167 info = -1;
00168 LAPACKE_xerbla( "LAPACKE_dggsvd_work", info );
00169 }
00170 return info;
00171 }