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_dbdsqr( char *uplo, lapack_int *n, lapack_int *ncvt,
00055 lapack_int *nru, lapack_int *ncc,
00056 lapack_int *ldvt, lapack_int *ldu,
00057 lapack_int *ldc );
00058 static void init_d( lapack_int size, double *d );
00059 static void init_e( lapack_int size, double *e );
00060 static void init_vt( lapack_int size, double *vt );
00061 static void init_u( lapack_int size, double *u );
00062 static void init_c( lapack_int size, double *c );
00063 static void init_work( lapack_int size, double *work );
00064 static int compare_dbdsqr( double *d, double *d_i, double *e, double *e_i,
00065 double *vt, double *vt_i, double *u, double *u_i,
00066 double *c, double *c_i, lapack_int info,
00067 lapack_int info_i, lapack_int ldc, lapack_int ldu,
00068 lapack_int ldvt, lapack_int n, lapack_int ncc,
00069 lapack_int ncvt, lapack_int nru );
00070
00071 int main(void)
00072 {
00073
00074 char uplo, uplo_i;
00075 lapack_int n, n_i;
00076 lapack_int ncvt, ncvt_i;
00077 lapack_int nru, nru_i;
00078 lapack_int ncc, ncc_i;
00079 lapack_int ldvt, ldvt_i;
00080 lapack_int ldvt_r;
00081 lapack_int ldu, ldu_i;
00082 lapack_int ldu_r;
00083 lapack_int ldc, ldc_i;
00084 lapack_int ldc_r;
00085 lapack_int info, info_i;
00086 lapack_int i;
00087 int failed;
00088
00089
00090 double *d = NULL, *d_i = NULL;
00091 double *e = NULL, *e_i = NULL;
00092 double *vt = NULL, *vt_i = NULL;
00093 double *u = NULL, *u_i = NULL;
00094 double *c = NULL, *c_i = NULL;
00095 double *work = NULL, *work_i = NULL;
00096 double *d_save = NULL;
00097 double *e_save = NULL;
00098 double *vt_save = NULL;
00099 double *u_save = NULL;
00100 double *c_save = NULL;
00101 double *vt_r = NULL;
00102 double *u_r = NULL;
00103 double *c_r = NULL;
00104
00105
00106 init_scalars_dbdsqr( &uplo, &n, &ncvt, &nru, &ncc, &ldvt, &ldu, &ldc );
00107 ldvt_r = ncvt+2;
00108 ldu_r = n+2;
00109 ldc_r = ncc+2;
00110 uplo_i = uplo;
00111 n_i = n;
00112 ncvt_i = ncvt;
00113 nru_i = nru;
00114 ncc_i = ncc;
00115 ldvt_i = ldvt;
00116 ldu_i = ldu;
00117 ldc_i = ldc;
00118
00119
00120 d = (double *)LAPACKE_malloc( n * sizeof(double) );
00121 e = (double *)LAPACKE_malloc( n * sizeof(double) );
00122 vt = (double *)LAPACKE_malloc( ldvt*ncvt * sizeof(double) );
00123 u = (double *)LAPACKE_malloc( ldu*n * sizeof(double) );
00124 c = (double *)LAPACKE_malloc( ldc*ncc * sizeof(double) );
00125 work = (double *)LAPACKE_malloc( 4*n * sizeof(double) );
00126
00127
00128 d_i = (double *)LAPACKE_malloc( n * sizeof(double) );
00129 e_i = (double *)LAPACKE_malloc( n * sizeof(double) );
00130 vt_i = (double *)LAPACKE_malloc( ldvt*ncvt * sizeof(double) );
00131 u_i = (double *)LAPACKE_malloc( ldu*n * sizeof(double) );
00132 c_i = (double *)LAPACKE_malloc( ldc*ncc * sizeof(double) );
00133 work_i = (double *)LAPACKE_malloc( 4*n * sizeof(double) );
00134
00135
00136 d_save = (double *)LAPACKE_malloc( n * sizeof(double) );
00137 e_save = (double *)LAPACKE_malloc( n * sizeof(double) );
00138 vt_save = (double *)LAPACKE_malloc( ldvt*ncvt * sizeof(double) );
00139 u_save = (double *)LAPACKE_malloc( ldu*n * sizeof(double) );
00140 c_save = (double *)LAPACKE_malloc( ldc*ncc * sizeof(double) );
00141
00142
00143 vt_r = (double *)LAPACKE_malloc( n*(ncvt+2) * sizeof(double) );
00144 u_r = (double *)LAPACKE_malloc( nru*(n+2) * sizeof(double) );
00145 c_r = (double *)LAPACKE_malloc( n*(ncc+2) * sizeof(double) );
00146
00147
00148 init_d( n, d );
00149 init_e( n, e );
00150 init_vt( ldvt*ncvt, vt );
00151 init_u( ldu*n, u );
00152 init_c( ldc*ncc, c );
00153 init_work( 4*n, work );
00154
00155
00156 for( i = 0; i < n; i++ ) {
00157 d_save[i] = d[i];
00158 }
00159 for( i = 0; i < n; i++ ) {
00160 e_save[i] = e[i];
00161 }
00162 for( i = 0; i < ldvt*ncvt; i++ ) {
00163 vt_save[i] = vt[i];
00164 }
00165 for( i = 0; i < ldu*n; i++ ) {
00166 u_save[i] = u[i];
00167 }
00168 for( i = 0; i < ldc*ncc; i++ ) {
00169 c_save[i] = c[i];
00170 }
00171
00172
00173 dbdsqr_( &uplo, &n, &ncvt, &nru, &ncc, d, e, vt, &ldvt, u, &ldu, c, &ldc,
00174 work, &info );
00175
00176
00177
00178 for( i = 0; i < n; i++ ) {
00179 d_i[i] = d_save[i];
00180 }
00181 for( i = 0; i < n; i++ ) {
00182 e_i[i] = e_save[i];
00183 }
00184 for( i = 0; i < ldvt*ncvt; i++ ) {
00185 vt_i[i] = vt_save[i];
00186 }
00187 for( i = 0; i < ldu*n; i++ ) {
00188 u_i[i] = u_save[i];
00189 }
00190 for( i = 0; i < ldc*ncc; i++ ) {
00191 c_i[i] = c_save[i];
00192 }
00193 for( i = 0; i < 4*n; i++ ) {
00194 work_i[i] = work[i];
00195 }
00196 info_i = LAPACKE_dbdsqr_work( LAPACK_COL_MAJOR, uplo_i, n_i, ncvt_i, nru_i,
00197 ncc_i, d_i, e_i, vt_i, ldvt_i, u_i, ldu_i,
00198 c_i, ldc_i, work_i );
00199
00200 failed = compare_dbdsqr( d, d_i, e, e_i, vt, vt_i, u, u_i, c, c_i, info,
00201 info_i, ldc, ldu, ldvt, n, ncc, ncvt, nru );
00202 if( failed == 0 ) {
00203 printf( "PASSED: column-major middle-level interface to dbdsqr\n" );
00204 } else {
00205 printf( "FAILED: column-major middle-level interface to dbdsqr\n" );
00206 }
00207
00208
00209
00210 for( i = 0; i < n; i++ ) {
00211 d_i[i] = d_save[i];
00212 }
00213 for( i = 0; i < n; i++ ) {
00214 e_i[i] = e_save[i];
00215 }
00216 for( i = 0; i < ldvt*ncvt; i++ ) {
00217 vt_i[i] = vt_save[i];
00218 }
00219 for( i = 0; i < ldu*n; i++ ) {
00220 u_i[i] = u_save[i];
00221 }
00222 for( i = 0; i < ldc*ncc; i++ ) {
00223 c_i[i] = c_save[i];
00224 }
00225 for( i = 0; i < 4*n; i++ ) {
00226 work_i[i] = work[i];
00227 }
00228 info_i = LAPACKE_dbdsqr( LAPACK_COL_MAJOR, uplo_i, n_i, ncvt_i, nru_i,
00229 ncc_i, d_i, e_i, vt_i, ldvt_i, u_i, ldu_i, c_i,
00230 ldc_i );
00231
00232 failed = compare_dbdsqr( d, d_i, e, e_i, vt, vt_i, u, u_i, c, c_i, info,
00233 info_i, ldc, ldu, ldvt, n, ncc, ncvt, nru );
00234 if( failed == 0 ) {
00235 printf( "PASSED: column-major high-level interface to dbdsqr\n" );
00236 } else {
00237 printf( "FAILED: column-major high-level interface to dbdsqr\n" );
00238 }
00239
00240
00241
00242 for( i = 0; i < n; i++ ) {
00243 d_i[i] = d_save[i];
00244 }
00245 for( i = 0; i < n; i++ ) {
00246 e_i[i] = e_save[i];
00247 }
00248 for( i = 0; i < ldvt*ncvt; i++ ) {
00249 vt_i[i] = vt_save[i];
00250 }
00251 for( i = 0; i < ldu*n; i++ ) {
00252 u_i[i] = u_save[i];
00253 }
00254 for( i = 0; i < ldc*ncc; i++ ) {
00255 c_i[i] = c_save[i];
00256 }
00257 for( i = 0; i < 4*n; i++ ) {
00258 work_i[i] = work[i];
00259 }
00260
00261 if( ncvt != 0 ) {
00262 LAPACKE_dge_trans( LAPACK_COL_MAJOR, n, ncvt, vt_i, ldvt, vt_r,
00263 ncvt+2 );
00264 }
00265 if( nru != 0 ) {
00266 LAPACKE_dge_trans( LAPACK_COL_MAJOR, nru, n, u_i, ldu, u_r, n+2 );
00267 }
00268 if( ncc != 0 ) {
00269 LAPACKE_dge_trans( LAPACK_COL_MAJOR, n, ncc, c_i, ldc, c_r, ncc+2 );
00270 }
00271 info_i = LAPACKE_dbdsqr_work( LAPACK_ROW_MAJOR, uplo_i, n_i, ncvt_i, nru_i,
00272 ncc_i, d_i, e_i, vt_r, ldvt_r, u_r, ldu_r,
00273 c_r, ldc_r, work_i );
00274
00275 if( ncvt != 0 ) {
00276 LAPACKE_dge_trans( LAPACK_ROW_MAJOR, n, ncvt, vt_r, ncvt+2, vt_i,
00277 ldvt );
00278 }
00279 if( nru != 0 ) {
00280 LAPACKE_dge_trans( LAPACK_ROW_MAJOR, nru, n, u_r, n+2, u_i, ldu );
00281 }
00282 if( ncc != 0 ) {
00283 LAPACKE_dge_trans( LAPACK_ROW_MAJOR, n, ncc, c_r, ncc+2, c_i, ldc );
00284 }
00285
00286 failed = compare_dbdsqr( d, d_i, e, e_i, vt, vt_i, u, u_i, c, c_i, info,
00287 info_i, ldc, ldu, ldvt, n, ncc, ncvt, nru );
00288 if( failed == 0 ) {
00289 printf( "PASSED: row-major middle-level interface to dbdsqr\n" );
00290 } else {
00291 printf( "FAILED: row-major middle-level interface to dbdsqr\n" );
00292 }
00293
00294
00295
00296 for( i = 0; i < n; i++ ) {
00297 d_i[i] = d_save[i];
00298 }
00299 for( i = 0; i < n; i++ ) {
00300 e_i[i] = e_save[i];
00301 }
00302 for( i = 0; i < ldvt*ncvt; i++ ) {
00303 vt_i[i] = vt_save[i];
00304 }
00305 for( i = 0; i < ldu*n; i++ ) {
00306 u_i[i] = u_save[i];
00307 }
00308 for( i = 0; i < ldc*ncc; i++ ) {
00309 c_i[i] = c_save[i];
00310 }
00311 for( i = 0; i < 4*n; i++ ) {
00312 work_i[i] = work[i];
00313 }
00314
00315
00316 if( ncvt != 0 ) {
00317 LAPACKE_dge_trans( LAPACK_COL_MAJOR, n, ncvt, vt_i, ldvt, vt_r,
00318 ncvt+2 );
00319 }
00320 if( nru != 0 ) {
00321 LAPACKE_dge_trans( LAPACK_COL_MAJOR, nru, n, u_i, ldu, u_r, n+2 );
00322 }
00323 if( ncc != 0 ) {
00324 LAPACKE_dge_trans( LAPACK_COL_MAJOR, n, ncc, c_i, ldc, c_r, ncc+2 );
00325 }
00326 info_i = LAPACKE_dbdsqr( LAPACK_ROW_MAJOR, uplo_i, n_i, ncvt_i, nru_i,
00327 ncc_i, d_i, e_i, vt_r, ldvt_r, u_r, ldu_r, c_r,
00328 ldc_r );
00329
00330 if( ncvt != 0 ) {
00331 LAPACKE_dge_trans( LAPACK_ROW_MAJOR, n, ncvt, vt_r, ncvt+2, vt_i,
00332 ldvt );
00333 }
00334 if( nru != 0 ) {
00335 LAPACKE_dge_trans( LAPACK_ROW_MAJOR, nru, n, u_r, n+2, u_i, ldu );
00336 }
00337 if( ncc != 0 ) {
00338 LAPACKE_dge_trans( LAPACK_ROW_MAJOR, n, ncc, c_r, ncc+2, c_i, ldc );
00339 }
00340
00341 failed = compare_dbdsqr( d, d_i, e, e_i, vt, vt_i, u, u_i, c, c_i, info,
00342 info_i, ldc, ldu, ldvt, n, ncc, ncvt, nru );
00343 if( failed == 0 ) {
00344 printf( "PASSED: row-major high-level interface to dbdsqr\n" );
00345 } else {
00346 printf( "FAILED: row-major high-level interface to dbdsqr\n" );
00347 }
00348
00349
00350 if( d != NULL ) {
00351 LAPACKE_free( d );
00352 }
00353 if( d_i != NULL ) {
00354 LAPACKE_free( d_i );
00355 }
00356 if( d_save != NULL ) {
00357 LAPACKE_free( d_save );
00358 }
00359 if( e != NULL ) {
00360 LAPACKE_free( e );
00361 }
00362 if( e_i != NULL ) {
00363 LAPACKE_free( e_i );
00364 }
00365 if( e_save != NULL ) {
00366 LAPACKE_free( e_save );
00367 }
00368 if( vt != NULL ) {
00369 LAPACKE_free( vt );
00370 }
00371 if( vt_i != NULL ) {
00372 LAPACKE_free( vt_i );
00373 }
00374 if( vt_r != NULL ) {
00375 LAPACKE_free( vt_r );
00376 }
00377 if( vt_save != NULL ) {
00378 LAPACKE_free( vt_save );
00379 }
00380 if( u != NULL ) {
00381 LAPACKE_free( u );
00382 }
00383 if( u_i != NULL ) {
00384 LAPACKE_free( u_i );
00385 }
00386 if( u_r != NULL ) {
00387 LAPACKE_free( u_r );
00388 }
00389 if( u_save != NULL ) {
00390 LAPACKE_free( u_save );
00391 }
00392 if( c != NULL ) {
00393 LAPACKE_free( c );
00394 }
00395 if( c_i != NULL ) {
00396 LAPACKE_free( c_i );
00397 }
00398 if( c_r != NULL ) {
00399 LAPACKE_free( c_r );
00400 }
00401 if( c_save != NULL ) {
00402 LAPACKE_free( c_save );
00403 }
00404 if( work != NULL ) {
00405 LAPACKE_free( work );
00406 }
00407 if( work_i != NULL ) {
00408 LAPACKE_free( work_i );
00409 }
00410
00411 return 0;
00412 }
00413
00414
00415 static void init_scalars_dbdsqr( char *uplo, lapack_int *n, lapack_int *ncvt,
00416 lapack_int *nru, lapack_int *ncc,
00417 lapack_int *ldvt, lapack_int *ldu,
00418 lapack_int *ldc )
00419 {
00420 *uplo = 'U';
00421 *n = 4;
00422 *ncvt = 4;
00423 *nru = 6;
00424 *ncc = 0;
00425 *ldvt = 8;
00426 *ldu = 8;
00427 *ldc = 8;
00428
00429 return;
00430 }
00431
00432
00433 static void init_d( lapack_int size, double *d ) {
00434 lapack_int i;
00435 for( i = 0; i < size; i++ ) {
00436 d[i] = 0;
00437 }
00438 d[0] = 3.61767881382524030e+000;
00439 d[1] = 2.41605534042332650e+000;
00440 d[2] = -1.92131096370999990e+000;
00441 d[3] = -1.42650154173085970e+000;
00442 }
00443 static void init_e( lapack_int size, double *e ) {
00444 lapack_int i;
00445 for( i = 0; i < size; i++ ) {
00446 e[i] = 0;
00447 }
00448 e[0] = 1.25871140471055650e+000;
00449 e[1] = 1.52615488994111790e+000;
00450 e[2] = -1.18946667701938940e+000;
00451 e[3] = 0.00000000000000000e+000;
00452 }
00453 static void init_vt( lapack_int size, double *vt ) {
00454 lapack_int i;
00455 for( i = 0; i < size; i++ ) {
00456 vt[i] = 0;
00457 }
00458 vt[0] = 1.00000000000000000e+000;
00459 vt[8] = 0.00000000000000000e+000;
00460 vt[16] = 0.00000000000000000e+000;
00461 vt[24] = 0.00000000000000000e+000;
00462 vt[1] = 0.00000000000000000e+000;
00463 vt[9] = -4.42198259457560860e-001;
00464 vt[17] = 6.73201611227231210e-001;
00465 vt[25] = 5.92672160619817870e-001;
00466 vt[2] = 0.00000000000000000e+000;
00467 vt[10] = -3.78772676930914820e-001;
00468 vt[18] = -7.39137581173476500e-001;
00469 vt[26] = 5.56962202763895850e-001;
00470 vt[3] = 0.00000000000000000e+000;
00471 vt[11] = 8.13014119522713070e-001;
00472 vt[19] = 2.17996958254463300e-002;
00473 vt[27] = 5.81835728293326530e-001;
00474 }
00475 static void init_u( lapack_int size, double *u ) {
00476 lapack_int i;
00477 for( i = 0; i < size; i++ ) {
00478 u[i] = 0;
00479 }
00480 u[0] = -1.57559592582321220e-001;
00481 u[8] = 2.69015046597834760e-001;
00482 u[16] = -2.61162619932311710e-001;
00483 u[24] = 8.51275109388800240e-001;
00484 u[1] = -5.33491252076982340e-001;
00485 u[9] = -5.31061312483676050e-001;
00486 u[17] = 2.92175054706849240e-001;
00487 u[25] = 1.84357485566072340e-002;
00488 u[2] = 6.35766777086559380e-001;
00489 u[10] = -3.49548734634038350e-001;
00490 u[18] = 2.49994593134254260e-002;
00491 u[26] = -2.09861574166003710e-002;
00492 u[3] = -5.33491252076982340e-001;
00493 u[11] = -3.47452857149511960e-003;
00494 u[19] = -1.53684614779524580e-001;
00495 u[27] = -2.59155252656677880e-001;
00496 u[4] = 4.14630506795582170e-002;
00497 u[12] = -5.57214452216133420e-001;
00498 u[20] = 2.91694775979778830e-001;
00499 u[28] = 4.52277601020145530e-001;
00500 u[5] = -5.52840675727442940e-003;
00501 u[13] = -4.61433040749688440e-001;
00502 u[21] = -8.58544144023679840e-001;
00503 u[29] = -5.32345671819077910e-002;
00504 }
00505 static void init_c( lapack_int size, double *c ) {
00506 lapack_int i;
00507 for( i = 0; i < size; i++ ) {
00508 c[i] = 0;
00509 }
00510 }
00511 static void init_work( lapack_int size, double *work ) {
00512 lapack_int i;
00513 for( i = 0; i < size; i++ ) {
00514 work[i] = 0;
00515 }
00516 }
00517
00518
00519
00520 static int compare_dbdsqr( double *d, double *d_i, double *e, double *e_i,
00521 double *vt, double *vt_i, double *u, double *u_i,
00522 double *c, double *c_i, lapack_int info,
00523 lapack_int info_i, lapack_int ldc, lapack_int ldu,
00524 lapack_int ldvt, lapack_int n, lapack_int ncc,
00525 lapack_int ncvt, lapack_int nru )
00526 {
00527 lapack_int i;
00528 int failed = 0;
00529 for( i = 0; i < n; i++ ) {
00530 failed += compare_doubles(d[i],d_i[i]);
00531 }
00532 for( i = 0; i < n; i++ ) {
00533 failed += compare_doubles(e[i],e_i[i]);
00534 }
00535 if( ncvt != 0 ) {
00536 for( i = 0; i < ldvt*ncvt; i++ ) {
00537 failed += compare_doubles(vt[i],vt_i[i]);
00538 }
00539 }
00540 if( nru != 0 ) {
00541 for( i = 0; i < ldu*n; i++ ) {
00542 failed += compare_doubles(u[i],u_i[i]);
00543 }
00544 }
00545 if( ncc != 0 ) {
00546 for( i = 0; i < ldc*ncc; i++ ) {
00547 failed += compare_doubles(c[i],c_i[i]);
00548 }
00549 }
00550 failed += (info == info_i) ? 0 : 1;
00551 if( info != 0 || info_i != 0 ) {
00552 printf( "info=%d, info_i=%d\n",(int)info,(int)info_i );
00553 }
00554
00555 return failed;
00556 }