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_stgexc_work( int matrix_order, lapack_logical wantq,
00038 lapack_logical wantz, lapack_int n, float* a,
00039 lapack_int lda, float* b, lapack_int ldb,
00040 float* q, lapack_int ldq, float* z,
00041 lapack_int ldz, lapack_int* ifst,
00042 lapack_int* ilst, float* work,
00043 lapack_int lwork )
00044 {
00045 lapack_int info = 0;
00046 if( matrix_order == LAPACK_COL_MAJOR ) {
00047
00048 LAPACK_stgexc( &wantq, &wantz, &n, a, &lda, b, &ldb, q, &ldq, z, &ldz,
00049 ifst, ilst, work, &lwork, &info );
00050 if( info < 0 ) {
00051 info = info - 1;
00052 }
00053 } else if( matrix_order == LAPACK_ROW_MAJOR ) {
00054 lapack_int lda_t = MAX(1,n);
00055 lapack_int ldb_t = MAX(1,n);
00056 lapack_int ldq_t = MAX(1,n);
00057 lapack_int ldz_t = MAX(1,n);
00058 float* a_t = NULL;
00059 float* b_t = NULL;
00060 float* q_t = NULL;
00061 float* z_t = NULL;
00062
00063 if( lda < n ) {
00064 info = -6;
00065 LAPACKE_xerbla( "LAPACKE_stgexc_work", info );
00066 return info;
00067 }
00068 if( ldb < n ) {
00069 info = -8;
00070 LAPACKE_xerbla( "LAPACKE_stgexc_work", info );
00071 return info;
00072 }
00073 if( ldq < n ) {
00074 info = -10;
00075 LAPACKE_xerbla( "LAPACKE_stgexc_work", info );
00076 return info;
00077 }
00078 if( ldz < n ) {
00079 info = -12;
00080 LAPACKE_xerbla( "LAPACKE_stgexc_work", info );
00081 return info;
00082 }
00083
00084 if( lwork == -1 ) {
00085 LAPACK_stgexc( &wantq, &wantz, &n, a, &lda_t, b, &ldb_t, q, &ldq_t,
00086 z, &ldz_t, ifst, ilst, work, &lwork, &info );
00087 return (info < 0) ? (info - 1) : info;
00088 }
00089
00090 a_t = (float*)LAPACKE_malloc( sizeof(float) * lda_t * MAX(1,n) );
00091 if( a_t == NULL ) {
00092 info = LAPACK_TRANSPOSE_MEMORY_ERROR;
00093 goto exit_level_0;
00094 }
00095 b_t = (float*)LAPACKE_malloc( sizeof(float) * ldb_t * MAX(1,n) );
00096 if( b_t == NULL ) {
00097 info = LAPACK_TRANSPOSE_MEMORY_ERROR;
00098 goto exit_level_1;
00099 }
00100 if( wantq ) {
00101 q_t = (float*)LAPACKE_malloc( sizeof(float) * ldq_t * MAX(1,n) );
00102 if( q_t == NULL ) {
00103 info = LAPACK_TRANSPOSE_MEMORY_ERROR;
00104 goto exit_level_2;
00105 }
00106 }
00107 if( wantz ) {
00108 z_t = (float*)LAPACKE_malloc( sizeof(float) * ldz_t * MAX(1,n) );
00109 if( z_t == NULL ) {
00110 info = LAPACK_TRANSPOSE_MEMORY_ERROR;
00111 goto exit_level_3;
00112 }
00113 }
00114
00115 LAPACKE_sge_trans( matrix_order, n, n, a, lda, a_t, lda_t );
00116 LAPACKE_sge_trans( matrix_order, n, n, b, ldb, b_t, ldb_t );
00117 if( wantq ) {
00118 LAPACKE_sge_trans( matrix_order, n, n, q, ldq, q_t, ldq_t );
00119 }
00120 if( wantz ) {
00121 LAPACKE_sge_trans( matrix_order, n, n, z, ldz, z_t, ldz_t );
00122 }
00123
00124 LAPACK_stgexc( &wantq, &wantz, &n, a_t, &lda_t, b_t, &ldb_t, q_t,
00125 &ldq_t, z_t, &ldz_t, ifst, ilst, work, &lwork, &info );
00126 if( info < 0 ) {
00127 info = info - 1;
00128 }
00129
00130 LAPACKE_sge_trans( LAPACK_COL_MAJOR, n, n, a_t, lda_t, a, lda );
00131 LAPACKE_sge_trans( LAPACK_COL_MAJOR, n, n, b_t, ldb_t, b, ldb );
00132 if( wantq ) {
00133 LAPACKE_sge_trans( LAPACK_COL_MAJOR, n, n, q_t, ldq_t, q, ldq );
00134 }
00135 if( wantz ) {
00136 LAPACKE_sge_trans( LAPACK_COL_MAJOR, n, n, z_t, ldz_t, z, ldz );
00137 }
00138
00139 if( wantz ) {
00140 LAPACKE_free( z_t );
00141 }
00142 exit_level_3:
00143 if( wantq ) {
00144 LAPACKE_free( q_t );
00145 }
00146 exit_level_2:
00147 LAPACKE_free( b_t );
00148 exit_level_1:
00149 LAPACKE_free( a_t );
00150 exit_level_0:
00151 if( info == LAPACK_TRANSPOSE_MEMORY_ERROR ) {
00152 LAPACKE_xerbla( "LAPACKE_stgexc_work", info );
00153 }
00154 } else {
00155 info = -1;
00156 LAPACKE_xerbla( "LAPACKE_stgexc_work", info );
00157 }
00158 return info;
00159 }