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_dsbtrd( char *vect, char *uplo, lapack_int *n,
00055 lapack_int *kd, lapack_int *ldab,
00056 lapack_int *ldq );
00057 static void init_ab( lapack_int size, double *ab );
00058 static void init_d( lapack_int size, double *d );
00059 static void init_e( lapack_int size, double *e );
00060 static void init_q( lapack_int size, double *q );
00061 static void init_work( lapack_int size, double *work );
00062 static int compare_dsbtrd( double *ab, double *ab_i, double *d, double *d_i,
00063 double *e, double *e_i, double *q, double *q_i,
00064 lapack_int info, lapack_int info_i, lapack_int ldab,
00065 lapack_int ldq, lapack_int n, char vect );
00066
00067 int main(void)
00068 {
00069
00070 char vect, vect_i;
00071 char uplo, uplo_i;
00072 lapack_int n, n_i;
00073 lapack_int kd, kd_i;
00074 lapack_int ldab, ldab_i;
00075 lapack_int ldab_r;
00076 lapack_int ldq, ldq_i;
00077 lapack_int ldq_r;
00078 lapack_int info, info_i;
00079 lapack_int i;
00080 int failed;
00081
00082
00083 double *ab = NULL, *ab_i = NULL;
00084 double *d = NULL, *d_i = NULL;
00085 double *e = NULL, *e_i = NULL;
00086 double *q = NULL, *q_i = NULL;
00087 double *work = NULL, *work_i = NULL;
00088 double *ab_save = NULL;
00089 double *d_save = NULL;
00090 double *e_save = NULL;
00091 double *q_save = NULL;
00092 double *ab_r = NULL;
00093 double *q_r = NULL;
00094
00095
00096 init_scalars_dsbtrd( &vect, &uplo, &n, &kd, &ldab, &ldq );
00097 ldab_r = n+2;
00098 ldq_r = n+2;
00099 vect_i = vect;
00100 uplo_i = uplo;
00101 n_i = n;
00102 kd_i = kd;
00103 ldab_i = ldab;
00104 ldq_i = ldq;
00105
00106
00107 ab = (double *)LAPACKE_malloc( ldab*n * sizeof(double) );
00108 d = (double *)LAPACKE_malloc( n * sizeof(double) );
00109 e = (double *)LAPACKE_malloc( (n-1) * sizeof(double) );
00110 q = (double *)LAPACKE_malloc( ldq*n * sizeof(double) );
00111 work = (double *)LAPACKE_malloc( n * sizeof(double) );
00112
00113
00114 ab_i = (double *)LAPACKE_malloc( ldab*n * sizeof(double) );
00115 d_i = (double *)LAPACKE_malloc( n * sizeof(double) );
00116 e_i = (double *)LAPACKE_malloc( (n-1) * sizeof(double) );
00117 q_i = (double *)LAPACKE_malloc( ldq*n * sizeof(double) );
00118 work_i = (double *)LAPACKE_malloc( n * sizeof(double) );
00119
00120
00121 ab_save = (double *)LAPACKE_malloc( ldab*n * sizeof(double) );
00122 d_save = (double *)LAPACKE_malloc( n * sizeof(double) );
00123 e_save = (double *)LAPACKE_malloc( (n-1) * sizeof(double) );
00124 q_save = (double *)LAPACKE_malloc( ldq*n * sizeof(double) );
00125
00126
00127 ab_r = (double *)LAPACKE_malloc( (kd+1)*(n+2) * sizeof(double) );
00128 q_r = (double *)LAPACKE_malloc( n*(n+2) * sizeof(double) );
00129
00130
00131 init_ab( ldab*n, ab );
00132 init_d( n, d );
00133 init_e( (n-1), e );
00134 init_q( ldq*n, q );
00135 init_work( n, work );
00136
00137
00138 for( i = 0; i < ldab*n; i++ ) {
00139 ab_save[i] = ab[i];
00140 }
00141 for( i = 0; i < n; i++ ) {
00142 d_save[i] = d[i];
00143 }
00144 for( i = 0; i < (n-1); i++ ) {
00145 e_save[i] = e[i];
00146 }
00147 for( i = 0; i < ldq*n; i++ ) {
00148 q_save[i] = q[i];
00149 }
00150
00151
00152 dsbtrd_( &vect, &uplo, &n, &kd, ab, &ldab, d, e, q, &ldq, work, &info );
00153
00154
00155
00156 for( i = 0; i < ldab*n; i++ ) {
00157 ab_i[i] = ab_save[i];
00158 }
00159 for( i = 0; i < n; i++ ) {
00160 d_i[i] = d_save[i];
00161 }
00162 for( i = 0; i < (n-1); i++ ) {
00163 e_i[i] = e_save[i];
00164 }
00165 for( i = 0; i < ldq*n; i++ ) {
00166 q_i[i] = q_save[i];
00167 }
00168 for( i = 0; i < n; i++ ) {
00169 work_i[i] = work[i];
00170 }
00171 info_i = LAPACKE_dsbtrd_work( LAPACK_COL_MAJOR, vect_i, uplo_i, n_i, kd_i,
00172 ab_i, ldab_i, d_i, e_i, q_i, ldq_i, work_i );
00173
00174 failed = compare_dsbtrd( ab, ab_i, d, d_i, e, e_i, q, q_i, info, info_i,
00175 ldab, ldq, n, vect );
00176 if( failed == 0 ) {
00177 printf( "PASSED: column-major middle-level interface to dsbtrd\n" );
00178 } else {
00179 printf( "FAILED: column-major middle-level interface to dsbtrd\n" );
00180 }
00181
00182
00183
00184 for( i = 0; i < ldab*n; i++ ) {
00185 ab_i[i] = ab_save[i];
00186 }
00187 for( i = 0; i < n; i++ ) {
00188 d_i[i] = d_save[i];
00189 }
00190 for( i = 0; i < (n-1); i++ ) {
00191 e_i[i] = e_save[i];
00192 }
00193 for( i = 0; i < ldq*n; i++ ) {
00194 q_i[i] = q_save[i];
00195 }
00196 for( i = 0; i < n; i++ ) {
00197 work_i[i] = work[i];
00198 }
00199 info_i = LAPACKE_dsbtrd( LAPACK_COL_MAJOR, vect_i, uplo_i, n_i, kd_i, ab_i,
00200 ldab_i, d_i, e_i, q_i, ldq_i );
00201
00202 failed = compare_dsbtrd( ab, ab_i, d, d_i, e, e_i, q, q_i, info, info_i,
00203 ldab, ldq, n, vect );
00204 if( failed == 0 ) {
00205 printf( "PASSED: column-major high-level interface to dsbtrd\n" );
00206 } else {
00207 printf( "FAILED: column-major high-level interface to dsbtrd\n" );
00208 }
00209
00210
00211
00212 for( i = 0; i < ldab*n; i++ ) {
00213 ab_i[i] = ab_save[i];
00214 }
00215 for( i = 0; i < n; i++ ) {
00216 d_i[i] = d_save[i];
00217 }
00218 for( i = 0; i < (n-1); i++ ) {
00219 e_i[i] = e_save[i];
00220 }
00221 for( i = 0; i < ldq*n; i++ ) {
00222 q_i[i] = q_save[i];
00223 }
00224 for( i = 0; i < n; i++ ) {
00225 work_i[i] = work[i];
00226 }
00227
00228 LAPACKE_dge_trans( LAPACK_COL_MAJOR, kd+1, n, ab_i, ldab, ab_r, n+2 );
00229 if( LAPACKE_lsame( vect, 'u' ) || LAPACKE_lsame( vect, 'v' ) ) {
00230 LAPACKE_dge_trans( LAPACK_COL_MAJOR, n, n, q_i, ldq, q_r, n+2 );
00231 }
00232 info_i = LAPACKE_dsbtrd_work( LAPACK_ROW_MAJOR, vect_i, uplo_i, n_i, kd_i,
00233 ab_r, ldab_r, d_i, e_i, q_r, ldq_r, work_i );
00234
00235 LAPACKE_dge_trans( LAPACK_ROW_MAJOR, kd+1, n, ab_r, n+2, ab_i, ldab );
00236 if( LAPACKE_lsame( vect, 'u' ) || LAPACKE_lsame( vect, 'v' ) ) {
00237 LAPACKE_dge_trans( LAPACK_ROW_MAJOR, n, n, q_r, n+2, q_i, ldq );
00238 }
00239
00240 failed = compare_dsbtrd( ab, ab_i, d, d_i, e, e_i, q, q_i, info, info_i,
00241 ldab, ldq, n, vect );
00242 if( failed == 0 ) {
00243 printf( "PASSED: row-major middle-level interface to dsbtrd\n" );
00244 } else {
00245 printf( "FAILED: row-major middle-level interface to dsbtrd\n" );
00246 }
00247
00248
00249
00250 for( i = 0; i < ldab*n; i++ ) {
00251 ab_i[i] = ab_save[i];
00252 }
00253 for( i = 0; i < n; i++ ) {
00254 d_i[i] = d_save[i];
00255 }
00256 for( i = 0; i < (n-1); i++ ) {
00257 e_i[i] = e_save[i];
00258 }
00259 for( i = 0; i < ldq*n; i++ ) {
00260 q_i[i] = q_save[i];
00261 }
00262 for( i = 0; i < n; i++ ) {
00263 work_i[i] = work[i];
00264 }
00265
00266
00267 LAPACKE_dge_trans( LAPACK_COL_MAJOR, kd+1, n, ab_i, ldab, ab_r, n+2 );
00268 if( LAPACKE_lsame( vect, 'u' ) || LAPACKE_lsame( vect, 'v' ) ) {
00269 LAPACKE_dge_trans( LAPACK_COL_MAJOR, n, n, q_i, ldq, q_r, n+2 );
00270 }
00271 info_i = LAPACKE_dsbtrd( LAPACK_ROW_MAJOR, vect_i, uplo_i, n_i, kd_i, ab_r,
00272 ldab_r, d_i, e_i, q_r, ldq_r );
00273
00274 LAPACKE_dge_trans( LAPACK_ROW_MAJOR, kd+1, n, ab_r, n+2, ab_i, ldab );
00275 if( LAPACKE_lsame( vect, 'u' ) || LAPACKE_lsame( vect, 'v' ) ) {
00276 LAPACKE_dge_trans( LAPACK_ROW_MAJOR, n, n, q_r, n+2, q_i, ldq );
00277 }
00278
00279 failed = compare_dsbtrd( ab, ab_i, d, d_i, e, e_i, q, q_i, info, info_i,
00280 ldab, ldq, n, vect );
00281 if( failed == 0 ) {
00282 printf( "PASSED: row-major high-level interface to dsbtrd\n" );
00283 } else {
00284 printf( "FAILED: row-major high-level interface to dsbtrd\n" );
00285 }
00286
00287
00288 if( ab != NULL ) {
00289 LAPACKE_free( ab );
00290 }
00291 if( ab_i != NULL ) {
00292 LAPACKE_free( ab_i );
00293 }
00294 if( ab_r != NULL ) {
00295 LAPACKE_free( ab_r );
00296 }
00297 if( ab_save != NULL ) {
00298 LAPACKE_free( ab_save );
00299 }
00300 if( d != NULL ) {
00301 LAPACKE_free( d );
00302 }
00303 if( d_i != NULL ) {
00304 LAPACKE_free( d_i );
00305 }
00306 if( d_save != NULL ) {
00307 LAPACKE_free( d_save );
00308 }
00309 if( e != NULL ) {
00310 LAPACKE_free( e );
00311 }
00312 if( e_i != NULL ) {
00313 LAPACKE_free( e_i );
00314 }
00315 if( e_save != NULL ) {
00316 LAPACKE_free( e_save );
00317 }
00318 if( q != NULL ) {
00319 LAPACKE_free( q );
00320 }
00321 if( q_i != NULL ) {
00322 LAPACKE_free( q_i );
00323 }
00324 if( q_r != NULL ) {
00325 LAPACKE_free( q_r );
00326 }
00327 if( q_save != NULL ) {
00328 LAPACKE_free( q_save );
00329 }
00330 if( work != NULL ) {
00331 LAPACKE_free( work );
00332 }
00333 if( work_i != NULL ) {
00334 LAPACKE_free( work_i );
00335 }
00336
00337 return 0;
00338 }
00339
00340
00341 static void init_scalars_dsbtrd( char *vect, char *uplo, lapack_int *n,
00342 lapack_int *kd, lapack_int *ldab,
00343 lapack_int *ldq )
00344 {
00345 *vect = 'V';
00346 *uplo = 'L';
00347 *n = 4;
00348 *kd = 2;
00349 *ldab = 9;
00350 *ldq = 8;
00351
00352 return;
00353 }
00354
00355
00356 static void init_ab( lapack_int size, double *ab ) {
00357 lapack_int i;
00358 for( i = 0; i < size; i++ ) {
00359 ab[i] = 0;
00360 }
00361 ab[0] = 4.99000000000000020e+000;
00362 ab[9] = 1.05000000000000000e+000;
00363 ab[18] = -2.31000000000000010e+000;
00364 ab[27] = -4.29999999999999990e-001;
00365 ab[1] = 4.00000000000000010e-002;
00366 ab[10] = -7.90000000000000040e-001;
00367 ab[19] = -1.30000000000000000e+000;
00368 ab[28] = 0.00000000000000000e+000;
00369 ab[2] = 2.20000000000000000e-001;
00370 ab[11] = 1.04000000000000000e+000;
00371 ab[20] = 0.00000000000000000e+000;
00372 ab[29] = 0.00000000000000000e+000;
00373 }
00374 static void init_d( lapack_int size, double *d ) {
00375 lapack_int i;
00376 for( i = 0; i < size; i++ ) {
00377 d[i] = 0;
00378 }
00379 }
00380 static void init_e( lapack_int size, double *e ) {
00381 lapack_int i;
00382 for( i = 0; i < size; i++ ) {
00383 e[i] = 0;
00384 }
00385 }
00386 static void init_q( lapack_int size, double *q ) {
00387 lapack_int i;
00388 for( i = 0; i < size; i++ ) {
00389 q[i] = 0;
00390 }
00391 q[0] = 0.00000000000000000e+000;
00392 q[8] = 0.00000000000000000e+000;
00393 q[16] = 0.00000000000000000e+000;
00394 q[24] = 0.00000000000000000e+000;
00395 q[1] = 0.00000000000000000e+000;
00396 q[9] = 0.00000000000000000e+000;
00397 q[17] = 0.00000000000000000e+000;
00398 q[25] = 0.00000000000000000e+000;
00399 q[2] = 0.00000000000000000e+000;
00400 q[10] = 0.00000000000000000e+000;
00401 q[18] = 0.00000000000000000e+000;
00402 q[26] = 0.00000000000000000e+000;
00403 q[3] = 0.00000000000000000e+000;
00404 q[11] = 0.00000000000000000e+000;
00405 q[19] = 0.00000000000000000e+000;
00406 q[27] = 0.00000000000000000e+000;
00407 }
00408 static void init_work( lapack_int size, double *work ) {
00409 lapack_int i;
00410 for( i = 0; i < size; i++ ) {
00411 work[i] = 0;
00412 }
00413 }
00414
00415
00416
00417 static int compare_dsbtrd( double *ab, double *ab_i, double *d, double *d_i,
00418 double *e, double *e_i, double *q, double *q_i,
00419 lapack_int info, lapack_int info_i, lapack_int ldab,
00420 lapack_int ldq, lapack_int n, char vect )
00421 {
00422 lapack_int i;
00423 int failed = 0;
00424 for( i = 0; i < ldab*n; i++ ) {
00425 failed += compare_doubles(ab[i],ab_i[i]);
00426 }
00427 for( i = 0; i < n; i++ ) {
00428 failed += compare_doubles(d[i],d_i[i]);
00429 }
00430 for( i = 0; i < (n-1); i++ ) {
00431 failed += compare_doubles(e[i],e_i[i]);
00432 }
00433 if( LAPACKE_lsame( vect, 'u' ) || LAPACKE_lsame( vect, 'v' ) ) {
00434 for( i = 0; i < ldq*n; i++ ) {
00435 failed += compare_doubles(q[i],q_i[i]);
00436 }
00437 }
00438 failed += (info == info_i) ? 0 : 1;
00439 if( info != 0 || info_i != 0 ) {
00440 printf( "info=%d, info_i=%d\n",(int)info,(int)info_i );
00441 }
00442
00443 return failed;
00444 }