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_zpprfs( char *uplo, lapack_int *n, lapack_int *nrhs,
00055 lapack_int *ldb, lapack_int *ldx );
00056 static void init_ap( lapack_int size, lapack_complex_double *ap );
00057 static void init_afp( lapack_int size, lapack_complex_double *afp );
00058 static void init_b( lapack_int size, lapack_complex_double *b );
00059 static void init_x( lapack_int size, lapack_complex_double *x );
00060 static void init_ferr( lapack_int size, double *ferr );
00061 static void init_berr( lapack_int size, double *berr );
00062 static void init_work( lapack_int size, lapack_complex_double *work );
00063 static void init_rwork( lapack_int size, double *rwork );
00064 static int compare_zpprfs( lapack_complex_double *x, lapack_complex_double *x_i,
00065 double *ferr, double *ferr_i, double *berr,
00066 double *berr_i, lapack_int info, lapack_int info_i,
00067 lapack_int ldx, lapack_int nrhs );
00068
00069 int main(void)
00070 {
00071
00072 char uplo, uplo_i;
00073 lapack_int n, n_i;
00074 lapack_int nrhs, nrhs_i;
00075 lapack_int ldb, ldb_i;
00076 lapack_int ldb_r;
00077 lapack_int ldx, ldx_i;
00078 lapack_int ldx_r;
00079 lapack_int info, info_i;
00080 lapack_int i;
00081 int failed;
00082
00083
00084 lapack_complex_double *ap = NULL, *ap_i = NULL;
00085 lapack_complex_double *afp = NULL, *afp_i = NULL;
00086 lapack_complex_double *b = NULL, *b_i = NULL;
00087 lapack_complex_double *x = NULL, *x_i = NULL;
00088 double *ferr = NULL, *ferr_i = NULL;
00089 double *berr = NULL, *berr_i = NULL;
00090 lapack_complex_double *work = NULL, *work_i = NULL;
00091 double *rwork = NULL, *rwork_i = NULL;
00092 lapack_complex_double *x_save = NULL;
00093 double *ferr_save = NULL;
00094 double *berr_save = NULL;
00095 lapack_complex_double *ap_r = NULL;
00096 lapack_complex_double *afp_r = NULL;
00097 lapack_complex_double *b_r = NULL;
00098 lapack_complex_double *x_r = NULL;
00099
00100
00101 init_scalars_zpprfs( &uplo, &n, &nrhs, &ldb, &ldx );
00102 ldb_r = nrhs+2;
00103 ldx_r = nrhs+2;
00104 uplo_i = uplo;
00105 n_i = n;
00106 nrhs_i = nrhs;
00107 ldb_i = ldb;
00108 ldx_i = ldx;
00109
00110
00111 ap = (lapack_complex_double *)
00112 LAPACKE_malloc( ((n*(n+1)/2)) * sizeof(lapack_complex_double) );
00113 afp = (lapack_complex_double *)
00114 LAPACKE_malloc( ((n*(n+1)/2)) * sizeof(lapack_complex_double) );
00115 b = (lapack_complex_double *)
00116 LAPACKE_malloc( ldb*nrhs * sizeof(lapack_complex_double) );
00117 x = (lapack_complex_double *)
00118 LAPACKE_malloc( ldx*nrhs * sizeof(lapack_complex_double) );
00119 ferr = (double *)LAPACKE_malloc( nrhs * sizeof(double) );
00120 berr = (double *)LAPACKE_malloc( nrhs * sizeof(double) );
00121 work = (lapack_complex_double *)
00122 LAPACKE_malloc( 2*n * sizeof(lapack_complex_double) );
00123 rwork = (double *)LAPACKE_malloc( n * sizeof(double) );
00124
00125
00126 ap_i = (lapack_complex_double *)
00127 LAPACKE_malloc( ((n*(n+1)/2)) * sizeof(lapack_complex_double) );
00128 afp_i = (lapack_complex_double *)
00129 LAPACKE_malloc( ((n*(n+1)/2)) * sizeof(lapack_complex_double) );
00130 b_i = (lapack_complex_double *)
00131 LAPACKE_malloc( ldb*nrhs * sizeof(lapack_complex_double) );
00132 x_i = (lapack_complex_double *)
00133 LAPACKE_malloc( ldx*nrhs * sizeof(lapack_complex_double) );
00134 ferr_i = (double *)LAPACKE_malloc( nrhs * sizeof(double) );
00135 berr_i = (double *)LAPACKE_malloc( nrhs * sizeof(double) );
00136 work_i = (lapack_complex_double *)
00137 LAPACKE_malloc( 2*n * sizeof(lapack_complex_double) );
00138 rwork_i = (double *)LAPACKE_malloc( n * sizeof(double) );
00139
00140
00141 x_save = (lapack_complex_double *)
00142 LAPACKE_malloc( ldx*nrhs * sizeof(lapack_complex_double) );
00143 ferr_save = (double *)LAPACKE_malloc( nrhs * sizeof(double) );
00144 berr_save = (double *)LAPACKE_malloc( nrhs * sizeof(double) );
00145
00146
00147 ap_r = (lapack_complex_double *)
00148 LAPACKE_malloc( n*(n+1)/2 * sizeof(lapack_complex_double) );
00149 afp_r = (lapack_complex_double *)
00150 LAPACKE_malloc( n*(n+1)/2 * sizeof(lapack_complex_double) );
00151 b_r = (lapack_complex_double *)
00152 LAPACKE_malloc( n*(nrhs+2) * sizeof(lapack_complex_double) );
00153 x_r = (lapack_complex_double *)
00154 LAPACKE_malloc( n*(nrhs+2) * sizeof(lapack_complex_double) );
00155
00156
00157 init_ap( (n*(n+1)/2), ap );
00158 init_afp( (n*(n+1)/2), afp );
00159 init_b( ldb*nrhs, b );
00160 init_x( ldx*nrhs, x );
00161 init_ferr( nrhs, ferr );
00162 init_berr( nrhs, berr );
00163 init_work( 2*n, work );
00164 init_rwork( n, rwork );
00165
00166
00167 for( i = 0; i < ldx*nrhs; i++ ) {
00168 x_save[i] = x[i];
00169 }
00170 for( i = 0; i < nrhs; i++ ) {
00171 ferr_save[i] = ferr[i];
00172 }
00173 for( i = 0; i < nrhs; i++ ) {
00174 berr_save[i] = berr[i];
00175 }
00176
00177
00178 zpprfs_( &uplo, &n, &nrhs, ap, afp, b, &ldb, x, &ldx, ferr, berr, work,
00179 rwork, &info );
00180
00181
00182
00183 for( i = 0; i < (n*(n+1)/2); i++ ) {
00184 ap_i[i] = ap[i];
00185 }
00186 for( i = 0; i < (n*(n+1)/2); i++ ) {
00187 afp_i[i] = afp[i];
00188 }
00189 for( i = 0; i < ldb*nrhs; i++ ) {
00190 b_i[i] = b[i];
00191 }
00192 for( i = 0; i < ldx*nrhs; i++ ) {
00193 x_i[i] = x_save[i];
00194 }
00195 for( i = 0; i < nrhs; i++ ) {
00196 ferr_i[i] = ferr_save[i];
00197 }
00198 for( i = 0; i < nrhs; i++ ) {
00199 berr_i[i] = berr_save[i];
00200 }
00201 for( i = 0; i < 2*n; i++ ) {
00202 work_i[i] = work[i];
00203 }
00204 for( i = 0; i < n; i++ ) {
00205 rwork_i[i] = rwork[i];
00206 }
00207 info_i = LAPACKE_zpprfs_work( LAPACK_COL_MAJOR, uplo_i, n_i, nrhs_i, ap_i,
00208 afp_i, b_i, ldb_i, x_i, ldx_i, ferr_i, berr_i,
00209 work_i, rwork_i );
00210
00211 failed = compare_zpprfs( x, x_i, ferr, ferr_i, berr, berr_i, info, info_i,
00212 ldx, nrhs );
00213 if( failed == 0 ) {
00214 printf( "PASSED: column-major middle-level interface to zpprfs\n" );
00215 } else {
00216 printf( "FAILED: column-major middle-level interface to zpprfs\n" );
00217 }
00218
00219
00220
00221 for( i = 0; i < (n*(n+1)/2); i++ ) {
00222 ap_i[i] = ap[i];
00223 }
00224 for( i = 0; i < (n*(n+1)/2); i++ ) {
00225 afp_i[i] = afp[i];
00226 }
00227 for( i = 0; i < ldb*nrhs; i++ ) {
00228 b_i[i] = b[i];
00229 }
00230 for( i = 0; i < ldx*nrhs; i++ ) {
00231 x_i[i] = x_save[i];
00232 }
00233 for( i = 0; i < nrhs; i++ ) {
00234 ferr_i[i] = ferr_save[i];
00235 }
00236 for( i = 0; i < nrhs; i++ ) {
00237 berr_i[i] = berr_save[i];
00238 }
00239 for( i = 0; i < 2*n; i++ ) {
00240 work_i[i] = work[i];
00241 }
00242 for( i = 0; i < n; i++ ) {
00243 rwork_i[i] = rwork[i];
00244 }
00245 info_i = LAPACKE_zpprfs( LAPACK_COL_MAJOR, uplo_i, n_i, nrhs_i, ap_i, afp_i,
00246 b_i, ldb_i, x_i, ldx_i, ferr_i, berr_i );
00247
00248 failed = compare_zpprfs( x, x_i, ferr, ferr_i, berr, berr_i, info, info_i,
00249 ldx, nrhs );
00250 if( failed == 0 ) {
00251 printf( "PASSED: column-major high-level interface to zpprfs\n" );
00252 } else {
00253 printf( "FAILED: column-major high-level interface to zpprfs\n" );
00254 }
00255
00256
00257
00258 for( i = 0; i < (n*(n+1)/2); i++ ) {
00259 ap_i[i] = ap[i];
00260 }
00261 for( i = 0; i < (n*(n+1)/2); i++ ) {
00262 afp_i[i] = afp[i];
00263 }
00264 for( i = 0; i < ldb*nrhs; i++ ) {
00265 b_i[i] = b[i];
00266 }
00267 for( i = 0; i < ldx*nrhs; i++ ) {
00268 x_i[i] = x_save[i];
00269 }
00270 for( i = 0; i < nrhs; i++ ) {
00271 ferr_i[i] = ferr_save[i];
00272 }
00273 for( i = 0; i < nrhs; i++ ) {
00274 berr_i[i] = berr_save[i];
00275 }
00276 for( i = 0; i < 2*n; i++ ) {
00277 work_i[i] = work[i];
00278 }
00279 for( i = 0; i < n; i++ ) {
00280 rwork_i[i] = rwork[i];
00281 }
00282
00283 LAPACKE_zpp_trans( LAPACK_COL_MAJOR, uplo, n, ap_i, ap_r );
00284 LAPACKE_zpp_trans( LAPACK_COL_MAJOR, uplo, n, afp_i, afp_r );
00285 LAPACKE_zge_trans( LAPACK_COL_MAJOR, n, nrhs, b_i, ldb, b_r, nrhs+2 );
00286 LAPACKE_zge_trans( LAPACK_COL_MAJOR, n, nrhs, x_i, ldx, x_r, nrhs+2 );
00287 info_i = LAPACKE_zpprfs_work( LAPACK_ROW_MAJOR, uplo_i, n_i, nrhs_i, ap_r,
00288 afp_r, b_r, ldb_r, x_r, ldx_r, ferr_i, berr_i,
00289 work_i, rwork_i );
00290
00291 LAPACKE_zge_trans( LAPACK_ROW_MAJOR, n, nrhs, x_r, nrhs+2, x_i, ldx );
00292
00293 failed = compare_zpprfs( x, x_i, ferr, ferr_i, berr, berr_i, info, info_i,
00294 ldx, nrhs );
00295 if( failed == 0 ) {
00296 printf( "PASSED: row-major middle-level interface to zpprfs\n" );
00297 } else {
00298 printf( "FAILED: row-major middle-level interface to zpprfs\n" );
00299 }
00300
00301
00302
00303 for( i = 0; i < (n*(n+1)/2); i++ ) {
00304 ap_i[i] = ap[i];
00305 }
00306 for( i = 0; i < (n*(n+1)/2); i++ ) {
00307 afp_i[i] = afp[i];
00308 }
00309 for( i = 0; i < ldb*nrhs; i++ ) {
00310 b_i[i] = b[i];
00311 }
00312 for( i = 0; i < ldx*nrhs; i++ ) {
00313 x_i[i] = x_save[i];
00314 }
00315 for( i = 0; i < nrhs; i++ ) {
00316 ferr_i[i] = ferr_save[i];
00317 }
00318 for( i = 0; i < nrhs; i++ ) {
00319 berr_i[i] = berr_save[i];
00320 }
00321 for( i = 0; i < 2*n; i++ ) {
00322 work_i[i] = work[i];
00323 }
00324 for( i = 0; i < n; i++ ) {
00325 rwork_i[i] = rwork[i];
00326 }
00327
00328
00329 LAPACKE_zpp_trans( LAPACK_COL_MAJOR, uplo, n, ap_i, ap_r );
00330 LAPACKE_zpp_trans( LAPACK_COL_MAJOR, uplo, n, afp_i, afp_r );
00331 LAPACKE_zge_trans( LAPACK_COL_MAJOR, n, nrhs, b_i, ldb, b_r, nrhs+2 );
00332 LAPACKE_zge_trans( LAPACK_COL_MAJOR, n, nrhs, x_i, ldx, x_r, nrhs+2 );
00333 info_i = LAPACKE_zpprfs( LAPACK_ROW_MAJOR, uplo_i, n_i, nrhs_i, ap_r, afp_r,
00334 b_r, ldb_r, x_r, ldx_r, ferr_i, berr_i );
00335
00336 LAPACKE_zge_trans( LAPACK_ROW_MAJOR, n, nrhs, x_r, nrhs+2, x_i, ldx );
00337
00338 failed = compare_zpprfs( x, x_i, ferr, ferr_i, berr, berr_i, info, info_i,
00339 ldx, nrhs );
00340 if( failed == 0 ) {
00341 printf( "PASSED: row-major high-level interface to zpprfs\n" );
00342 } else {
00343 printf( "FAILED: row-major high-level interface to zpprfs\n" );
00344 }
00345
00346
00347 if( ap != NULL ) {
00348 LAPACKE_free( ap );
00349 }
00350 if( ap_i != NULL ) {
00351 LAPACKE_free( ap_i );
00352 }
00353 if( ap_r != NULL ) {
00354 LAPACKE_free( ap_r );
00355 }
00356 if( afp != NULL ) {
00357 LAPACKE_free( afp );
00358 }
00359 if( afp_i != NULL ) {
00360 LAPACKE_free( afp_i );
00361 }
00362 if( afp_r != NULL ) {
00363 LAPACKE_free( afp_r );
00364 }
00365 if( b != NULL ) {
00366 LAPACKE_free( b );
00367 }
00368 if( b_i != NULL ) {
00369 LAPACKE_free( b_i );
00370 }
00371 if( b_r != NULL ) {
00372 LAPACKE_free( b_r );
00373 }
00374 if( x != NULL ) {
00375 LAPACKE_free( x );
00376 }
00377 if( x_i != NULL ) {
00378 LAPACKE_free( x_i );
00379 }
00380 if( x_r != NULL ) {
00381 LAPACKE_free( x_r );
00382 }
00383 if( x_save != NULL ) {
00384 LAPACKE_free( x_save );
00385 }
00386 if( ferr != NULL ) {
00387 LAPACKE_free( ferr );
00388 }
00389 if( ferr_i != NULL ) {
00390 LAPACKE_free( ferr_i );
00391 }
00392 if( ferr_save != NULL ) {
00393 LAPACKE_free( ferr_save );
00394 }
00395 if( berr != NULL ) {
00396 LAPACKE_free( berr );
00397 }
00398 if( berr_i != NULL ) {
00399 LAPACKE_free( berr_i );
00400 }
00401 if( berr_save != NULL ) {
00402 LAPACKE_free( berr_save );
00403 }
00404 if( work != NULL ) {
00405 LAPACKE_free( work );
00406 }
00407 if( work_i != NULL ) {
00408 LAPACKE_free( work_i );
00409 }
00410 if( rwork != NULL ) {
00411 LAPACKE_free( rwork );
00412 }
00413 if( rwork_i != NULL ) {
00414 LAPACKE_free( rwork_i );
00415 }
00416
00417 return 0;
00418 }
00419
00420
00421 static void init_scalars_zpprfs( char *uplo, lapack_int *n, lapack_int *nrhs,
00422 lapack_int *ldb, lapack_int *ldx )
00423 {
00424 *uplo = 'L';
00425 *n = 4;
00426 *nrhs = 2;
00427 *ldb = 8;
00428 *ldx = 8;
00429
00430 return;
00431 }
00432
00433
00434 static void init_ap( lapack_int size, lapack_complex_double *ap ) {
00435 lapack_int i;
00436 for( i = 0; i < size; i++ ) {
00437 ap[i] = lapack_make_complex_double( 0.0, 0.0 );
00438 }
00439 ap[0] = lapack_make_complex_double( 3.23000000000000000e+000,
00440 0.00000000000000000e+000 );
00441 ap[1] = lapack_make_complex_double( 1.51000000000000000e+000,
00442 1.91999999999999990e+000 );
00443 ap[2] = lapack_make_complex_double( 1.89999999999999990e+000,
00444 -8.39999999999999970e-001 );
00445 ap[3] = lapack_make_complex_double( 4.19999999999999980e-001,
00446 -2.50000000000000000e+000 );
00447 ap[4] = lapack_make_complex_double( 3.58000000000000010e+000,
00448 0.00000000000000000e+000 );
00449 ap[5] = lapack_make_complex_double( -2.30000000000000010e-001,
00450 -1.11000000000000010e+000 );
00451 ap[6] = lapack_make_complex_double( -1.17999999999999990e+000,
00452 -1.37000000000000010e+000 );
00453 ap[7] = lapack_make_complex_double( 4.08999999999999990e+000,
00454 0.00000000000000000e+000 );
00455 ap[8] = lapack_make_complex_double( 2.33000000000000010e+000,
00456 1.40000000000000010e-001 );
00457 ap[9] = lapack_make_complex_double( 4.29000000000000000e+000,
00458 0.00000000000000000e+000 );
00459 }
00460 static void init_afp( lapack_int size, lapack_complex_double *afp ) {
00461 lapack_int i;
00462 for( i = 0; i < size; i++ ) {
00463 afp[i] = lapack_make_complex_double( 0.0, 0.0 );
00464 }
00465 afp[0] = lapack_make_complex_double( 1.79722007556114290e+000,
00466 0.00000000000000000e+000 );
00467 afp[1] = lapack_make_complex_double( 8.40186474952732460e-001,
00468 1.06831657742334190e+000 );
00469 afp[2] = lapack_make_complex_double( 1.05718827974184860e+000,
00470 -4.67388502622712030e-001 );
00471 afp[3] = lapack_make_complex_double( 2.33694251311356020e-001,
00472 -1.39103721018664310e+000 );
00473 afp[4] = lapack_make_complex_double( 1.31635343950968520e+000,
00474 0.00000000000000000e+000 );
00475 afp[5] = lapack_make_complex_double( -4.70174947010632950e-001,
00476 3.13065815599946400e-001 );
00477 afp[6] = lapack_make_complex_double( 8.33525092394419580e-002,
00478 3.67607144303747430e-002 );
00479 afp[7] = lapack_make_complex_double( 1.56039297713712430e+000,
00480 0.00000000000000000e+000 );
00481 afp[8] = lapack_make_complex_double( 9.35961733792340160e-001,
00482 9.89969219281573890e-001 );
00483 afp[9] = lapack_make_complex_double( 6.60333297365588770e-001,
00484 0.00000000000000000e+000 );
00485 }
00486 static void init_b( lapack_int size, lapack_complex_double *b ) {
00487 lapack_int i;
00488 for( i = 0; i < size; i++ ) {
00489 b[i] = lapack_make_complex_double( 0.0, 0.0 );
00490 }
00491 b[0] = lapack_make_complex_double( 3.93000000000000020e+000,
00492 -6.13999999999999970e+000 );
00493 b[8] = lapack_make_complex_double( 1.48000000000000000e+000,
00494 6.58000000000000010e+000 );
00495 b[1] = lapack_make_complex_double( 6.16999999999999990e+000,
00496 9.41999999999999990e+000 );
00497 b[9] = lapack_make_complex_double( 4.65000000000000040e+000,
00498 -4.75000000000000000e+000 );
00499 b[2] = lapack_make_complex_double( -7.16999999999999990e+000,
00500 -2.18299999999999980e+001 );
00501 b[10] = lapack_make_complex_double( -4.91000000000000010e+000,
00502 2.29000000000000000e+000 );
00503 b[3] = lapack_make_complex_double( 1.99000000000000000e+000,
00504 -1.43800000000000010e+001 );
00505 b[11] = lapack_make_complex_double( 7.63999999999999970e+000,
00506 -1.07899999999999990e+001 );
00507 }
00508 static void init_x( lapack_int size, lapack_complex_double *x ) {
00509 lapack_int i;
00510 for( i = 0; i < size; i++ ) {
00511 x[i] = lapack_make_complex_double( 0.0, 0.0 );
00512 }
00513 x[0] = lapack_make_complex_double( 1.00000000000000470e+000,
00514 -1.00000000000000890e+000 );
00515 x[8] = lapack_make_complex_double( -1.00000000000001090e+000,
00516 2.00000000000000490e+000 );
00517 x[1] = lapack_make_complex_double( -4.30138081129436570e-015,
00518 3.00000000000000220e+000 );
00519 x[9] = lapack_make_complex_double( 3.00000000000000440e+000,
00520 -3.99999999999999780e+000 );
00521 x[2] = lapack_make_complex_double( -4.00000000000000530e+000,
00522 -4.99999999999999640e+000 );
00523 x[10] = lapack_make_complex_double( -1.99999999999999180e+000,
00524 3.00000000000000040e+000 );
00525 x[3] = lapack_make_complex_double( 2.00000000000000580e+000,
00526 1.00000000000000110e+000 );
00527 x[11] = lapack_make_complex_double( 3.99999999999999380e+000,
00528 -5.00000000000000440e+000 );
00529 }
00530 static void init_ferr( lapack_int size, double *ferr ) {
00531 lapack_int i;
00532 for( i = 0; i < size; i++ ) {
00533 ferr[i] = 0;
00534 }
00535 }
00536 static void init_berr( lapack_int size, double *berr ) {
00537 lapack_int i;
00538 for( i = 0; i < size; i++ ) {
00539 berr[i] = 0;
00540 }
00541 }
00542 static void init_work( lapack_int size, lapack_complex_double *work ) {
00543 lapack_int i;
00544 for( i = 0; i < size; i++ ) {
00545 work[i] = lapack_make_complex_double( 0.0, 0.0 );
00546 }
00547 }
00548 static void init_rwork( lapack_int size, double *rwork ) {
00549 lapack_int i;
00550 for( i = 0; i < size; i++ ) {
00551 rwork[i] = 0;
00552 }
00553 }
00554
00555
00556
00557 static int compare_zpprfs( lapack_complex_double *x, lapack_complex_double *x_i,
00558 double *ferr, double *ferr_i, double *berr,
00559 double *berr_i, lapack_int info, lapack_int info_i,
00560 lapack_int ldx, lapack_int nrhs )
00561 {
00562 lapack_int i;
00563 int failed = 0;
00564 for( i = 0; i < ldx*nrhs; i++ ) {
00565 failed += compare_complex_doubles(x[i],x_i[i]);
00566 }
00567 for( i = 0; i < nrhs; i++ ) {
00568 failed += compare_doubles(ferr[i],ferr_i[i]);
00569 }
00570 for( i = 0; i < nrhs; i++ ) {
00571 failed += compare_doubles(berr[i],berr_i[i]);
00572 }
00573 failed += (info == info_i) ? 0 : 1;
00574 if( info != 0 || info_i != 0 ) {
00575 printf( "info=%d, info_i=%d\n",(int)info,(int)info_i );
00576 }
00577
00578 return failed;
00579 }