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
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044
00045
00046
00047
00048
00049 #include <stdio.h>
00050 #include "lapacke.h"
00051 #include "lapacke_utils.h"
00052 #include "test_utils.h"
00053
00054 static void init_scalars_shseqr( char *job, char *compz, lapack_int *n,
00055 lapack_int *ilo, lapack_int *ihi,
00056 lapack_int *ldh, lapack_int *ldz,
00057 lapack_int *lwork );
00058 static void init_h( lapack_int size, float *h );
00059 static void init_wr( lapack_int size, float *wr );
00060 static void init_wi( lapack_int size, float *wi );
00061 static void init_z( lapack_int size, float *z );
00062 static void init_work( lapack_int size, float *work );
00063 static int compare_shseqr( float *h, float *h_i, float *wr, float *wr_i,
00064 float *wi, float *wi_i, float *z, float *z_i,
00065 lapack_int info, lapack_int info_i, char compz,
00066 lapack_int ldh, lapack_int ldz, lapack_int n );
00067
00068 int main(void)
00069 {
00070
00071 char job, job_i;
00072 char compz, compz_i;
00073 lapack_int n, n_i;
00074 lapack_int ilo, ilo_i;
00075 lapack_int ihi, ihi_i;
00076 lapack_int ldh, ldh_i;
00077 lapack_int ldh_r;
00078 lapack_int ldz, ldz_i;
00079 lapack_int ldz_r;
00080 lapack_int lwork, lwork_i;
00081 lapack_int info, info_i;
00082 lapack_int i;
00083 int failed;
00084
00085
00086 float *h = NULL, *h_i = NULL;
00087 float *wr = NULL, *wr_i = NULL;
00088 float *wi = NULL, *wi_i = NULL;
00089 float *z = NULL, *z_i = NULL;
00090 float *work = NULL, *work_i = NULL;
00091 float *h_save = NULL;
00092 float *wr_save = NULL;
00093 float *wi_save = NULL;
00094 float *z_save = NULL;
00095 float *h_r = NULL;
00096 float *z_r = NULL;
00097
00098
00099 init_scalars_shseqr( &job, &compz, &n, &ilo, &ihi, &ldh, &ldz, &lwork );
00100 ldh_r = n+2;
00101 ldz_r = n+2;
00102 job_i = job;
00103 compz_i = compz;
00104 n_i = n;
00105 ilo_i = ilo;
00106 ihi_i = ihi;
00107 ldh_i = ldh;
00108 ldz_i = ldz;
00109 lwork_i = lwork;
00110
00111
00112 h = (float *)LAPACKE_malloc( ldh*n * sizeof(float) );
00113 wr = (float *)LAPACKE_malloc( n * sizeof(float) );
00114 wi = (float *)LAPACKE_malloc( n * sizeof(float) );
00115 z = (float *)LAPACKE_malloc( ldz*n * sizeof(float) );
00116 work = (float *)LAPACKE_malloc( lwork * sizeof(float) );
00117
00118
00119 h_i = (float *)LAPACKE_malloc( ldh*n * sizeof(float) );
00120 wr_i = (float *)LAPACKE_malloc( n * sizeof(float) );
00121 wi_i = (float *)LAPACKE_malloc( n * sizeof(float) );
00122 z_i = (float *)LAPACKE_malloc( ldz*n * sizeof(float) );
00123 work_i = (float *)LAPACKE_malloc( lwork * sizeof(float) );
00124
00125
00126 h_save = (float *)LAPACKE_malloc( ldh*n * sizeof(float) );
00127 wr_save = (float *)LAPACKE_malloc( n * sizeof(float) );
00128 wi_save = (float *)LAPACKE_malloc( n * sizeof(float) );
00129 z_save = (float *)LAPACKE_malloc( ldz*n * sizeof(float) );
00130
00131
00132 h_r = (float *)LAPACKE_malloc( n*(n+2) * sizeof(float) );
00133 z_r = (float *)LAPACKE_malloc( n*(n+2) * sizeof(float) );
00134
00135
00136 init_h( ldh*n, h );
00137 init_wr( n, wr );
00138 init_wi( n, wi );
00139 init_z( ldz*n, z );
00140 init_work( lwork, work );
00141
00142
00143 for( i = 0; i < ldh*n; i++ ) {
00144 h_save[i] = h[i];
00145 }
00146 for( i = 0; i < n; i++ ) {
00147 wr_save[i] = wr[i];
00148 }
00149 for( i = 0; i < n; i++ ) {
00150 wi_save[i] = wi[i];
00151 }
00152 for( i = 0; i < ldz*n; i++ ) {
00153 z_save[i] = z[i];
00154 }
00155
00156
00157 shseqr_( &job, &compz, &n, &ilo, &ihi, h, &ldh, wr, wi, z, &ldz, work,
00158 &lwork, &info );
00159
00160
00161
00162 for( i = 0; i < ldh*n; i++ ) {
00163 h_i[i] = h_save[i];
00164 }
00165 for( i = 0; i < n; i++ ) {
00166 wr_i[i] = wr_save[i];
00167 }
00168 for( i = 0; i < n; i++ ) {
00169 wi_i[i] = wi_save[i];
00170 }
00171 for( i = 0; i < ldz*n; i++ ) {
00172 z_i[i] = z_save[i];
00173 }
00174 for( i = 0; i < lwork; i++ ) {
00175 work_i[i] = work[i];
00176 }
00177 info_i = LAPACKE_shseqr_work( LAPACK_COL_MAJOR, job_i, compz_i, n_i, ilo_i,
00178 ihi_i, h_i, ldh_i, wr_i, wi_i, z_i, ldz_i,
00179 work_i, lwork_i );
00180
00181 failed = compare_shseqr( h, h_i, wr, wr_i, wi, wi_i, z, z_i, info, info_i,
00182 compz, ldh, ldz, n );
00183 if( failed == 0 ) {
00184 printf( "PASSED: column-major middle-level interface to shseqr\n" );
00185 } else {
00186 printf( "FAILED: column-major middle-level interface to shseqr\n" );
00187 }
00188
00189
00190
00191 for( i = 0; i < ldh*n; i++ ) {
00192 h_i[i] = h_save[i];
00193 }
00194 for( i = 0; i < n; i++ ) {
00195 wr_i[i] = wr_save[i];
00196 }
00197 for( i = 0; i < n; i++ ) {
00198 wi_i[i] = wi_save[i];
00199 }
00200 for( i = 0; i < ldz*n; i++ ) {
00201 z_i[i] = z_save[i];
00202 }
00203 for( i = 0; i < lwork; i++ ) {
00204 work_i[i] = work[i];
00205 }
00206 info_i = LAPACKE_shseqr( LAPACK_COL_MAJOR, job_i, compz_i, n_i, ilo_i,
00207 ihi_i, h_i, ldh_i, wr_i, wi_i, z_i, ldz_i );
00208
00209 failed = compare_shseqr( h, h_i, wr, wr_i, wi, wi_i, z, z_i, info, info_i,
00210 compz, ldh, ldz, n );
00211 if( failed == 0 ) {
00212 printf( "PASSED: column-major high-level interface to shseqr\n" );
00213 } else {
00214 printf( "FAILED: column-major high-level interface to shseqr\n" );
00215 }
00216
00217
00218
00219 for( i = 0; i < ldh*n; i++ ) {
00220 h_i[i] = h_save[i];
00221 }
00222 for( i = 0; i < n; i++ ) {
00223 wr_i[i] = wr_save[i];
00224 }
00225 for( i = 0; i < n; i++ ) {
00226 wi_i[i] = wi_save[i];
00227 }
00228 for( i = 0; i < ldz*n; i++ ) {
00229 z_i[i] = z_save[i];
00230 }
00231 for( i = 0; i < lwork; i++ ) {
00232 work_i[i] = work[i];
00233 }
00234
00235 LAPACKE_sge_trans( LAPACK_COL_MAJOR, n, n, h_i, ldh, h_r, n+2 );
00236 if( LAPACKE_lsame( compz, 'i' ) || LAPACKE_lsame( compz, 'v' ) ) {
00237 LAPACKE_sge_trans( LAPACK_COL_MAJOR, n, n, z_i, ldz, z_r, n+2 );
00238 }
00239 info_i = LAPACKE_shseqr_work( LAPACK_ROW_MAJOR, job_i, compz_i, n_i, ilo_i,
00240 ihi_i, h_r, ldh_r, wr_i, wi_i, z_r, ldz_r,
00241 work_i, lwork_i );
00242
00243 LAPACKE_sge_trans( LAPACK_ROW_MAJOR, n, n, h_r, n+2, h_i, ldh );
00244 if( LAPACKE_lsame( compz, 'i' ) || LAPACKE_lsame( compz, 'v' ) ) {
00245 LAPACKE_sge_trans( LAPACK_ROW_MAJOR, n, n, z_r, n+2, z_i, ldz );
00246 }
00247
00248 failed = compare_shseqr( h, h_i, wr, wr_i, wi, wi_i, z, z_i, info, info_i,
00249 compz, ldh, ldz, n );
00250 if( failed == 0 ) {
00251 printf( "PASSED: row-major middle-level interface to shseqr\n" );
00252 } else {
00253 printf( "FAILED: row-major middle-level interface to shseqr\n" );
00254 }
00255
00256
00257
00258 for( i = 0; i < ldh*n; i++ ) {
00259 h_i[i] = h_save[i];
00260 }
00261 for( i = 0; i < n; i++ ) {
00262 wr_i[i] = wr_save[i];
00263 }
00264 for( i = 0; i < n; i++ ) {
00265 wi_i[i] = wi_save[i];
00266 }
00267 for( i = 0; i < ldz*n; i++ ) {
00268 z_i[i] = z_save[i];
00269 }
00270 for( i = 0; i < lwork; i++ ) {
00271 work_i[i] = work[i];
00272 }
00273
00274
00275 LAPACKE_sge_trans( LAPACK_COL_MAJOR, n, n, h_i, ldh, h_r, n+2 );
00276 if( LAPACKE_lsame( compz, 'i' ) || LAPACKE_lsame( compz, 'v' ) ) {
00277 LAPACKE_sge_trans( LAPACK_COL_MAJOR, n, n, z_i, ldz, z_r, n+2 );
00278 }
00279 info_i = LAPACKE_shseqr( LAPACK_ROW_MAJOR, job_i, compz_i, n_i, ilo_i,
00280 ihi_i, h_r, ldh_r, wr_i, wi_i, z_r, ldz_r );
00281
00282 LAPACKE_sge_trans( LAPACK_ROW_MAJOR, n, n, h_r, n+2, h_i, ldh );
00283 if( LAPACKE_lsame( compz, 'i' ) || LAPACKE_lsame( compz, 'v' ) ) {
00284 LAPACKE_sge_trans( LAPACK_ROW_MAJOR, n, n, z_r, n+2, z_i, ldz );
00285 }
00286
00287 failed = compare_shseqr( h, h_i, wr, wr_i, wi, wi_i, z, z_i, info, info_i,
00288 compz, ldh, ldz, n );
00289 if( failed == 0 ) {
00290 printf( "PASSED: row-major high-level interface to shseqr\n" );
00291 } else {
00292 printf( "FAILED: row-major high-level interface to shseqr\n" );
00293 }
00294
00295
00296 if( h != NULL ) {
00297 LAPACKE_free( h );
00298 }
00299 if( h_i != NULL ) {
00300 LAPACKE_free( h_i );
00301 }
00302 if( h_r != NULL ) {
00303 LAPACKE_free( h_r );
00304 }
00305 if( h_save != NULL ) {
00306 LAPACKE_free( h_save );
00307 }
00308 if( wr != NULL ) {
00309 LAPACKE_free( wr );
00310 }
00311 if( wr_i != NULL ) {
00312 LAPACKE_free( wr_i );
00313 }
00314 if( wr_save != NULL ) {
00315 LAPACKE_free( wr_save );
00316 }
00317 if( wi != NULL ) {
00318 LAPACKE_free( wi );
00319 }
00320 if( wi_i != NULL ) {
00321 LAPACKE_free( wi_i );
00322 }
00323 if( wi_save != NULL ) {
00324 LAPACKE_free( wi_save );
00325 }
00326 if( z != NULL ) {
00327 LAPACKE_free( z );
00328 }
00329 if( z_i != NULL ) {
00330 LAPACKE_free( z_i );
00331 }
00332 if( z_r != NULL ) {
00333 LAPACKE_free( z_r );
00334 }
00335 if( z_save != NULL ) {
00336 LAPACKE_free( z_save );
00337 }
00338 if( work != NULL ) {
00339 LAPACKE_free( work );
00340 }
00341 if( work_i != NULL ) {
00342 LAPACKE_free( work_i );
00343 }
00344
00345 return 0;
00346 }
00347
00348
00349 static void init_scalars_shseqr( char *job, char *compz, lapack_int *n,
00350 lapack_int *ilo, lapack_int *ihi,
00351 lapack_int *ldh, lapack_int *ldz,
00352 lapack_int *lwork )
00353 {
00354 *job = 'S';
00355 *compz = 'V';
00356 *n = 4;
00357 *ilo = 1;
00358 *ihi = 4;
00359 *ldh = 8;
00360 *ldz = 8;
00361 *lwork = 448;
00362
00363 return;
00364 }
00365
00366
00367 static void init_h( lapack_int size, float *h ) {
00368 lapack_int i;
00369 for( i = 0; i < size; i++ ) {
00370 h[i] = 0;
00371 }
00372 h[0] = 3.499999940e-001;
00373 h[8] = -1.159524173e-001;
00374 h[16] = -3.886010647e-001;
00375 h[24] = -2.941840589e-001;
00376 h[1] = -5.140038729e-001;
00377 h[9] = 1.224867105e-001;
00378 h[17] = 1.003598347e-001;
00379 h[25] = 1.125619039e-001;
00380 h[2] = -7.284721732e-001;
00381 h[10] = 6.442636251e-001;
00382 h[18] = -1.357001662e-001;
00383 h[26] = -9.768159688e-002;
00384 h[3] = 4.139046371e-001;
00385 h[11] = -1.665445715e-001;
00386 h[19] = 4.262443781e-001;
00387 h[27] = 1.632134318e-001;
00388 }
00389 static void init_wr( lapack_int size, float *wr ) {
00390 lapack_int i;
00391 for( i = 0; i < size; i++ ) {
00392 wr[i] = 0;
00393 }
00394 }
00395 static void init_wi( lapack_int size, float *wi ) {
00396 lapack_int i;
00397 for( i = 0; i < size; i++ ) {
00398 wi[i] = 0;
00399 }
00400 }
00401 static void init_z( lapack_int size, float *z ) {
00402 lapack_int i;
00403 for( i = 0; i < size; i++ ) {
00404 z[i] = 0;
00405 }
00406 z[0] = 1.000000000e+000;
00407 z[8] = 0.000000000e+000;
00408 z[16] = 0.000000000e+000;
00409 z[24] = 0.000000000e+000;
00410 z[1] = 0.000000000e+000;
00411 z[9] = -1.750959158e-001;
00412 z[17] = -9.674538374e-001;
00413 z[25] = -1.826872081e-001;
00414 z[2] = 0.000000000e+000;
00415 z[10] = 8.560246825e-001;
00416 z[18] = -2.412597835e-001;
00417 z[26] = 4.571821094e-001;
00418 z[3] = 0.000000000e+000;
00419 z[11] = -4.863776565e-001;
00420 z[19] = -7.633405924e-002;
00421 z[27] = 8.704078794e-001;
00422 }
00423 static void init_work( lapack_int size, float *work ) {
00424 lapack_int i;
00425 for( i = 0; i < size; i++ ) {
00426 work[i] = 0;
00427 }
00428 }
00429
00430
00431
00432 static int compare_shseqr( float *h, float *h_i, float *wr, float *wr_i,
00433 float *wi, float *wi_i, float *z, float *z_i,
00434 lapack_int info, lapack_int info_i, char compz,
00435 lapack_int ldh, lapack_int ldz, lapack_int n )
00436 {
00437 lapack_int i;
00438 int failed = 0;
00439 for( i = 0; i < ldh*n; i++ ) {
00440 failed += compare_floats(h[i],h_i[i]);
00441 }
00442 for( i = 0; i < n; i++ ) {
00443 failed += compare_floats(wr[i],wr_i[i]);
00444 }
00445 for( i = 0; i < n; i++ ) {
00446 failed += compare_floats(wi[i],wi_i[i]);
00447 }
00448 if( LAPACKE_lsame( compz, 'i' ) || LAPACKE_lsame( compz, 'v' ) ) {
00449 for( i = 0; i < ldz*n; i++ ) {
00450 failed += compare_floats(z[i],z_i[i]);
00451 }
00452 }
00453 failed += (info == info_i) ? 0 : 1;
00454 if( info != 0 || info_i != 0 ) {
00455 printf( "info=%d, info_i=%d\n",(int)info,(int)info_i );
00456 }
00457
00458 return failed;
00459 }