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_stbrfs( char *uplo, char *trans, char *diag,
00055 lapack_int *n, lapack_int *kd,
00056 lapack_int *nrhs, lapack_int *ldab,
00057 lapack_int *ldb, lapack_int *ldx );
00058 static void init_ab( lapack_int size, float *ab );
00059 static void init_b( lapack_int size, float *b );
00060 static void init_x( lapack_int size, float *x );
00061 static void init_ferr( lapack_int size, float *ferr );
00062 static void init_berr( lapack_int size, float *berr );
00063 static void init_work( lapack_int size, float *work );
00064 static void init_iwork( lapack_int size, lapack_int *iwork );
00065 static int compare_stbrfs( float *ferr, float *ferr_i, float *berr,
00066 float *berr_i, lapack_int info, lapack_int info_i,
00067 lapack_int nrhs );
00068
00069 int main(void)
00070 {
00071
00072 char uplo, uplo_i;
00073 char trans, trans_i;
00074 char diag, diag_i;
00075 lapack_int n, n_i;
00076 lapack_int kd, kd_i;
00077 lapack_int nrhs, nrhs_i;
00078 lapack_int ldab, ldab_i;
00079 lapack_int ldab_r;
00080 lapack_int ldb, ldb_i;
00081 lapack_int ldb_r;
00082 lapack_int ldx, ldx_i;
00083 lapack_int ldx_r;
00084 lapack_int info, info_i;
00085 lapack_int i;
00086 int failed;
00087
00088
00089 float *ab = NULL, *ab_i = NULL;
00090 float *b = NULL, *b_i = NULL;
00091 float *x = NULL, *x_i = NULL;
00092 float *ferr = NULL, *ferr_i = NULL;
00093 float *berr = NULL, *berr_i = NULL;
00094 float *work = NULL, *work_i = NULL;
00095 lapack_int *iwork = NULL, *iwork_i = NULL;
00096 float *ferr_save = NULL;
00097 float *berr_save = NULL;
00098 float *ab_r = NULL;
00099 float *b_r = NULL;
00100 float *x_r = NULL;
00101
00102
00103 init_scalars_stbrfs( &uplo, &trans, &diag, &n, &kd, &nrhs, &ldab, &ldb,
00104 &ldx );
00105 ldab_r = n+2;
00106 ldb_r = nrhs+2;
00107 ldx_r = nrhs+2;
00108 uplo_i = uplo;
00109 trans_i = trans;
00110 diag_i = diag;
00111 n_i = n;
00112 kd_i = kd;
00113 nrhs_i = nrhs;
00114 ldab_i = ldab;
00115 ldb_i = ldb;
00116 ldx_i = ldx;
00117
00118
00119 ab = (float *)LAPACKE_malloc( ldab*n * sizeof(float) );
00120 b = (float *)LAPACKE_malloc( ldb*nrhs * sizeof(float) );
00121 x = (float *)LAPACKE_malloc( ldx*nrhs * sizeof(float) );
00122 ferr = (float *)LAPACKE_malloc( nrhs * sizeof(float) );
00123 berr = (float *)LAPACKE_malloc( nrhs * sizeof(float) );
00124 work = (float *)LAPACKE_malloc( 3*n * sizeof(float) );
00125 iwork = (lapack_int *)LAPACKE_malloc( n * sizeof(lapack_int) );
00126
00127
00128 ab_i = (float *)LAPACKE_malloc( ldab*n * sizeof(float) );
00129 b_i = (float *)LAPACKE_malloc( ldb*nrhs * sizeof(float) );
00130 x_i = (float *)LAPACKE_malloc( ldx*nrhs * sizeof(float) );
00131 ferr_i = (float *)LAPACKE_malloc( nrhs * sizeof(float) );
00132 berr_i = (float *)LAPACKE_malloc( nrhs * sizeof(float) );
00133 work_i = (float *)LAPACKE_malloc( 3*n * sizeof(float) );
00134 iwork_i = (lapack_int *)LAPACKE_malloc( n * sizeof(lapack_int) );
00135
00136
00137 ferr_save = (float *)LAPACKE_malloc( nrhs * sizeof(float) );
00138 berr_save = (float *)LAPACKE_malloc( nrhs * sizeof(float) );
00139
00140
00141 ab_r = (float *)LAPACKE_malloc( (kd+1)*(n+2) * sizeof(float) );
00142 b_r = (float *)LAPACKE_malloc( n*(nrhs+2) * sizeof(float) );
00143 x_r = (float *)LAPACKE_malloc( n*(nrhs+2) * sizeof(float) );
00144
00145
00146 init_ab( ldab*n, ab );
00147 init_b( ldb*nrhs, b );
00148 init_x( ldx*nrhs, x );
00149 init_ferr( nrhs, ferr );
00150 init_berr( nrhs, berr );
00151 init_work( 3*n, work );
00152 init_iwork( n, iwork );
00153
00154
00155 for( i = 0; i < nrhs; i++ ) {
00156 ferr_save[i] = ferr[i];
00157 }
00158 for( i = 0; i < nrhs; i++ ) {
00159 berr_save[i] = berr[i];
00160 }
00161
00162
00163 stbrfs_( &uplo, &trans, &diag, &n, &kd, &nrhs, ab, &ldab, b, &ldb, x, &ldx,
00164 ferr, berr, work, iwork, &info );
00165
00166
00167
00168 for( i = 0; i < ldab*n; i++ ) {
00169 ab_i[i] = ab[i];
00170 }
00171 for( i = 0; i < ldb*nrhs; i++ ) {
00172 b_i[i] = b[i];
00173 }
00174 for( i = 0; i < ldx*nrhs; i++ ) {
00175 x_i[i] = x[i];
00176 }
00177 for( i = 0; i < nrhs; i++ ) {
00178 ferr_i[i] = ferr_save[i];
00179 }
00180 for( i = 0; i < nrhs; i++ ) {
00181 berr_i[i] = berr_save[i];
00182 }
00183 for( i = 0; i < 3*n; i++ ) {
00184 work_i[i] = work[i];
00185 }
00186 for( i = 0; i < n; i++ ) {
00187 iwork_i[i] = iwork[i];
00188 }
00189 info_i = LAPACKE_stbrfs_work( LAPACK_COL_MAJOR, uplo_i, trans_i, diag_i,
00190 n_i, kd_i, nrhs_i, ab_i, ldab_i, b_i, ldb_i,
00191 x_i, ldx_i, ferr_i, berr_i, work_i, iwork_i );
00192
00193 failed = compare_stbrfs( ferr, ferr_i, berr, berr_i, info, info_i, nrhs );
00194 if( failed == 0 ) {
00195 printf( "PASSED: column-major middle-level interface to stbrfs\n" );
00196 } else {
00197 printf( "FAILED: column-major middle-level interface to stbrfs\n" );
00198 }
00199
00200
00201
00202 for( i = 0; i < ldab*n; i++ ) {
00203 ab_i[i] = ab[i];
00204 }
00205 for( i = 0; i < ldb*nrhs; i++ ) {
00206 b_i[i] = b[i];
00207 }
00208 for( i = 0; i < ldx*nrhs; i++ ) {
00209 x_i[i] = x[i];
00210 }
00211 for( i = 0; i < nrhs; i++ ) {
00212 ferr_i[i] = ferr_save[i];
00213 }
00214 for( i = 0; i < nrhs; i++ ) {
00215 berr_i[i] = berr_save[i];
00216 }
00217 for( i = 0; i < 3*n; i++ ) {
00218 work_i[i] = work[i];
00219 }
00220 for( i = 0; i < n; i++ ) {
00221 iwork_i[i] = iwork[i];
00222 }
00223 info_i = LAPACKE_stbrfs( LAPACK_COL_MAJOR, uplo_i, trans_i, diag_i, n_i,
00224 kd_i, nrhs_i, ab_i, ldab_i, b_i, ldb_i, x_i, ldx_i,
00225 ferr_i, berr_i );
00226
00227 failed = compare_stbrfs( ferr, ferr_i, berr, berr_i, info, info_i, nrhs );
00228 if( failed == 0 ) {
00229 printf( "PASSED: column-major high-level interface to stbrfs\n" );
00230 } else {
00231 printf( "FAILED: column-major high-level interface to stbrfs\n" );
00232 }
00233
00234
00235
00236 for( i = 0; i < ldab*n; i++ ) {
00237 ab_i[i] = ab[i];
00238 }
00239 for( i = 0; i < ldb*nrhs; i++ ) {
00240 b_i[i] = b[i];
00241 }
00242 for( i = 0; i < ldx*nrhs; i++ ) {
00243 x_i[i] = x[i];
00244 }
00245 for( i = 0; i < nrhs; i++ ) {
00246 ferr_i[i] = ferr_save[i];
00247 }
00248 for( i = 0; i < nrhs; i++ ) {
00249 berr_i[i] = berr_save[i];
00250 }
00251 for( i = 0; i < 3*n; i++ ) {
00252 work_i[i] = work[i];
00253 }
00254 for( i = 0; i < n; i++ ) {
00255 iwork_i[i] = iwork[i];
00256 }
00257
00258 LAPACKE_sge_trans( LAPACK_COL_MAJOR, kd+1, n, ab_i, ldab, ab_r, n+2 );
00259 LAPACKE_sge_trans( LAPACK_COL_MAJOR, n, nrhs, b_i, ldb, b_r, nrhs+2 );
00260 LAPACKE_sge_trans( LAPACK_COL_MAJOR, n, nrhs, x_i, ldx, x_r, nrhs+2 );
00261 info_i = LAPACKE_stbrfs_work( LAPACK_ROW_MAJOR, uplo_i, trans_i, diag_i,
00262 n_i, kd_i, nrhs_i, ab_r, ldab_r, b_r, ldb_r,
00263 x_r, ldx_r, ferr_i, berr_i, work_i, iwork_i );
00264
00265 failed = compare_stbrfs( ferr, ferr_i, berr, berr_i, info, info_i, nrhs );
00266 if( failed == 0 ) {
00267 printf( "PASSED: row-major middle-level interface to stbrfs\n" );
00268 } else {
00269 printf( "FAILED: row-major middle-level interface to stbrfs\n" );
00270 }
00271
00272
00273
00274 for( i = 0; i < ldab*n; i++ ) {
00275 ab_i[i] = ab[i];
00276 }
00277 for( i = 0; i < ldb*nrhs; i++ ) {
00278 b_i[i] = b[i];
00279 }
00280 for( i = 0; i < ldx*nrhs; i++ ) {
00281 x_i[i] = x[i];
00282 }
00283 for( i = 0; i < nrhs; i++ ) {
00284 ferr_i[i] = ferr_save[i];
00285 }
00286 for( i = 0; i < nrhs; i++ ) {
00287 berr_i[i] = berr_save[i];
00288 }
00289 for( i = 0; i < 3*n; i++ ) {
00290 work_i[i] = work[i];
00291 }
00292 for( i = 0; i < n; i++ ) {
00293 iwork_i[i] = iwork[i];
00294 }
00295
00296
00297 LAPACKE_sge_trans( LAPACK_COL_MAJOR, kd+1, n, ab_i, ldab, ab_r, n+2 );
00298 LAPACKE_sge_trans( LAPACK_COL_MAJOR, n, nrhs, b_i, ldb, b_r, nrhs+2 );
00299 LAPACKE_sge_trans( LAPACK_COL_MAJOR, n, nrhs, x_i, ldx, x_r, nrhs+2 );
00300 info_i = LAPACKE_stbrfs( LAPACK_ROW_MAJOR, uplo_i, trans_i, diag_i, n_i,
00301 kd_i, nrhs_i, ab_r, ldab_r, b_r, ldb_r, x_r, ldx_r,
00302 ferr_i, berr_i );
00303
00304 failed = compare_stbrfs( ferr, ferr_i, berr, berr_i, info, info_i, nrhs );
00305 if( failed == 0 ) {
00306 printf( "PASSED: row-major high-level interface to stbrfs\n" );
00307 } else {
00308 printf( "FAILED: row-major high-level interface to stbrfs\n" );
00309 }
00310
00311
00312 if( ab != NULL ) {
00313 LAPACKE_free( ab );
00314 }
00315 if( ab_i != NULL ) {
00316 LAPACKE_free( ab_i );
00317 }
00318 if( ab_r != NULL ) {
00319 LAPACKE_free( ab_r );
00320 }
00321 if( b != NULL ) {
00322 LAPACKE_free( b );
00323 }
00324 if( b_i != NULL ) {
00325 LAPACKE_free( b_i );
00326 }
00327 if( b_r != NULL ) {
00328 LAPACKE_free( b_r );
00329 }
00330 if( x != NULL ) {
00331 LAPACKE_free( x );
00332 }
00333 if( x_i != NULL ) {
00334 LAPACKE_free( x_i );
00335 }
00336 if( x_r != NULL ) {
00337 LAPACKE_free( x_r );
00338 }
00339 if( ferr != NULL ) {
00340 LAPACKE_free( ferr );
00341 }
00342 if( ferr_i != NULL ) {
00343 LAPACKE_free( ferr_i );
00344 }
00345 if( ferr_save != NULL ) {
00346 LAPACKE_free( ferr_save );
00347 }
00348 if( berr != NULL ) {
00349 LAPACKE_free( berr );
00350 }
00351 if( berr_i != NULL ) {
00352 LAPACKE_free( berr_i );
00353 }
00354 if( berr_save != NULL ) {
00355 LAPACKE_free( berr_save );
00356 }
00357 if( work != NULL ) {
00358 LAPACKE_free( work );
00359 }
00360 if( work_i != NULL ) {
00361 LAPACKE_free( work_i );
00362 }
00363 if( iwork != NULL ) {
00364 LAPACKE_free( iwork );
00365 }
00366 if( iwork_i != NULL ) {
00367 LAPACKE_free( iwork_i );
00368 }
00369
00370 return 0;
00371 }
00372
00373
00374 static void init_scalars_stbrfs( char *uplo, char *trans, char *diag,
00375 lapack_int *n, lapack_int *kd,
00376 lapack_int *nrhs, lapack_int *ldab,
00377 lapack_int *ldb, lapack_int *ldx )
00378 {
00379 *uplo = 'L';
00380 *trans = 'N';
00381 *diag = 'N';
00382 *n = 4;
00383 *kd = 1;
00384 *nrhs = 2;
00385 *ldab = 9;
00386 *ldb = 8;
00387 *ldx = 8;
00388
00389 return;
00390 }
00391
00392
00393 static void init_ab( lapack_int size, float *ab ) {
00394 lapack_int i;
00395 for( i = 0; i < size; i++ ) {
00396 ab[i] = 0;
00397 }
00398 ab[0] = -4.159999847e+000;
00399 ab[9] = 4.780000210e+000;
00400 ab[18] = 6.320000172e+000;
00401 ab[27] = 1.599999964e-001;
00402 ab[1] = -2.250000000e+000;
00403 ab[10] = 5.860000134e+000;
00404 ab[19] = -4.820000172e+000;
00405 ab[28] = 0.000000000e+000;
00406 }
00407 static void init_b( lapack_int size, float *b ) {
00408 lapack_int i;
00409 for( i = 0; i < size; i++ ) {
00410 b[i] = 0;
00411 }
00412 b[0] = -1.663999939e+001;
00413 b[8] = -4.159999847e+000;
00414 b[1] = -1.377999973e+001;
00415 b[9] = -1.659000015e+001;
00416 b[2] = 1.310000038e+001;
00417 b[10] = -4.940000057e+000;
00418 b[3] = -1.414000034e+001;
00419 b[11] = -9.960000038e+000;
00420 }
00421 static void init_x( lapack_int size, float *x ) {
00422 lapack_int i;
00423 for( i = 0; i < size; i++ ) {
00424 x[i] = 0;
00425 }
00426 x[0] = 4.000000000e+000;
00427 x[8] = 1.000000000e+000;
00428 x[1] = -9.999998808e-001;
00429 x[9] = -3.000000000e+000;
00430 x[2] = 2.999999762e+000;
00431 x[10] = 1.999999881e+000;
00432 x[3] = 1.999992132e+000;
00433 x[11] = -2.000004053e+000;
00434 }
00435 static void init_ferr( lapack_int size, float *ferr ) {
00436 lapack_int i;
00437 for( i = 0; i < size; i++ ) {
00438 ferr[i] = 0;
00439 }
00440 }
00441 static void init_berr( lapack_int size, float *berr ) {
00442 lapack_int i;
00443 for( i = 0; i < size; i++ ) {
00444 berr[i] = 0;
00445 }
00446 }
00447 static void init_work( lapack_int size, float *work ) {
00448 lapack_int i;
00449 for( i = 0; i < size; i++ ) {
00450 work[i] = 0;
00451 }
00452 }
00453 static void init_iwork( lapack_int size, lapack_int *iwork ) {
00454 lapack_int i;
00455 for( i = 0; i < size; i++ ) {
00456 iwork[i] = 0;
00457 }
00458 }
00459
00460
00461
00462 static int compare_stbrfs( float *ferr, float *ferr_i, float *berr,
00463 float *berr_i, lapack_int info, lapack_int info_i,
00464 lapack_int nrhs )
00465 {
00466 lapack_int i;
00467 int failed = 0;
00468 for( i = 0; i < nrhs; i++ ) {
00469 failed += compare_floats(ferr[i],ferr_i[i]);
00470 }
00471 for( i = 0; i < nrhs; i++ ) {
00472 failed += compare_floats(berr[i],berr_i[i]);
00473 }
00474 failed += (info == info_i) ? 0 : 1;
00475 if( info != 0 || info_i != 0 ) {
00476 printf( "info=%d, info_i=%d\n",(int)info,(int)info_i );
00477 }
00478
00479 return failed;
00480 }