00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035 #ifndef TEMPLATE_BLAS_SYRK_HEADER
00036 #define TEMPLATE_BLAS_SYRK_HEADER
00037
00038
00039 template<class Treal>
00040 int template_blas_syrk(const char *uplo, const char *trans, const integer *n, const integer *k,
00041 const Treal *alpha, const Treal *a, const integer *lda, const Treal *beta,
00042 Treal *c__, const integer *ldc)
00043 {
00044
00045 integer a_dim1, a_offset, c_dim1, c_offset, i__1, i__2, i__3;
00046
00047 integer info;
00048 Treal temp;
00049 integer i__, j, l;
00050 integer nrowa;
00051 logical upper;
00052 #define a_ref(a_1,a_2) a[(a_2)*a_dim1 + a_1]
00053 #define c___ref(a_1,a_2) c__[(a_2)*c_dim1 + a_1]
00054
00055
00056
00057
00058
00059
00060
00061
00062
00063
00064
00065
00066
00067
00068
00069
00070
00071
00072
00073
00074
00075
00076
00077
00078
00079
00080
00081
00082
00083
00084
00085
00086
00087
00088
00089
00090
00091
00092
00093
00094
00095
00096
00097
00098
00099
00100
00101
00102
00103
00104
00105
00106
00107
00108
00109
00110
00111
00112
00113
00114
00115
00116
00117
00118
00119
00120
00121
00122
00123
00124
00125
00126
00127
00128
00129
00130
00131
00132
00133
00134
00135
00136 a_dim1 = *lda;
00137 a_offset = 1 + a_dim1 * 1;
00138 a -= a_offset;
00139 c_dim1 = *ldc;
00140 c_offset = 1 + c_dim1 * 1;
00141 c__ -= c_offset;
00142
00143 if (template_blas_lsame(trans, "N")) {
00144 nrowa = *n;
00145 } else {
00146 nrowa = *k;
00147 }
00148 upper = template_blas_lsame(uplo, "U");
00149 info = 0;
00150 if (! upper && ! template_blas_lsame(uplo, "L")) {
00151 info = 1;
00152 } else if (! template_blas_lsame(trans, "N") && ! template_blas_lsame(trans,
00153 "T") && ! template_blas_lsame(trans, "C")) {
00154 info = 2;
00155 } else if (*n < 0) {
00156 info = 3;
00157 } else if (*k < 0) {
00158 info = 4;
00159 } else if (*lda < maxMACRO(1,nrowa)) {
00160 info = 7;
00161 } else if (*ldc < maxMACRO(1,*n)) {
00162 info = 10;
00163 }
00164 if (info != 0) {
00165 template_blas_erbla("DSYRK ", &info);
00166 return 0;
00167 }
00168
00169 if (*n == 0 || ( (*alpha == 0. || *k == 0) && *beta == 1. ) ) {
00170 return 0;
00171 }
00172
00173 if (*alpha == 0.) {
00174 if (upper) {
00175 if (*beta == 0.) {
00176 i__1 = *n;
00177 for (j = 1; j <= i__1; ++j) {
00178 i__2 = j;
00179 for (i__ = 1; i__ <= i__2; ++i__) {
00180 c___ref(i__, j) = 0.;
00181
00182 }
00183
00184 }
00185 } else {
00186 i__1 = *n;
00187 for (j = 1; j <= i__1; ++j) {
00188 i__2 = j;
00189 for (i__ = 1; i__ <= i__2; ++i__) {
00190 c___ref(i__, j) = *beta * c___ref(i__, j);
00191
00192 }
00193
00194 }
00195 }
00196 } else {
00197 if (*beta == 0.) {
00198 i__1 = *n;
00199 for (j = 1; j <= i__1; ++j) {
00200 i__2 = *n;
00201 for (i__ = j; i__ <= i__2; ++i__) {
00202 c___ref(i__, j) = 0.;
00203
00204 }
00205
00206 }
00207 } else {
00208 i__1 = *n;
00209 for (j = 1; j <= i__1; ++j) {
00210 i__2 = *n;
00211 for (i__ = j; i__ <= i__2; ++i__) {
00212 c___ref(i__, j) = *beta * c___ref(i__, j);
00213
00214 }
00215
00216 }
00217 }
00218 }
00219 return 0;
00220 }
00221
00222 if (template_blas_lsame(trans, "N")) {
00223
00224 if (upper) {
00225 i__1 = *n;
00226 for (j = 1; j <= i__1; ++j) {
00227 if (*beta == 0.) {
00228 i__2 = j;
00229 for (i__ = 1; i__ <= i__2; ++i__) {
00230 c___ref(i__, j) = 0.;
00231
00232 }
00233 } else if (*beta != 1.) {
00234 i__2 = j;
00235 for (i__ = 1; i__ <= i__2; ++i__) {
00236 c___ref(i__, j) = *beta * c___ref(i__, j);
00237
00238 }
00239 }
00240 i__2 = *k;
00241 for (l = 1; l <= i__2; ++l) {
00242 if (a_ref(j, l) != 0.) {
00243 temp = *alpha * a_ref(j, l);
00244 i__3 = j;
00245 for (i__ = 1; i__ <= i__3; ++i__) {
00246 c___ref(i__, j) = c___ref(i__, j) + temp * a_ref(
00247 i__, l);
00248
00249 }
00250 }
00251
00252 }
00253
00254 }
00255 } else {
00256 i__1 = *n;
00257 for (j = 1; j <= i__1; ++j) {
00258 if (*beta == 0.) {
00259 i__2 = *n;
00260 for (i__ = j; i__ <= i__2; ++i__) {
00261 c___ref(i__, j) = 0.;
00262
00263 }
00264 } else if (*beta != 1.) {
00265 i__2 = *n;
00266 for (i__ = j; i__ <= i__2; ++i__) {
00267 c___ref(i__, j) = *beta * c___ref(i__, j);
00268
00269 }
00270 }
00271 i__2 = *k;
00272 for (l = 1; l <= i__2; ++l) {
00273 if (a_ref(j, l) != 0.) {
00274 temp = *alpha * a_ref(j, l);
00275 i__3 = *n;
00276 for (i__ = j; i__ <= i__3; ++i__) {
00277 c___ref(i__, j) = c___ref(i__, j) + temp * a_ref(
00278 i__, l);
00279
00280 }
00281 }
00282
00283 }
00284
00285 }
00286 }
00287 } else {
00288
00289 if (upper) {
00290 i__1 = *n;
00291 for (j = 1; j <= i__1; ++j) {
00292 i__2 = j;
00293 for (i__ = 1; i__ <= i__2; ++i__) {
00294 temp = 0.;
00295 i__3 = *k;
00296 for (l = 1; l <= i__3; ++l) {
00297 temp += a_ref(l, i__) * a_ref(l, j);
00298
00299 }
00300 if (*beta == 0.) {
00301 c___ref(i__, j) = *alpha * temp;
00302 } else {
00303 c___ref(i__, j) = *alpha * temp + *beta * c___ref(i__,
00304 j);
00305 }
00306
00307 }
00308
00309 }
00310 } else {
00311 i__1 = *n;
00312 for (j = 1; j <= i__1; ++j) {
00313 i__2 = *n;
00314 for (i__ = j; i__ <= i__2; ++i__) {
00315 temp = 0.;
00316 i__3 = *k;
00317 for (l = 1; l <= i__3; ++l) {
00318 temp += a_ref(l, i__) * a_ref(l, j);
00319
00320 }
00321 if (*beta == 0.) {
00322 c___ref(i__, j) = *alpha * temp;
00323 } else {
00324 c___ref(i__, j) = *alpha * temp + *beta * c___ref(i__,
00325 j);
00326 }
00327
00328 }
00329
00330 }
00331 }
00332 }
00333 return 0;
00334
00335 }
00336 #undef c___ref
00337 #undef a_ref
00338
00339 #endif