4 #include <types/simple.h>
9 F77_FUNC(dgemm,DGEMM)(const char *transa,
23 const char tra=toupper(*transa);
24 const char trb=toupper(*transb);
27 int nrowa,ncola,nrowb;
36 double alpha = *alpha__;
37 double beta = *beta__;
52 if(m==0 || n==0 || (( fabs(alpha)<GMX_DOUBLE_MIN || k==0) && fabs(beta-1.0)<GMX_DOUBLE_EPS))
55 if(fabs(alpha)<GMX_DOUBLE_MIN) {
56 if(fabs(beta)<GMX_DOUBLE_MIN) {
73 if(fabs(beta)<GMX_DOUBLE_MIN) {
76 } else if(fabs(beta-1.0)>GMX_DOUBLE_EPS) {
81 if( fabs(b[ j*(ldb) + l ])>GMX_DOUBLE_MIN) {
82 temp = alpha * b[ j*(ldb) + l ];
84 c[j*(ldc)+i] += temp * a[l*(lda)+i];
89 /* transpose A, but not B */
94 temp += a[i*(lda)+l] * b[j*(ldb)+l];
95 if(fabs(beta)<GMX_DOUBLE_MIN)
96 c[j*(ldc)+i] = alpha * temp;
98 c[j*(ldc)+i] = alpha * temp + beta * c[j*(ldc)+i];
106 /* transpose B, but not A */
109 if(fabs(beta)<GMX_DOUBLE_MIN) {
112 } else if(fabs(beta-1.0)>GMX_DOUBLE_EPS) {
114 c[j*(ldc)+i] *= beta;
117 if( fabs(b[ l*(ldb) + j ])>GMX_DOUBLE_MIN) {
118 temp = alpha * b[ l*(ldb) + j ];
120 c[j*(ldc)+i] += temp * a[l*(lda)+i];
126 /* Transpose both A and B */
131 temp += a[i*(lda)+l] * b[l*(ldb)+j];
132 if(fabs(beta)<GMX_DOUBLE_MIN)
133 c[j*(ldc)+i] = alpha * temp;
135 c[j*(ldc)+i] = alpha * temp + beta * c[j*(ldc)+i];