1 /* 2 * Written by D.P. Manley, Digital Equipment Corporation. 3 * Prefixed "C_" to BLAS routines and their declarations. 4 * 5 * Modified by T. H. Do, 1/23/98, SGI/CRAY Research. 6 */ 7 #include <stdlib.h> 8 #include "cblas.h" 9 #include "cblas_test.h" 10 F77_dgemv(int * order,char * transp,int * m,int * n,double * alpha,double * a,int * lda,double * x,int * incx,double * beta,double * y,int * incy)11 void F77_dgemv(int *order, char *transp, int *m, int *n, double *alpha, 12 double *a, int *lda, double *x, int *incx, double *beta, 13 double *y, int *incy ) { 14 15 double *A; 16 int i,j,LDA; 17 enum CBLAS_TRANSPOSE trans; 18 19 get_transpose_type(transp, &trans); 20 if (*order == TEST_ROW_MJR) { 21 LDA = *n+1; 22 A = ( double* )malloc( (*m)*LDA*sizeof( double ) ); 23 for( i=0; i<*m; i++ ) 24 for( j=0; j<*n; j++ ) 25 A[ LDA*i+j ]=a[ (*lda)*j+i ]; 26 cblas_dgemv( CblasRowMajor, trans, 27 *m, *n, *alpha, A, LDA, x, *incx, *beta, y, *incy ); 28 free(A); 29 } 30 else if (*order == TEST_COL_MJR) 31 cblas_dgemv( CblasColMajor, trans, 32 *m, *n, *alpha, a, *lda, x, *incx, *beta, y, *incy ); 33 else 34 cblas_dgemv( UNDEFINED, trans, 35 *m, *n, *alpha, a, *lda, x, *incx, *beta, y, *incy ); 36 } 37 F77_dger(int * order,int * m,int * n,double * alpha,double * x,int * incx,double * y,int * incy,double * a,int * lda)38 void F77_dger(int *order, int *m, int *n, double *alpha, double *x, int *incx, 39 double *y, int *incy, double *a, int *lda ) { 40 41 double *A; 42 int i,j,LDA; 43 44 if (*order == TEST_ROW_MJR) { 45 LDA = *n+1; 46 A = ( double* )malloc( (*m)*LDA*sizeof( double ) ); 47 48 for( i=0; i<*m; i++ ) { 49 for( j=0; j<*n; j++ ) 50 A[ LDA*i+j ]=a[ (*lda)*j+i ]; 51 } 52 53 cblas_dger(CblasRowMajor, *m, *n, *alpha, x, *incx, y, *incy, A, LDA ); 54 for( i=0; i<*m; i++ ) 55 for( j=0; j<*n; j++ ) 56 a[ (*lda)*j+i ]=A[ LDA*i+j ]; 57 free(A); 58 } 59 else 60 cblas_dger( CblasColMajor, *m, *n, *alpha, x, *incx, y, *incy, a, *lda ); 61 } 62 F77_dtrmv(int * order,char * uplow,char * transp,char * diagn,int * n,double * a,int * lda,double * x,int * incx)63 void F77_dtrmv(int *order, char *uplow, char *transp, char *diagn, 64 int *n, double *a, int *lda, double *x, int *incx) { 65 double *A; 66 int i,j,LDA; 67 enum CBLAS_TRANSPOSE trans; 68 enum CBLAS_UPLO uplo; 69 enum CBLAS_DIAG diag; 70 71 get_transpose_type(transp,&trans); 72 get_uplo_type(uplow,&uplo); 73 get_diag_type(diagn,&diag); 74 75 if (*order == TEST_ROW_MJR) { 76 LDA = *n+1; 77 A = ( double* )malloc( (*n)*LDA*sizeof( double ) ); 78 for( i=0; i<*n; i++ ) 79 for( j=0; j<*n; j++ ) 80 A[ LDA*i+j ]=a[ (*lda)*j+i ]; 81 cblas_dtrmv(CblasRowMajor, uplo, trans, diag, *n, A, LDA, x, *incx); 82 free(A); 83 } 84 else if (*order == TEST_COL_MJR) 85 cblas_dtrmv(CblasColMajor, uplo, trans, diag, *n, a, *lda, x, *incx); 86 else { 87 cblas_dtrmv(UNDEFINED, uplo, trans, diag, *n, a, *lda, x, *incx); 88 } 89 } 90 F77_dtrsv(int * order,char * uplow,char * transp,char * diagn,int * n,double * a,int * lda,double * x,int * incx)91 void F77_dtrsv(int *order, char *uplow, char *transp, char *diagn, 92 int *n, double *a, int *lda, double *x, int *incx ) { 93 double *A; 94 int i,j,LDA; 95 enum CBLAS_TRANSPOSE trans; 96 enum CBLAS_UPLO uplo; 97 enum CBLAS_DIAG diag; 98 99 get_transpose_type(transp,&trans); 100 get_uplo_type(uplow,&uplo); 101 get_diag_type(diagn,&diag); 102 103 if (*order == TEST_ROW_MJR) { 104 LDA = *n+1; 105 A = ( double* )malloc( (*n)*LDA*sizeof( double ) ); 106 for( i=0; i<*n; i++ ) 107 for( j=0; j<*n; j++ ) 108 A[ LDA*i+j ]=a[ (*lda)*j+i ]; 109 cblas_dtrsv(CblasRowMajor, uplo, trans, diag, *n, A, LDA, x, *incx ); 110 free(A); 111 } 112 else 113 cblas_dtrsv(CblasColMajor, uplo, trans, diag, *n, a, *lda, x, *incx ); 114 } F77_dsymv(int * order,char * uplow,int * n,double * alpha,double * a,int * lda,double * x,int * incx,double * beta,double * y,int * incy)115 void F77_dsymv(int *order, char *uplow, int *n, double *alpha, double *a, 116 int *lda, double *x, int *incx, double *beta, double *y, 117 int *incy) { 118 double *A; 119 int i,j,LDA; 120 enum CBLAS_UPLO uplo; 121 122 get_uplo_type(uplow,&uplo); 123 124 if (*order == TEST_ROW_MJR) { 125 LDA = *n+1; 126 A = ( double* )malloc( (*n)*LDA*sizeof( double ) ); 127 for( i=0; i<*n; i++ ) 128 for( j=0; j<*n; j++ ) 129 A[ LDA*i+j ]=a[ (*lda)*j+i ]; 130 cblas_dsymv(CblasRowMajor, uplo, *n, *alpha, A, LDA, x, *incx, 131 *beta, y, *incy ); 132 free(A); 133 } 134 else 135 cblas_dsymv(CblasColMajor, uplo, *n, *alpha, a, *lda, x, *incx, 136 *beta, y, *incy ); 137 } 138 F77_dsyr(int * order,char * uplow,int * n,double * alpha,double * x,int * incx,double * a,int * lda)139 void F77_dsyr(int *order, char *uplow, int *n, double *alpha, double *x, 140 int *incx, double *a, int *lda) { 141 double *A; 142 int i,j,LDA; 143 enum CBLAS_UPLO uplo; 144 145 get_uplo_type(uplow,&uplo); 146 147 if (*order == TEST_ROW_MJR) { 148 LDA = *n+1; 149 A = ( double* )malloc( (*n)*LDA*sizeof( double ) ); 150 for( i=0; i<*n; i++ ) 151 for( j=0; j<*n; j++ ) 152 A[ LDA*i+j ]=a[ (*lda)*j+i ]; 153 cblas_dsyr(CblasRowMajor, uplo, *n, *alpha, x, *incx, A, LDA); 154 for( i=0; i<*n; i++ ) 155 for( j=0; j<*n; j++ ) 156 a[ (*lda)*j+i ]=A[ LDA*i+j ]; 157 free(A); 158 } 159 else 160 cblas_dsyr(CblasColMajor, uplo, *n, *alpha, x, *incx, a, *lda); 161 } 162 F77_dsyr2(int * order,char * uplow,int * n,double * alpha,double * x,int * incx,double * y,int * incy,double * a,int * lda)163 void F77_dsyr2(int *order, char *uplow, int *n, double *alpha, double *x, 164 int *incx, double *y, int *incy, double *a, int *lda) { 165 double *A; 166 int i,j,LDA; 167 enum CBLAS_UPLO uplo; 168 169 get_uplo_type(uplow,&uplo); 170 171 if (*order == TEST_ROW_MJR) { 172 LDA = *n+1; 173 A = ( double* )malloc( (*n)*LDA*sizeof( double ) ); 174 for( i=0; i<*n; i++ ) 175 for( j=0; j<*n; j++ ) 176 A[ LDA*i+j ]=a[ (*lda)*j+i ]; 177 cblas_dsyr2(CblasRowMajor, uplo, *n, *alpha, x, *incx, y, *incy, A, LDA); 178 for( i=0; i<*n; i++ ) 179 for( j=0; j<*n; j++ ) 180 a[ (*lda)*j+i ]=A[ LDA*i+j ]; 181 free(A); 182 } 183 else 184 cblas_dsyr2(CblasColMajor, uplo, *n, *alpha, x, *incx, y, *incy, a, *lda); 185 } 186 F77_dgbmv(int * order,char * transp,int * m,int * n,int * kl,int * ku,double * alpha,double * a,int * lda,double * x,int * incx,double * beta,double * y,int * incy)187 void F77_dgbmv(int *order, char *transp, int *m, int *n, int *kl, int *ku, 188 double *alpha, double *a, int *lda, double *x, int *incx, 189 double *beta, double *y, int *incy ) { 190 191 double *A; 192 int i,irow,j,jcol,LDA; 193 enum CBLAS_TRANSPOSE trans; 194 195 get_transpose_type(transp, &trans); 196 197 if (*order == TEST_ROW_MJR) { 198 LDA = *ku+*kl+2; 199 A = ( double* )malloc( (*n+*kl)*LDA*sizeof( double ) ); 200 for( i=0; i<*ku; i++ ){ 201 irow=*ku+*kl-i; 202 jcol=(*ku)-i; 203 for( j=jcol; j<*n; j++ ) 204 A[ LDA*(j-jcol)+irow ]=a[ (*lda)*j+i ]; 205 } 206 i=*ku; 207 irow=*ku+*kl-i; 208 for( j=0; j<*n; j++ ) 209 A[ LDA*j+irow ]=a[ (*lda)*j+i ]; 210 for( i=*ku+1; i<*ku+*kl+1; i++ ){ 211 irow=*ku+*kl-i; 212 jcol=i-(*ku); 213 for( j=jcol; j<(*n+*kl); j++ ) 214 A[ LDA*j+irow ]=a[ (*lda)*(j-jcol)+i ]; 215 } 216 cblas_dgbmv( CblasRowMajor, trans, *m, *n, *kl, *ku, *alpha, 217 A, LDA, x, *incx, *beta, y, *incy ); 218 free(A); 219 } 220 else 221 cblas_dgbmv( CblasColMajor, trans, *m, *n, *kl, *ku, *alpha, 222 a, *lda, x, *incx, *beta, y, *incy ); 223 } 224 F77_dtbmv(int * order,char * uplow,char * transp,char * diagn,int * n,int * k,double * a,int * lda,double * x,int * incx)225 void F77_dtbmv(int *order, char *uplow, char *transp, char *diagn, 226 int *n, int *k, double *a, int *lda, double *x, int *incx) { 227 double *A; 228 int irow, jcol, i, j, LDA; 229 enum CBLAS_TRANSPOSE trans; 230 enum CBLAS_UPLO uplo; 231 enum CBLAS_DIAG diag; 232 233 get_transpose_type(transp,&trans); 234 get_uplo_type(uplow,&uplo); 235 get_diag_type(diagn,&diag); 236 237 if (*order == TEST_ROW_MJR) { 238 LDA = *k+1; 239 A = ( double* )malloc( (*n+*k)*LDA*sizeof( double ) ); 240 if (uplo == CblasUpper) { 241 for( i=0; i<*k; i++ ){ 242 irow=*k-i; 243 jcol=(*k)-i; 244 for( j=jcol; j<*n; j++ ) 245 A[ LDA*(j-jcol)+irow ]=a[ (*lda)*j+i ]; 246 } 247 i=*k; 248 irow=*k-i; 249 for( j=0; j<*n; j++ ) 250 A[ LDA*j+irow ]=a[ (*lda)*j+i ]; 251 } 252 else { 253 i=0; 254 irow=*k-i; 255 for( j=0; j<*n; j++ ) 256 A[ LDA*j+irow ]=a[ (*lda)*j+i ]; 257 for( i=1; i<*k+1; i++ ){ 258 irow=*k-i; 259 jcol=i; 260 for( j=jcol; j<(*n+*k); j++ ) 261 A[ LDA*j+irow ]=a[ (*lda)*(j-jcol)+i ]; 262 } 263 } 264 cblas_dtbmv(CblasRowMajor, uplo, trans, diag, *n, *k, A, LDA, x, *incx); 265 free(A); 266 } 267 else 268 cblas_dtbmv(CblasColMajor, uplo, trans, diag, *n, *k, a, *lda, x, *incx); 269 } 270 F77_dtbsv(int * order,char * uplow,char * transp,char * diagn,int * n,int * k,double * a,int * lda,double * x,int * incx)271 void F77_dtbsv(int *order, char *uplow, char *transp, char *diagn, 272 int *n, int *k, double *a, int *lda, double *x, int *incx) { 273 double *A; 274 int irow, jcol, i, j, LDA; 275 enum CBLAS_TRANSPOSE trans; 276 enum CBLAS_UPLO uplo; 277 enum CBLAS_DIAG diag; 278 279 get_transpose_type(transp,&trans); 280 get_uplo_type(uplow,&uplo); 281 get_diag_type(diagn,&diag); 282 283 if (*order == TEST_ROW_MJR) { 284 LDA = *k+1; 285 A = ( double* )malloc( (*n+*k)*LDA*sizeof( double ) ); 286 if (uplo == CblasUpper) { 287 for( i=0; i<*k; i++ ){ 288 irow=*k-i; 289 jcol=(*k)-i; 290 for( j=jcol; j<*n; j++ ) 291 A[ LDA*(j-jcol)+irow ]=a[ (*lda)*j+i ]; 292 } 293 i=*k; 294 irow=*k-i; 295 for( j=0; j<*n; j++ ) 296 A[ LDA*j+irow ]=a[ (*lda)*j+i ]; 297 } 298 else { 299 i=0; 300 irow=*k-i; 301 for( j=0; j<*n; j++ ) 302 A[ LDA*j+irow ]=a[ (*lda)*j+i ]; 303 for( i=1; i<*k+1; i++ ){ 304 irow=*k-i; 305 jcol=i; 306 for( j=jcol; j<(*n+*k); j++ ) 307 A[ LDA*j+irow ]=a[ (*lda)*(j-jcol)+i ]; 308 } 309 } 310 cblas_dtbsv(CblasRowMajor, uplo, trans, diag, *n, *k, A, LDA, x, *incx); 311 free(A); 312 } 313 else 314 cblas_dtbsv(CblasColMajor, uplo, trans, diag, *n, *k, a, *lda, x, *incx); 315 } 316 F77_dsbmv(int * order,char * uplow,int * n,int * k,double * alpha,double * a,int * lda,double * x,int * incx,double * beta,double * y,int * incy)317 void F77_dsbmv(int *order, char *uplow, int *n, int *k, double *alpha, 318 double *a, int *lda, double *x, int *incx, double *beta, 319 double *y, int *incy) { 320 double *A; 321 int i,j,irow,jcol,LDA; 322 enum CBLAS_UPLO uplo; 323 324 get_uplo_type(uplow,&uplo); 325 326 if (*order == TEST_ROW_MJR) { 327 LDA = *k+1; 328 A = ( double* )malloc( (*n+*k)*LDA*sizeof( double ) ); 329 if (uplo == CblasUpper) { 330 for( i=0; i<*k; i++ ){ 331 irow=*k-i; 332 jcol=(*k)-i; 333 for( j=jcol; j<*n; j++ ) 334 A[ LDA*(j-jcol)+irow ]=a[ (*lda)*j+i ]; 335 } 336 i=*k; 337 irow=*k-i; 338 for( j=0; j<*n; j++ ) 339 A[ LDA*j+irow ]=a[ (*lda)*j+i ]; 340 } 341 else { 342 i=0; 343 irow=*k-i; 344 for( j=0; j<*n; j++ ) 345 A[ LDA*j+irow ]=a[ (*lda)*j+i ]; 346 for( i=1; i<*k+1; i++ ){ 347 irow=*k-i; 348 jcol=i; 349 for( j=jcol; j<(*n+*k); j++ ) 350 A[ LDA*j+irow ]=a[ (*lda)*(j-jcol)+i ]; 351 } 352 } 353 cblas_dsbmv(CblasRowMajor, uplo, *n, *k, *alpha, A, LDA, x, *incx, 354 *beta, y, *incy ); 355 free(A); 356 } 357 else 358 cblas_dsbmv(CblasColMajor, uplo, *n, *k, *alpha, a, *lda, x, *incx, 359 *beta, y, *incy ); 360 } 361 F77_dspmv(int * order,char * uplow,int * n,double * alpha,double * ap,double * x,int * incx,double * beta,double * y,int * incy)362 void F77_dspmv(int *order, char *uplow, int *n, double *alpha, double *ap, 363 double *x, int *incx, double *beta, double *y, int *incy) { 364 double *A,*AP; 365 int i,j,k,LDA; 366 enum CBLAS_UPLO uplo; 367 368 get_uplo_type(uplow,&uplo); 369 370 if (*order == TEST_ROW_MJR) { 371 LDA = *n; 372 A = ( double* )malloc( LDA*LDA*sizeof( double ) ); 373 AP = ( double* )malloc( (((LDA+1)*LDA)/2)*sizeof( double ) ); 374 if (uplo == CblasUpper) { 375 for( j=0, k=0; j<*n; j++ ) 376 for( i=0; i<j+1; i++, k++ ) 377 A[ LDA*i+j ]=ap[ k ]; 378 for( i=0, k=0; i<*n; i++ ) 379 for( j=i; j<*n; j++, k++ ) 380 AP[ k ]=A[ LDA*i+j ]; 381 } 382 else { 383 for( j=0, k=0; j<*n; j++ ) 384 for( i=j; i<*n; i++, k++ ) 385 A[ LDA*i+j ]=ap[ k ]; 386 for( i=0, k=0; i<*n; i++ ) 387 for( j=0; j<i+1; j++, k++ ) 388 AP[ k ]=A[ LDA*i+j ]; 389 } 390 cblas_dspmv( CblasRowMajor, uplo, *n, *alpha, AP, x, *incx, *beta, y, 391 *incy ); 392 free(A); 393 free(AP); 394 } 395 else 396 cblas_dspmv( CblasColMajor, uplo, *n, *alpha, ap, x, *incx, *beta, y, 397 *incy ); 398 } 399 F77_dtpmv(int * order,char * uplow,char * transp,char * diagn,int * n,double * ap,double * x,int * incx)400 void F77_dtpmv(int *order, char *uplow, char *transp, char *diagn, 401 int *n, double *ap, double *x, int *incx) { 402 double *A, *AP; 403 int i, j, k, LDA; 404 enum CBLAS_TRANSPOSE trans; 405 enum CBLAS_UPLO uplo; 406 enum CBLAS_DIAG diag; 407 408 get_transpose_type(transp,&trans); 409 get_uplo_type(uplow,&uplo); 410 get_diag_type(diagn,&diag); 411 412 if (*order == TEST_ROW_MJR) { 413 LDA = *n; 414 A = ( double* )malloc( LDA*LDA*sizeof( double ) ); 415 AP = ( double* )malloc( (((LDA+1)*LDA)/2)*sizeof( double ) ); 416 if (uplo == CblasUpper) { 417 for( j=0, k=0; j<*n; j++ ) 418 for( i=0; i<j+1; i++, k++ ) 419 A[ LDA*i+j ]=ap[ k ]; 420 for( i=0, k=0; i<*n; i++ ) 421 for( j=i; j<*n; j++, k++ ) 422 AP[ k ]=A[ LDA*i+j ]; 423 } 424 else { 425 for( j=0, k=0; j<*n; j++ ) 426 for( i=j; i<*n; i++, k++ ) 427 A[ LDA*i+j ]=ap[ k ]; 428 for( i=0, k=0; i<*n; i++ ) 429 for( j=0; j<i+1; j++, k++ ) 430 AP[ k ]=A[ LDA*i+j ]; 431 } 432 cblas_dtpmv( CblasRowMajor, uplo, trans, diag, *n, AP, x, *incx ); 433 free(A); 434 free(AP); 435 } 436 else 437 cblas_dtpmv( CblasColMajor, uplo, trans, diag, *n, ap, x, *incx ); 438 } 439 F77_dtpsv(int * order,char * uplow,char * transp,char * diagn,int * n,double * ap,double * x,int * incx)440 void F77_dtpsv(int *order, char *uplow, char *transp, char *diagn, 441 int *n, double *ap, double *x, int *incx) { 442 double *A, *AP; 443 int i, j, k, LDA; 444 enum CBLAS_TRANSPOSE trans; 445 enum CBLAS_UPLO uplo; 446 enum CBLAS_DIAG diag; 447 448 get_transpose_type(transp,&trans); 449 get_uplo_type(uplow,&uplo); 450 get_diag_type(diagn,&diag); 451 452 if (*order == TEST_ROW_MJR) { 453 LDA = *n; 454 A = ( double* )malloc( LDA*LDA*sizeof( double ) ); 455 AP = ( double* )malloc( (((LDA+1)*LDA)/2)*sizeof( double ) ); 456 if (uplo == CblasUpper) { 457 for( j=0, k=0; j<*n; j++ ) 458 for( i=0; i<j+1; i++, k++ ) 459 A[ LDA*i+j ]=ap[ k ]; 460 for( i=0, k=0; i<*n; i++ ) 461 for( j=i; j<*n; j++, k++ ) 462 AP[ k ]=A[ LDA*i+j ]; 463 464 } 465 else { 466 for( j=0, k=0; j<*n; j++ ) 467 for( i=j; i<*n; i++, k++ ) 468 A[ LDA*i+j ]=ap[ k ]; 469 for( i=0, k=0; i<*n; i++ ) 470 for( j=0; j<i+1; j++, k++ ) 471 AP[ k ]=A[ LDA*i+j ]; 472 } 473 cblas_dtpsv( CblasRowMajor, uplo, trans, diag, *n, AP, x, *incx ); 474 free(A); 475 free(AP); 476 } 477 else 478 cblas_dtpsv( CblasColMajor, uplo, trans, diag, *n, ap, x, *incx ); 479 } 480 F77_dspr(int * order,char * uplow,int * n,double * alpha,double * x,int * incx,double * ap)481 void F77_dspr(int *order, char *uplow, int *n, double *alpha, double *x, 482 int *incx, double *ap ){ 483 double *A, *AP; 484 int i,j,k,LDA; 485 enum CBLAS_UPLO uplo; 486 487 get_uplo_type(uplow,&uplo); 488 489 if (*order == TEST_ROW_MJR) { 490 LDA = *n; 491 A = ( double* )malloc( LDA*LDA*sizeof( double ) ); 492 AP = ( double* )malloc( (((LDA+1)*LDA)/2)*sizeof( double ) ); 493 if (uplo == CblasUpper) { 494 for( j=0, k=0; j<*n; j++ ) 495 for( i=0; i<j+1; i++, k++ ) 496 A[ LDA*i+j ]=ap[ k ]; 497 for( i=0, k=0; i<*n; i++ ) 498 for( j=i; j<*n; j++, k++ ) 499 AP[ k ]=A[ LDA*i+j ]; 500 } 501 else { 502 for( j=0, k=0; j<*n; j++ ) 503 for( i=j; i<*n; i++, k++ ) 504 A[ LDA*i+j ]=ap[ k ]; 505 for( i=0, k=0; i<*n; i++ ) 506 for( j=0; j<i+1; j++, k++ ) 507 AP[ k ]=A[ LDA*i+j ]; 508 } 509 cblas_dspr( CblasRowMajor, uplo, *n, *alpha, x, *incx, AP ); 510 if (uplo == CblasUpper) { 511 for( i=0, k=0; i<*n; i++ ) 512 for( j=i; j<*n; j++, k++ ) 513 A[ LDA*i+j ]=AP[ k ]; 514 for( j=0, k=0; j<*n; j++ ) 515 for( i=0; i<j+1; i++, k++ ) 516 ap[ k ]=A[ LDA*i+j ]; 517 } 518 else { 519 for( i=0, k=0; i<*n; i++ ) 520 for( j=0; j<i+1; j++, k++ ) 521 A[ LDA*i+j ]=AP[ k ]; 522 for( j=0, k=0; j<*n; j++ ) 523 for( i=j; i<*n; i++, k++ ) 524 ap[ k ]=A[ LDA*i+j ]; 525 } 526 free(A); 527 free(AP); 528 } 529 else 530 cblas_dspr( CblasColMajor, uplo, *n, *alpha, x, *incx, ap ); 531 } 532 F77_dspr2(int * order,char * uplow,int * n,double * alpha,double * x,int * incx,double * y,int * incy,double * ap)533 void F77_dspr2(int *order, char *uplow, int *n, double *alpha, double *x, 534 int *incx, double *y, int *incy, double *ap ){ 535 double *A, *AP; 536 int i,j,k,LDA; 537 enum CBLAS_UPLO uplo; 538 539 get_uplo_type(uplow,&uplo); 540 541 if (*order == TEST_ROW_MJR) { 542 LDA = *n; 543 A = ( double* )malloc( LDA*LDA*sizeof( double ) ); 544 AP = ( double* )malloc( (((LDA+1)*LDA)/2)*sizeof( double ) ); 545 if (uplo == CblasUpper) { 546 for( j=0, k=0; j<*n; j++ ) 547 for( i=0; i<j+1; i++, k++ ) 548 A[ LDA*i+j ]=ap[ k ]; 549 for( i=0, k=0; i<*n; i++ ) 550 for( j=i; j<*n; j++, k++ ) 551 AP[ k ]=A[ LDA*i+j ]; 552 } 553 else { 554 for( j=0, k=0; j<*n; j++ ) 555 for( i=j; i<*n; i++, k++ ) 556 A[ LDA*i+j ]=ap[ k ]; 557 for( i=0, k=0; i<*n; i++ ) 558 for( j=0; j<i+1; j++, k++ ) 559 AP[ k ]=A[ LDA*i+j ]; 560 } 561 cblas_dspr2( CblasRowMajor, uplo, *n, *alpha, x, *incx, y, *incy, AP ); 562 if (uplo == CblasUpper) { 563 for( i=0, k=0; i<*n; i++ ) 564 for( j=i; j<*n; j++, k++ ) 565 A[ LDA*i+j ]=AP[ k ]; 566 for( j=0, k=0; j<*n; j++ ) 567 for( i=0; i<j+1; i++, k++ ) 568 ap[ k ]=A[ LDA*i+j ]; 569 } 570 else { 571 for( i=0, k=0; i<*n; i++ ) 572 for( j=0; j<i+1; j++, k++ ) 573 A[ LDA*i+j ]=AP[ k ]; 574 for( j=0, k=0; j<*n; j++ ) 575 for( i=j; i<*n; i++, k++ ) 576 ap[ k ]=A[ LDA*i+j ]; 577 } 578 free(A); 579 free(AP); 580 } 581 else 582 cblas_dspr2( CblasColMajor, uplo, *n, *alpha, x, *incx, y, *incy, ap ); 583 } 584