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_zungtr( char *uplo, lapack_int *n, lapack_int *lda,
00055 lapack_int *lwork );
00056 static void init_a( lapack_int size, lapack_complex_double *a );
00057 static void init_tau( lapack_int size, lapack_complex_double *tau );
00058 static void init_work( lapack_int size, lapack_complex_double *work );
00059 static int compare_zungtr( lapack_complex_double *a, lapack_complex_double *a_i,
00060 lapack_int info, lapack_int info_i, lapack_int lda,
00061 lapack_int n );
00062
00063 int main(void)
00064 {
00065
00066 char uplo, uplo_i;
00067 lapack_int n, n_i;
00068 lapack_int lda, lda_i;
00069 lapack_int lda_r;
00070 lapack_int lwork, lwork_i;
00071 lapack_int info, info_i;
00072 lapack_int i;
00073 int failed;
00074
00075
00076 lapack_complex_double *a = NULL, *a_i = NULL;
00077 lapack_complex_double *tau = NULL, *tau_i = NULL;
00078 lapack_complex_double *work = NULL, *work_i = NULL;
00079 lapack_complex_double *a_save = NULL;
00080 lapack_complex_double *a_r = NULL;
00081
00082
00083 init_scalars_zungtr( &uplo, &n, &lda, &lwork );
00084 lda_r = n+2;
00085 uplo_i = uplo;
00086 n_i = n;
00087 lda_i = lda;
00088 lwork_i = lwork;
00089
00090
00091 a = (lapack_complex_double *)
00092 LAPACKE_malloc( lda*n * sizeof(lapack_complex_double) );
00093 tau = (lapack_complex_double *)
00094 LAPACKE_malloc( (n-1) * sizeof(lapack_complex_double) );
00095 work = (lapack_complex_double *)
00096 LAPACKE_malloc( lwork * sizeof(lapack_complex_double) );
00097
00098
00099 a_i = (lapack_complex_double *)
00100 LAPACKE_malloc( lda*n * sizeof(lapack_complex_double) );
00101 tau_i = (lapack_complex_double *)
00102 LAPACKE_malloc( (n-1) * sizeof(lapack_complex_double) );
00103 work_i = (lapack_complex_double *)
00104 LAPACKE_malloc( lwork * sizeof(lapack_complex_double) );
00105
00106
00107 a_save = (lapack_complex_double *)
00108 LAPACKE_malloc( lda*n * sizeof(lapack_complex_double) );
00109
00110
00111 a_r = (lapack_complex_double *)
00112 LAPACKE_malloc( n*(n+2) * sizeof(lapack_complex_double) );
00113
00114
00115 init_a( lda*n, a );
00116 init_tau( (n-1), tau );
00117 init_work( lwork, work );
00118
00119
00120 for( i = 0; i < lda*n; i++ ) {
00121 a_save[i] = a[i];
00122 }
00123
00124
00125 zungtr_( &uplo, &n, a, &lda, tau, work, &lwork, &info );
00126
00127
00128
00129 for( i = 0; i < lda*n; i++ ) {
00130 a_i[i] = a_save[i];
00131 }
00132 for( i = 0; i < (n-1); i++ ) {
00133 tau_i[i] = tau[i];
00134 }
00135 for( i = 0; i < lwork; i++ ) {
00136 work_i[i] = work[i];
00137 }
00138 info_i = LAPACKE_zungtr_work( LAPACK_COL_MAJOR, uplo_i, n_i, a_i, lda_i,
00139 tau_i, work_i, lwork_i );
00140
00141 failed = compare_zungtr( a, a_i, info, info_i, lda, n );
00142 if( failed == 0 ) {
00143 printf( "PASSED: column-major middle-level interface to zungtr\n" );
00144 } else {
00145 printf( "FAILED: column-major middle-level interface to zungtr\n" );
00146 }
00147
00148
00149
00150 for( i = 0; i < lda*n; i++ ) {
00151 a_i[i] = a_save[i];
00152 }
00153 for( i = 0; i < (n-1); i++ ) {
00154 tau_i[i] = tau[i];
00155 }
00156 for( i = 0; i < lwork; i++ ) {
00157 work_i[i] = work[i];
00158 }
00159 info_i = LAPACKE_zungtr( LAPACK_COL_MAJOR, uplo_i, n_i, a_i, lda_i, tau_i );
00160
00161 failed = compare_zungtr( a, a_i, info, info_i, lda, n );
00162 if( failed == 0 ) {
00163 printf( "PASSED: column-major high-level interface to zungtr\n" );
00164 } else {
00165 printf( "FAILED: column-major high-level interface to zungtr\n" );
00166 }
00167
00168
00169
00170 for( i = 0; i < lda*n; i++ ) {
00171 a_i[i] = a_save[i];
00172 }
00173 for( i = 0; i < (n-1); i++ ) {
00174 tau_i[i] = tau[i];
00175 }
00176 for( i = 0; i < lwork; i++ ) {
00177 work_i[i] = work[i];
00178 }
00179
00180 LAPACKE_zge_trans( LAPACK_COL_MAJOR, n, n, a_i, lda, a_r, n+2 );
00181 info_i = LAPACKE_zungtr_work( LAPACK_ROW_MAJOR, uplo_i, n_i, a_r, lda_r,
00182 tau_i, work_i, lwork_i );
00183
00184 LAPACKE_zge_trans( LAPACK_ROW_MAJOR, n, n, a_r, n+2, a_i, lda );
00185
00186 failed = compare_zungtr( a, a_i, info, info_i, lda, n );
00187 if( failed == 0 ) {
00188 printf( "PASSED: row-major middle-level interface to zungtr\n" );
00189 } else {
00190 printf( "FAILED: row-major middle-level interface to zungtr\n" );
00191 }
00192
00193
00194
00195 for( i = 0; i < lda*n; i++ ) {
00196 a_i[i] = a_save[i];
00197 }
00198 for( i = 0; i < (n-1); i++ ) {
00199 tau_i[i] = tau[i];
00200 }
00201 for( i = 0; i < lwork; i++ ) {
00202 work_i[i] = work[i];
00203 }
00204
00205
00206 LAPACKE_zge_trans( LAPACK_COL_MAJOR, n, n, a_i, lda, a_r, n+2 );
00207 info_i = LAPACKE_zungtr( LAPACK_ROW_MAJOR, uplo_i, n_i, a_r, lda_r, tau_i );
00208
00209 LAPACKE_zge_trans( LAPACK_ROW_MAJOR, n, n, a_r, n+2, a_i, lda );
00210
00211 failed = compare_zungtr( a, a_i, info, info_i, lda, n );
00212 if( failed == 0 ) {
00213 printf( "PASSED: row-major high-level interface to zungtr\n" );
00214 } else {
00215 printf( "FAILED: row-major high-level interface to zungtr\n" );
00216 }
00217
00218
00219 if( a != NULL ) {
00220 LAPACKE_free( a );
00221 }
00222 if( a_i != NULL ) {
00223 LAPACKE_free( a_i );
00224 }
00225 if( a_r != NULL ) {
00226 LAPACKE_free( a_r );
00227 }
00228 if( a_save != NULL ) {
00229 LAPACKE_free( a_save );
00230 }
00231 if( tau != NULL ) {
00232 LAPACKE_free( tau );
00233 }
00234 if( tau_i != NULL ) {
00235 LAPACKE_free( tau_i );
00236 }
00237 if( work != NULL ) {
00238 LAPACKE_free( work );
00239 }
00240 if( work_i != NULL ) {
00241 LAPACKE_free( work_i );
00242 }
00243
00244 return 0;
00245 }
00246
00247
00248 static void init_scalars_zungtr( char *uplo, lapack_int *n, lapack_int *lda,
00249 lapack_int *lwork )
00250 {
00251 *uplo = 'L';
00252 *n = 4;
00253 *lda = 8;
00254 *lwork = 512;
00255
00256 return;
00257 }
00258
00259
00260 static void init_a( lapack_int size, lapack_complex_double *a ) {
00261 lapack_int i;
00262 for( i = 0; i < size; i++ ) {
00263 a[i] = lapack_make_complex_double( 0.0, 0.0 );
00264 }
00265 a[0] = lapack_make_complex_double( -2.27999999999999980e+000,
00266 0.00000000000000000e+000 );
00267 a[8] = lapack_make_complex_double( 0.00000000000000000e+000,
00268 0.00000000000000000e+000 );
00269 a[16] = lapack_make_complex_double( 0.00000000000000000e+000,
00270 0.00000000000000000e+000 );
00271 a[24] = lapack_make_complex_double( 0.00000000000000000e+000,
00272 0.00000000000000000e+000 );
00273 a[1] = lapack_make_complex_double( -4.33845594653212970e+000,
00274 0.00000000000000000e+000 );
00275 a[9] = lapack_make_complex_double( -1.28456981649329280e-001,
00276 0.00000000000000000e+000 );
00277 a[17] = lapack_make_complex_double( 0.00000000000000000e+000,
00278 0.00000000000000000e+000 );
00279 a[25] = lapack_make_complex_double( 0.00000000000000000e+000,
00280 0.00000000000000000e+000 );
00281 a[2] = lapack_make_complex_double( 3.27860676092192380e-001,
00282 -1.25122609226443690e-001 );
00283 a[10] = lapack_make_complex_double( -2.02259457862261720e+000,
00284 0.00000000000000000e+000 );
00285 a[18] = lapack_make_complex_double( -1.66593253752407190e-001,
00286 0.00000000000000000e+000 );
00287 a[26] = lapack_make_complex_double( 0.00000000000000000e+000,
00288 0.00000000000000000e+000 );
00289 a[3] = lapack_make_complex_double( -1.41256563750694670e-001,
00290 -3.66636483973957040e-001 );
00291 a[11] = lapack_make_complex_double( -3.08321908008089010e-001,
00292 1.76322636472677850e-001 );
00293 a[19] = lapack_make_complex_double( -1.80232297833873440e+000,
00294 0.00000000000000000e+000 );
00295 a[27] = lapack_make_complex_double( -1.92494976459826360e+000,
00296 0.00000000000000000e+000 );
00297 }
00298 static void init_tau( lapack_int size, lapack_complex_double *tau ) {
00299 lapack_int i;
00300 for( i = 0; i < size; i++ ) {
00301 tau[i] = lapack_make_complex_double( 0.0, 0.0 );
00302 }
00303 tau[0] = lapack_make_complex_double( 1.41028421676675380e+000,
00304 4.67908404514893240e-001 );
00305 tau[1] = lapack_make_complex_double( 1.30242036943477490e+000,
00306 7.85332074252958030e-001 );
00307 tau[2] = lapack_make_complex_double( 1.09397371592308160e+000,
00308 -9.95574678623159850e-001 );
00309 }
00310 static void init_work( lapack_int size, lapack_complex_double *work ) {
00311 lapack_int i;
00312 for( i = 0; i < size; i++ ) {
00313 work[i] = lapack_make_complex_double( 0.0, 0.0 );
00314 }
00315 }
00316
00317
00318
00319 static int compare_zungtr( lapack_complex_double *a, lapack_complex_double *a_i,
00320 lapack_int info, lapack_int info_i, lapack_int lda,
00321 lapack_int n )
00322 {
00323 lapack_int i;
00324 int failed = 0;
00325 for( i = 0; i < lda*n; i++ ) {
00326 failed += compare_complex_doubles(a[i],a_i[i]);
00327 }
00328 failed += (info == info_i) ? 0 : 1;
00329 if( info != 0 || info_i != 0 ) {
00330 printf( "info=%d, info_i=%d\n",(int)info,(int)info_i );
00331 }
00332
00333 return failed;
00334 }