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
11void F77_sgemv(int *order, char *transp, int *m, int *n, float *alpha,
12	       float *a, int *lda, float *x, int *incx, float *beta,
13	       float *y, int *incy ) {
14
15  float *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   = ( float* )malloc( (*m)*LDA*sizeof( float ) );
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_sgemv( 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_sgemv( CblasColMajor, trans,
32		  *m, *n, *alpha, a, *lda, x, *incx, *beta, y, *incy );
33  else
34     cblas_sgemv( UNDEFINED, trans,
35		  *m, *n, *alpha, a, *lda, x, *incx, *beta, y, *incy );
36}
37
38void F77_sger(int *order, int *m, int *n, float *alpha, float *x, int *incx,
39	     float *y, int *incy, float *a, int *lda ) {
40
41  float *A;
42  int i,j,LDA;
43
44  if (*order == TEST_ROW_MJR) {
45     LDA = *n+1;
46     A   = ( float* )malloc( (*m)*LDA*sizeof( float ) );
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_sger(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_sger( CblasColMajor, *m, *n, *alpha, x, *incx, y, *incy, a, *lda );
61}
62
63void F77_strmv(int *order, char *uplow, char *transp, char *diagn,
64	      int *n, float *a, int *lda, float *x, int *incx) {
65  float *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   = ( float* )malloc( (*n)*LDA*sizeof( float ) );
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_strmv(CblasRowMajor, uplo, trans, diag, *n, A, LDA, x, *incx);
82     free(A);
83  }
84  else if (*order == TEST_COL_MJR)
85     cblas_strmv(CblasColMajor, uplo, trans, diag, *n, a, *lda, x, *incx);
86  else {
87     cblas_strmv(UNDEFINED, uplo, trans, diag, *n, a, *lda, x, *incx);
88  }
89}
90
91void F77_strsv(int *order, char *uplow, char *transp, char *diagn,
92	       int *n, float *a, int *lda, float *x, int *incx ) {
93  float *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   = ( float* )malloc( (*n)*LDA*sizeof( float ) );
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_strsv(CblasRowMajor, uplo, trans, diag, *n, A, LDA, x, *incx );
110     free(A);
111   }
112   else
113     cblas_strsv(CblasColMajor, uplo, trans, diag, *n, a, *lda, x, *incx );
114}
115void F77_ssymv(int *order, char *uplow, int *n, float *alpha, float *a,
116	      int *lda, float *x, int *incx, float *beta, float *y,
117	      int *incy) {
118  float *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   = ( float* )malloc( (*n)*LDA*sizeof( float ) );
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_ssymv(CblasRowMajor, uplo, *n, *alpha, A, LDA, x, *incx,
131		 *beta, y, *incy );
132     free(A);
133   }
134   else
135     cblas_ssymv(CblasColMajor, uplo, *n, *alpha, a, *lda, x, *incx,
136		 *beta, y, *incy );
137}
138
139void F77_ssyr(int *order, char *uplow, int *n, float *alpha, float *x,
140	     int *incx, float *a, int *lda) {
141  float *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   = ( float* )malloc( (*n)*LDA*sizeof( float ) );
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_ssyr(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_ssyr(CblasColMajor, uplo, *n, *alpha, x, *incx, a, *lda);
161}
162
163void F77_ssyr2(int *order, char *uplow, int *n, float *alpha, float *x,
164	     int *incx, float *y, int *incy, float *a, int *lda) {
165  float *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   = ( float* )malloc( (*n)*LDA*sizeof( float ) );
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_ssyr2(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_ssyr2(CblasColMajor, uplo, *n, *alpha, x, *incx, y, *incy, a, *lda);
185}
186
187void F77_sgbmv(int *order, char *transp, int *m, int *n, int *kl, int *ku,
188	       float *alpha, float *a, int *lda, float *x, int *incx,
189	       float *beta, float *y, int *incy ) {
190
191  float *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   = ( float* )malloc( (*n+*kl)*LDA*sizeof( float ) );
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_sgbmv( CblasRowMajor, trans, *m, *n, *kl, *ku, *alpha,
217		  A, LDA, x, *incx, *beta, y, *incy );
218     free(A);
219  }
220  else
221     cblas_sgbmv( CblasColMajor, trans, *m, *n, *kl, *ku, *alpha,
222		  a, *lda, x, *incx, *beta, y, *incy );
223}
224
225void F77_stbmv(int *order, char *uplow, char *transp, char *diagn,
226	      int *n, int *k, float *a, int *lda, float *x, int *incx) {
227  float *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 = ( float* )malloc( (*n+*k)*LDA*sizeof( float ) );
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_stbmv(CblasRowMajor, uplo, trans, diag, *n, *k, A, LDA, x, *incx);
265     free(A);
266   }
267   else
268     cblas_stbmv(CblasColMajor, uplo, trans, diag, *n, *k, a, *lda, x, *incx);
269}
270
271void F77_stbsv(int *order, char *uplow, char *transp, char *diagn,
272	      int *n, int *k, float *a, int *lda, float *x, int *incx) {
273  float *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 = ( float* )malloc( (*n+*k)*LDA*sizeof( float ) );
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_stbsv(CblasRowMajor, uplo, trans, diag, *n, *k, A, LDA, x, *incx);
311     free(A);
312  }
313  else
314     cblas_stbsv(CblasColMajor, uplo, trans, diag, *n, *k, a, *lda, x, *incx);
315}
316
317void F77_ssbmv(int *order, char *uplow, int *n, int *k, float *alpha,
318	      float *a, int *lda, float *x, int *incx, float *beta,
319	      float *y, int *incy) {
320  float *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   = ( float* )malloc( (*n+*k)*LDA*sizeof( float ) );
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_ssbmv(CblasRowMajor, uplo, *n, *k, *alpha, A, LDA, x, *incx,
354		 *beta, y, *incy );
355     free(A);
356   }
357   else
358     cblas_ssbmv(CblasColMajor, uplo, *n, *k, *alpha, a, *lda, x, *incx,
359		 *beta, y, *incy );
360}
361
362void F77_sspmv(int *order, char *uplow, int *n, float *alpha, float *ap,
363	      float *x, int *incx, float *beta, float *y, int *incy) {
364  float *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   = ( float* )malloc( LDA*LDA*sizeof( float ) );
373     AP  = ( float* )malloc( (((LDA+1)*LDA)/2)*sizeof( float ) );
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_sspmv( CblasRowMajor, uplo, *n, *alpha, AP, x, *incx, *beta, y,
391		  *incy );
392     free(A); free(AP);
393  }
394  else
395     cblas_sspmv( CblasColMajor, uplo, *n, *alpha, ap, x, *incx, *beta, y,
396		  *incy );
397}
398
399void F77_stpmv(int *order, char *uplow, char *transp, char *diagn,
400	      int *n, float *ap, float *x, int *incx) {
401  float *A, *AP;
402  int i, j, k, LDA;
403  enum CBLAS_TRANSPOSE trans;
404  enum CBLAS_UPLO uplo;
405  enum CBLAS_DIAG diag;
406
407  get_transpose_type(transp,&trans);
408  get_uplo_type(uplow,&uplo);
409  get_diag_type(diagn,&diag);
410
411  if (*order == TEST_ROW_MJR) {
412     LDA = *n;
413     A   = ( float* )malloc( LDA*LDA*sizeof( float ) );
414     AP  = ( float* )malloc( (((LDA+1)*LDA)/2)*sizeof( float ) );
415     if (uplo == CblasUpper) {
416        for( j=0, k=0; j<*n; j++ )
417           for( i=0; i<j+1; i++, k++ )
418              A[ LDA*i+j ]=ap[ k ];
419        for( i=0, k=0; i<*n; i++ )
420           for( j=i; j<*n; j++, k++ )
421              AP[ k ]=A[ LDA*i+j ];
422     }
423     else {
424        for( j=0, k=0; j<*n; j++ )
425           for( i=j; i<*n; i++, k++ )
426              A[ LDA*i+j ]=ap[ k ];
427        for( i=0, k=0; i<*n; i++ )
428           for( j=0; j<i+1; j++, k++ )
429              AP[ k ]=A[ LDA*i+j ];
430     }
431     cblas_stpmv( CblasRowMajor, uplo, trans, diag, *n, AP, x, *incx );
432     free(A); free(AP);
433  }
434  else
435     cblas_stpmv( CblasColMajor, uplo, trans, diag, *n, ap, x, *incx );
436}
437
438void F77_stpsv(int *order, char *uplow, char *transp, char *diagn,
439	      int *n, float *ap, float *x, int *incx) {
440  float *A, *AP;
441  int i, j, k, LDA;
442  enum CBLAS_TRANSPOSE trans;
443  enum CBLAS_UPLO uplo;
444  enum CBLAS_DIAG diag;
445
446  get_transpose_type(transp,&trans);
447  get_uplo_type(uplow,&uplo);
448  get_diag_type(diagn,&diag);
449
450  if (*order == TEST_ROW_MJR) {
451     LDA = *n;
452     A   = ( float* )malloc( LDA*LDA*sizeof( float ) );
453     AP  = ( float* )malloc( (((LDA+1)*LDA)/2)*sizeof( float ) );
454     if (uplo == CblasUpper) {
455        for( j=0, k=0; j<*n; j++ )
456           for( i=0; i<j+1; i++, k++ )
457              A[ LDA*i+j ]=ap[ k ];
458        for( i=0, k=0; i<*n; i++ )
459           for( j=i; j<*n; j++, k++ )
460              AP[ k ]=A[ LDA*i+j ];
461
462     }
463     else {
464        for( j=0, k=0; j<*n; j++ )
465           for( i=j; i<*n; i++, k++ )
466              A[ LDA*i+j ]=ap[ k ];
467        for( i=0, k=0; i<*n; i++ )
468           for( j=0; j<i+1; j++, k++ )
469              AP[ k ]=A[ LDA*i+j ];
470     }
471     cblas_stpsv( CblasRowMajor, uplo, trans, diag, *n, AP, x, *incx );
472     free(A); free(AP);
473  }
474  else
475     cblas_stpsv( CblasColMajor, uplo, trans, diag, *n, ap, x, *incx );
476}
477
478void F77_sspr(int *order, char *uplow, int *n, float *alpha, float *x,
479	     int *incx, float *ap ){
480  float *A, *AP;
481  int i,j,k,LDA;
482  enum CBLAS_UPLO uplo;
483
484  get_uplo_type(uplow,&uplo);
485
486  if (*order == TEST_ROW_MJR) {
487     LDA = *n;
488     A   = ( float* )malloc( LDA*LDA*sizeof( float ) );
489     AP  = ( float* )malloc( (((LDA+1)*LDA)/2)*sizeof( float ) );
490     if (uplo == CblasUpper) {
491        for( j=0, k=0; j<*n; j++ )
492           for( i=0; i<j+1; i++, k++ )
493              A[ LDA*i+j ]=ap[ k ];
494        for( i=0, k=0; i<*n; i++ )
495           for( j=i; j<*n; j++, k++ )
496              AP[ k ]=A[ LDA*i+j ];
497     }
498     else {
499        for( j=0, k=0; j<*n; j++ )
500           for( i=j; i<*n; i++, k++ )
501              A[ LDA*i+j ]=ap[ k ];
502        for( i=0, k=0; i<*n; i++ )
503           for( j=0; j<i+1; j++, k++ )
504              AP[ k ]=A[ LDA*i+j ];
505     }
506     cblas_sspr( CblasRowMajor, uplo, *n, *alpha, x, *incx, AP );
507     if (uplo == CblasUpper) {
508        for( i=0, k=0; i<*n; i++ )
509           for( j=i; j<*n; j++, k++ )
510              A[ LDA*i+j ]=AP[ k ];
511        for( j=0, k=0; j<*n; j++ )
512           for( i=0; i<j+1; i++, k++ )
513              ap[ k ]=A[ LDA*i+j ];
514     }
515     else {
516        for( i=0, k=0; i<*n; i++ )
517           for( j=0; j<i+1; j++, k++ )
518              A[ LDA*i+j ]=AP[ k ];
519        for( j=0, k=0; j<*n; j++ )
520           for( i=j; i<*n; i++, k++ )
521              ap[ k ]=A[ LDA*i+j ];
522     }
523     free(A); free(AP);
524  }
525  else
526     cblas_sspr( CblasColMajor, uplo, *n, *alpha, x, *incx, ap );
527}
528
529void F77_sspr2(int *order, char *uplow, int *n, float *alpha, float *x,
530	     int *incx, float *y, int *incy, float *ap ){
531  float *A, *AP;
532  int i,j,k,LDA;
533  enum CBLAS_UPLO uplo;
534
535  get_uplo_type(uplow,&uplo);
536
537  if (*order == TEST_ROW_MJR) {
538     LDA = *n;
539     A   = ( float* )malloc( LDA*LDA*sizeof( float ) );
540     AP  = ( float* )malloc( (((LDA+1)*LDA)/2)*sizeof( float ) );
541     if (uplo == CblasUpper) {
542        for( j=0, k=0; j<*n; j++ )
543           for( i=0; i<j+1; i++, k++ )
544              A[ LDA*i+j ]=ap[ k ];
545        for( i=0, k=0; i<*n; i++ )
546           for( j=i; j<*n; j++, k++ )
547              AP[ k ]=A[ LDA*i+j ];
548     }
549     else {
550        for( j=0, k=0; j<*n; j++ )
551           for( i=j; i<*n; i++, k++ )
552              A[ LDA*i+j ]=ap[ k ];
553        for( i=0, k=0; i<*n; i++ )
554           for( j=0; j<i+1; j++, k++ )
555              AP[ k ]=A[ LDA*i+j ];
556     }
557     cblas_sspr2( CblasRowMajor, uplo, *n, *alpha, x, *incx, y, *incy, AP );
558     if (uplo == CblasUpper) {
559        for( i=0, k=0; i<*n; i++ )
560           for( j=i; j<*n; j++, k++ )
561              A[ LDA*i+j ]=AP[ k ];
562        for( j=0, k=0; j<*n; j++ )
563           for( i=0; i<j+1; i++, k++ )
564              ap[ k ]=A[ LDA*i+j ];
565     }
566     else {
567        for( i=0, k=0; i<*n; i++ )
568           for( j=0; j<i+1; j++, k++ )
569              A[ LDA*i+j ]=AP[ k ];
570        for( j=0, k=0; j<*n; j++ )
571           for( i=j; i<*n; i++, k++ )
572              ap[ k ]=A[ LDA*i+j ];
573     }
574     free(A);
575     free(AP);
576  }
577  else
578     cblas_sspr2( CblasColMajor, uplo, *n, *alpha, x, *incx, y, *incy, ap );
579}
580