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, 2/19/98, SGI/CRAY Research.
6 */
7#include <stdio.h>
8#include <stdlib.h>
9#include "cblas.h"
10#include "cblas_test.h"
11
12void F77_sgemm(int *order, char *transpa, char *transpb, int *m, int *n,
13              int *k, float *alpha, float *a, int *lda, float *b, int *ldb,
14              float *beta, float *c, int *ldc ) {
15
16  float *A, *B, *C;
17  int i,j,LDA, LDB, LDC;
18  enum CBLAS_TRANSPOSE transa, transb;
19
20  get_transpose_type(transpa, &transa);
21  get_transpose_type(transpb, &transb);
22
23  if (*order == TEST_ROW_MJR) {
24     if (transa == CblasNoTrans) {
25        LDA = *k+1;
26        A = (float *)malloc( (*m)*LDA*sizeof( float ) );
27        for( i=0; i<*m; i++ )
28           for( j=0; j<*k; j++ )
29              A[i*LDA+j]=a[j*(*lda)+i];
30     }
31     else {
32        LDA = *m+1;
33        A   = ( float* )malloc( LDA*(*k)*sizeof( float ) );
34        for( i=0; i<*k; i++ )
35           for( j=0; j<*m; j++ )
36              A[i*LDA+j]=a[j*(*lda)+i];
37     }
38     if (transb == CblasNoTrans) {
39        LDB = *n+1;
40        B   = ( float* )malloc( (*k)*LDB*sizeof( float ) );
41        for( i=0; i<*k; i++ )
42           for( j=0; j<*n; j++ )
43              B[i*LDB+j]=b[j*(*ldb)+i];
44     }
45     else {
46        LDB = *k+1;
47        B   = ( float* )malloc( LDB*(*n)*sizeof( float ) );
48        for( i=0; i<*n; i++ )
49           for( j=0; j<*k; j++ )
50              B[i*LDB+j]=b[j*(*ldb)+i];
51     }
52     LDC = *n+1;
53     C   = ( float* )malloc( (*m)*LDC*sizeof( float ) );
54     for( j=0; j<*n; j++ )
55        for( i=0; i<*m; i++ )
56           C[i*LDC+j]=c[j*(*ldc)+i];
57     cblas_sgemm( CblasRowMajor, transa, transb, *m, *n, *k, *alpha, A, LDA,
58                  B, LDB, *beta, C, LDC );
59     for( j=0; j<*n; j++ )
60        for( i=0; i<*m; i++ )
61           c[j*(*ldc)+i]=C[i*LDC+j];
62     free(A);
63     free(B);
64     free(C);
65  }
66  else if (*order == TEST_COL_MJR)
67     cblas_sgemm( CblasColMajor, transa, transb, *m, *n, *k, *alpha, a, *lda,
68                  b, *ldb, *beta, c, *ldc );
69  else
70     cblas_sgemm( UNDEFINED, transa, transb, *m, *n, *k, *alpha, a, *lda,
71                  b, *ldb, *beta, c, *ldc );
72}
73void F77_ssymm(int *order, char *rtlf, char *uplow, int *m, int *n,
74              float *alpha, float *a, int *lda, float *b, int *ldb,
75              float *beta, float *c, int *ldc ) {
76
77  float *A, *B, *C;
78  int i,j,LDA, LDB, LDC;
79  enum CBLAS_UPLO uplo;
80  enum CBLAS_SIDE side;
81
82  get_uplo_type(uplow,&uplo);
83  get_side_type(rtlf,&side);
84
85  if (*order == TEST_ROW_MJR) {
86     if (side == CblasLeft) {
87        LDA = *m+1;
88        A   = ( float* )malloc( (*m)*LDA*sizeof( float ) );
89        for( i=0; i<*m; i++ )
90           for( j=0; j<*m; j++ )
91              A[i*LDA+j]=a[j*(*lda)+i];
92     }
93     else{
94        LDA = *n+1;
95        A   = ( float* )malloc( (*n)*LDA*sizeof( float ) );
96        for( i=0; i<*n; i++ )
97           for( j=0; j<*n; j++ )
98              A[i*LDA+j]=a[j*(*lda)+i];
99     }
100     LDB = *n+1;
101     B   = ( float* )malloc( (*m)*LDB*sizeof( float ) );
102     for( i=0; i<*m; i++ )
103        for( j=0; j<*n; j++ )
104           B[i*LDB+j]=b[j*(*ldb)+i];
105     LDC = *n+1;
106     C   = ( float* )malloc( (*m)*LDC*sizeof( float ) );
107     for( j=0; j<*n; j++ )
108        for( i=0; i<*m; i++ )
109           C[i*LDC+j]=c[j*(*ldc)+i];
110     cblas_ssymm( CblasRowMajor, side, uplo, *m, *n, *alpha, A, LDA, B, LDB,
111                  *beta, C, LDC );
112     for( j=0; j<*n; j++ )
113        for( i=0; i<*m; i++ )
114           c[j*(*ldc)+i]=C[i*LDC+j];
115     free(A);
116     free(B);
117     free(C);
118  }
119  else if (*order == TEST_COL_MJR)
120     cblas_ssymm( CblasColMajor, side, uplo, *m, *n, *alpha, a, *lda, b, *ldb,
121                  *beta, c, *ldc );
122  else
123     cblas_ssymm( UNDEFINED, side, uplo, *m, *n, *alpha, a, *lda, b, *ldb,
124                  *beta, c, *ldc );
125}
126
127void F77_ssyrk(int *order, char *uplow, char *transp, int *n, int *k,
128              float *alpha, float *a, int *lda,
129              float *beta, float *c, int *ldc ) {
130
131  int i,j,LDA,LDC;
132  float *A, *C;
133  enum CBLAS_UPLO uplo;
134  enum CBLAS_TRANSPOSE trans;
135
136  get_uplo_type(uplow,&uplo);
137  get_transpose_type(transp,&trans);
138
139  if (*order == TEST_ROW_MJR) {
140     if (trans == CblasNoTrans) {
141        LDA = *k+1;
142        A   = ( float* )malloc( (*n)*LDA*sizeof( float ) );
143        for( i=0; i<*n; i++ )
144           for( j=0; j<*k; j++ )
145              A[i*LDA+j]=a[j*(*lda)+i];
146     }
147     else{
148        LDA = *n+1;
149        A   = ( float* )malloc( (*k)*LDA*sizeof( float ) );
150        for( i=0; i<*k; i++ )
151           for( j=0; j<*n; j++ )
152              A[i*LDA+j]=a[j*(*lda)+i];
153     }
154     LDC = *n+1;
155     C   = ( float* )malloc( (*n)*LDC*sizeof( float ) );
156     for( i=0; i<*n; i++ )
157        for( j=0; j<*n; j++ )
158           C[i*LDC+j]=c[j*(*ldc)+i];
159     cblas_ssyrk(CblasRowMajor, uplo, trans, *n, *k, *alpha, A, LDA, *beta,
160	         C, LDC );
161     for( j=0; j<*n; j++ )
162        for( i=0; i<*n; i++ )
163           c[j*(*ldc)+i]=C[i*LDC+j];
164     free(A);
165     free(C);
166  }
167  else if (*order == TEST_COL_MJR)
168     cblas_ssyrk(CblasColMajor, uplo, trans, *n, *k, *alpha, a, *lda, *beta,
169	         c, *ldc );
170  else
171     cblas_ssyrk(UNDEFINED, uplo, trans, *n, *k, *alpha, a, *lda, *beta,
172	         c, *ldc );
173}
174
175void F77_ssyr2k(int *order, char *uplow, char *transp, int *n, int *k,
176               float *alpha, float *a, int *lda, float *b, int *ldb,
177               float *beta, float *c, int *ldc ) {
178  int i,j,LDA,LDB,LDC;
179  float *A, *B, *C;
180  enum CBLAS_UPLO uplo;
181  enum CBLAS_TRANSPOSE trans;
182
183  get_uplo_type(uplow,&uplo);
184  get_transpose_type(transp,&trans);
185
186  if (*order == TEST_ROW_MJR) {
187     if (trans == CblasNoTrans) {
188        LDA = *k+1;
189        LDB = *k+1;
190        A   = ( float* )malloc( (*n)*LDA*sizeof( float ) );
191        B   = ( float* )malloc( (*n)*LDB*sizeof( float ) );
192        for( i=0; i<*n; i++ )
193           for( j=0; j<*k; j++ ) {
194              A[i*LDA+j]=a[j*(*lda)+i];
195              B[i*LDB+j]=b[j*(*ldb)+i];
196           }
197     }
198     else {
199        LDA = *n+1;
200        LDB = *n+1;
201        A   = ( float* )malloc( LDA*(*k)*sizeof( float ) );
202        B   = ( float* )malloc( LDB*(*k)*sizeof( float ) );
203        for( i=0; i<*k; i++ )
204           for( j=0; j<*n; j++ ){
205              A[i*LDA+j]=a[j*(*lda)+i];
206              B[i*LDB+j]=b[j*(*ldb)+i];
207           }
208     }
209     LDC = *n+1;
210     C   = ( float* )malloc( (*n)*LDC*sizeof( float ) );
211     for( i=0; i<*n; i++ )
212        for( j=0; j<*n; j++ )
213           C[i*LDC+j]=c[j*(*ldc)+i];
214     cblas_ssyr2k(CblasRowMajor, uplo, trans, *n, *k, *alpha, A, LDA,
215		  B, LDB, *beta, C, LDC );
216     for( j=0; j<*n; j++ )
217        for( i=0; i<*n; i++ )
218           c[j*(*ldc)+i]=C[i*LDC+j];
219     free(A);
220     free(B);
221     free(C);
222  }
223  else if (*order == TEST_COL_MJR)
224     cblas_ssyr2k(CblasColMajor, uplo, trans, *n, *k, *alpha, a, *lda,
225		   b, *ldb, *beta, c, *ldc );
226  else
227     cblas_ssyr2k(UNDEFINED, uplo, trans, *n, *k, *alpha, a, *lda,
228		   b, *ldb, *beta, c, *ldc );
229}
230void F77_strmm(int *order, char *rtlf, char *uplow, char *transp, char *diagn,
231              int *m, int *n, float *alpha, float *a, int *lda, float *b,
232              int *ldb) {
233  int i,j,LDA,LDB;
234  float *A, *B;
235  enum CBLAS_SIDE side;
236  enum CBLAS_DIAG diag;
237  enum CBLAS_UPLO uplo;
238  enum CBLAS_TRANSPOSE trans;
239
240  get_uplo_type(uplow,&uplo);
241  get_transpose_type(transp,&trans);
242  get_diag_type(diagn,&diag);
243  get_side_type(rtlf,&side);
244
245  if (*order == TEST_ROW_MJR) {
246     if (side == CblasLeft) {
247        LDA = *m+1;
248        A   = ( float* )malloc( (*m)*LDA*sizeof( float ) );
249        for( i=0; i<*m; i++ )
250           for( j=0; j<*m; j++ )
251              A[i*LDA+j]=a[j*(*lda)+i];
252     }
253     else{
254        LDA = *n+1;
255        A   = ( float* )malloc( (*n)*LDA*sizeof( float ) );
256        for( i=0; i<*n; i++ )
257           for( j=0; j<*n; j++ )
258              A[i*LDA+j]=a[j*(*lda)+i];
259     }
260     LDB = *n+1;
261     B   = ( float* )malloc( (*m)*LDB*sizeof( float ) );
262     for( i=0; i<*m; i++ )
263        for( j=0; j<*n; j++ )
264           B[i*LDB+j]=b[j*(*ldb)+i];
265     cblas_strmm(CblasRowMajor, side, uplo, trans, diag, *m, *n, *alpha,
266		 A, LDA, B, LDB );
267     for( j=0; j<*n; j++ )
268        for( i=0; i<*m; i++ )
269           b[j*(*ldb)+i]=B[i*LDB+j];
270     free(A);
271     free(B);
272  }
273  else if (*order == TEST_COL_MJR)
274     cblas_strmm(CblasColMajor, side, uplo, trans, diag, *m, *n, *alpha,
275		   a, *lda, b, *ldb);
276  else
277     cblas_strmm(UNDEFINED, side, uplo, trans, diag, *m, *n, *alpha,
278		   a, *lda, b, *ldb);
279}
280
281void F77_strsm(int *order, char *rtlf, char *uplow, char *transp, char *diagn,
282              int *m, int *n, float *alpha, float *a, int *lda, float *b,
283              int *ldb) {
284  int i,j,LDA,LDB;
285  float *A, *B;
286  enum CBLAS_SIDE side;
287  enum CBLAS_DIAG diag;
288  enum CBLAS_UPLO uplo;
289  enum CBLAS_TRANSPOSE trans;
290
291  get_uplo_type(uplow,&uplo);
292  get_transpose_type(transp,&trans);
293  get_diag_type(diagn,&diag);
294  get_side_type(rtlf,&side);
295
296  if (*order == TEST_ROW_MJR) {
297     if (side == CblasLeft) {
298        LDA = *m+1;
299        A   = ( float* )malloc( (*m)*LDA*sizeof( float ) );
300        for( i=0; i<*m; i++ )
301           for( j=0; j<*m; j++ )
302              A[i*LDA+j]=a[j*(*lda)+i];
303     }
304     else{
305        LDA = *n+1;
306        A   = ( float* )malloc( (*n)*LDA*sizeof( float ) );
307        for( i=0; i<*n; i++ )
308           for( j=0; j<*n; j++ )
309              A[i*LDA+j]=a[j*(*lda)+i];
310     }
311     LDB = *n+1;
312     B   = ( float* )malloc( (*m)*LDB*sizeof( float ) );
313     for( i=0; i<*m; i++ )
314        for( j=0; j<*n; j++ )
315           B[i*LDB+j]=b[j*(*ldb)+i];
316     cblas_strsm(CblasRowMajor, side, uplo, trans, diag, *m, *n, *alpha,
317		 A, LDA, B, LDB );
318     for( j=0; j<*n; j++ )
319        for( i=0; i<*m; i++ )
320           b[j*(*ldb)+i]=B[i*LDB+j];
321     free(A);
322     free(B);
323  }
324  else if (*order == TEST_COL_MJR)
325     cblas_strsm(CblasColMajor, side, uplo, trans, diag, *m, *n, *alpha,
326		   a, *lda, b, *ldb);
327  else
328     cblas_strsm(UNDEFINED, side, uplo, trans, diag, *m, *n, *alpha,
329		   a, *lda, b, *ldb);
330}
331