• Home
  • Line#
  • Scopes#
  • Navigate#
  • Raw
  • Download
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