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