3 #include "gromacs/utility/real.h"
4 #include "../gmx_blas.h"
7 F77_FUNC(strmm,STRMM)(const char *side,
19 int a_dim1, a_offset, b_dim1, b_offset, i__1, i__2, i__3;
25 float alpha = *alpha__;
35 a_offset = 1 + a_dim1;
38 b_offset = 1 + b_dim1;
42 lside = (*side=='L' || *side=='l');
48 nounit = (*diag=='N' || *diag=='n');
49 upper = (*uplo=='U' || *uplo=='u');
56 if (fabs(alpha)<GMX_FLOAT_MIN) {
58 for (j = 1; j <= i__1; ++j) {
60 for (i__ = 1; i__ <= i__2; ++i__) {
61 b[i__ + j * b_dim1] = 0.;
67 if (*transa=='N' || *transa=='n') {
70 for (j = 1; j <= i__1; ++j) {
72 for (k = 1; k <= i__2; ++k) {
73 if ( fabs(b[k + j * b_dim1])>GMX_FLOAT_MIN) {
74 temp = alpha * b[k + j * b_dim1];
76 for (i__ = 1; i__ <= i__3; ++i__) {
77 b[i__ + j * b_dim1] += temp * a[i__ + k *
81 temp *= a[k + k * a_dim1];
83 b[k + j * b_dim1] = temp;
89 for (j = 1; j <= i__1; ++j) {
90 for (k = m; k >= 1; --k) {
91 if (fabs(b[k + j * b_dim1])>GMX_FLOAT_MIN) {
92 temp = alpha * b[k + j * b_dim1];
93 b[k + j * b_dim1] = temp;
95 b[k + j * b_dim1] *= a[k + k * a_dim1];
98 for (i__ = k + 1; i__ <= i__2; ++i__) {
99 b[i__ + j * b_dim1] += temp * a[i__ + k *
110 for (j = 1; j <= i__1; ++j) {
111 for (i__ = m; i__ >= 1; --i__) {
112 temp = b[i__ + j * b_dim1];
114 temp *= a[i__ + i__ * a_dim1];
117 for (k = 1; k <= i__2; ++k) {
118 temp += a[k + i__ * a_dim1] * b[k + j * b_dim1];
120 b[i__ + j * b_dim1] = alpha * temp;
125 for (j = 1; j <= i__1; ++j) {
127 for (i__ = 1; i__ <= i__2; ++i__) {
128 temp = b[i__ + j * b_dim1];
130 temp *= a[i__ + i__ * a_dim1];
133 for (k = i__ + 1; k <= i__3; ++k) {
134 temp += a[k + i__ * a_dim1] * b[k + j * b_dim1];
136 b[i__ + j * b_dim1] = alpha * temp;
142 if (*transa=='N' || *transa=='n') {
145 for (j = n; j >= 1; --j) {
148 temp *= a[j + j * a_dim1];
151 for (i__ = 1; i__ <= i__1; ++i__) {
152 b[i__ + j * b_dim1] = temp * b[i__ + j * b_dim1];
155 for (k = 1; k <= i__1; ++k) {
156 if ( fabs(a[k + j * a_dim1])>GMX_FLOAT_MIN) {
157 temp = alpha * a[k + j * a_dim1];
159 for (i__ = 1; i__ <= i__2; ++i__) {
160 b[i__ + j * b_dim1] += temp * b[i__ + k *
168 for (j = 1; j <= i__1; ++j) {
171 temp *= a[j + j * a_dim1];
174 for (i__ = 1; i__ <= i__2; ++i__) {
175 b[i__ + j * b_dim1] = temp * b[i__ + j * b_dim1];
178 for (k = j + 1; k <= i__2; ++k) {
179 if ( fabs(a[k + j * a_dim1])>GMX_FLOAT_MIN) {
180 temp = alpha * a[k + j * a_dim1];
182 for (i__ = 1; i__ <= i__3; ++i__) {
183 b[i__ + j * b_dim1] += temp * b[i__ + k *
194 for (k = 1; k <= i__1; ++k) {
196 for (j = 1; j <= i__2; ++j) {
197 if ( fabs(a[j + k * a_dim1])>GMX_FLOAT_MIN ) {
198 temp = alpha * a[j + k * a_dim1];
200 for (i__ = 1; i__ <= i__3; ++i__) {
201 b[i__ + j * b_dim1] += temp * b[i__ + k *
208 temp *= a[k + k * a_dim1];
210 if ( fabs(temp-1.0)>GMX_FLOAT_EPS) {
212 for (i__ = 1; i__ <= i__2; ++i__) {
213 b[i__ + k * b_dim1] = temp * b[i__ + k * b_dim1];
218 for (k = n; k >= 1; --k) {
220 for (j = k + 1; j <= i__1; ++j) {
221 if ( fabs(a[j + k * a_dim1])>GMX_FLOAT_MIN) {
222 temp = alpha * a[j + k * a_dim1];
224 for (i__ = 1; i__ <= i__2; ++i__) {
225 b[i__ + j * b_dim1] += temp * b[i__ + k *
232 temp *= a[k + k * a_dim1];
234 if ( fabs(temp-1.0)>GMX_FLOAT_EPS) {
236 for (i__ = 1; i__ <= i__1; ++i__) {
237 b[i__ + k * b_dim1] = temp * b[i__ + k * b_dim1];