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_LAPACK_LATRD_HEADER
00036 #define TEMPLATE_LAPACK_LATRD_HEADER
00037
00038
00039 template<class Treal>
00040 int template_lapack_latrd(const char *uplo, const integer *n, const integer *nb, Treal *
00041 a, const integer *lda, Treal *e, Treal *tau, Treal *w,
00042 const integer *ldw)
00043 {
00044
00045
00046
00047
00048
00049
00050
00051
00052
00053
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
00137
00138
00139
00140
00141
00142
00143
00144
00145
00146
00147
00148
00149
00150
00151
00152
00153
00154
00155
00156
00157
00158
00159
00160
00161
00162
00163
00164
00165
00166
00167
00168
00169
00170
00171
00172
00173
00174
00175
00176
00177
00178 Treal c_b5 = -1.;
00179 Treal c_b6 = 1.;
00180 integer c__1 = 1;
00181 Treal c_b16 = 0.;
00182
00183
00184 integer a_dim1, a_offset, w_dim1, w_offset, i__1, i__2, i__3;
00185
00186 integer i__;
00187 Treal alpha;
00188 integer iw;
00189 #define a_ref(a_1,a_2) a[(a_2)*a_dim1 + a_1]
00190 #define w_ref(a_1,a_2) w[(a_2)*w_dim1 + a_1]
00191
00192
00193 a_dim1 = *lda;
00194 a_offset = 1 + a_dim1 * 1;
00195 a -= a_offset;
00196 --e;
00197 --tau;
00198 w_dim1 = *ldw;
00199 w_offset = 1 + w_dim1 * 1;
00200 w -= w_offset;
00201
00202
00203 if (*n <= 0) {
00204 return 0;
00205 }
00206
00207 if (template_blas_lsame(uplo, "U")) {
00208
00209
00210
00211 i__1 = *n - *nb + 1;
00212 for (i__ = *n; i__ >= i__1; --i__) {
00213 iw = i__ - *n + *nb;
00214 if (i__ < *n) {
00215
00216
00217
00218 i__2 = *n - i__;
00219 template_blas_gemv("No transpose", &i__, &i__2, &c_b5, &a_ref(1, i__ + 1),
00220 lda, &w_ref(i__, iw + 1), ldw, &c_b6, &a_ref(1, i__),
00221 &c__1);
00222 i__2 = *n - i__;
00223 template_blas_gemv("No transpose", &i__, &i__2, &c_b5, &w_ref(1, iw + 1),
00224 ldw, &a_ref(i__, i__ + 1), lda, &c_b6, &a_ref(1, i__),
00225 &c__1);
00226 }
00227 if (i__ > 1) {
00228
00229
00230
00231
00232 i__2 = i__ - 1;
00233 template_lapack_larfg(&i__2, &a_ref(i__ - 1, i__), &a_ref(1, i__), &c__1, &
00234 tau[i__ - 1]);
00235 e[i__ - 1] = a_ref(i__ - 1, i__);
00236 a_ref(i__ - 1, i__) = 1.;
00237
00238
00239
00240 i__2 = i__ - 1;
00241 template_blas_symv("Upper", &i__2, &c_b6, &a[a_offset], lda, &a_ref(1,
00242 i__), &c__1, &c_b16, &w_ref(1, iw), &c__1);
00243 if (i__ < *n) {
00244 i__2 = i__ - 1;
00245 i__3 = *n - i__;
00246 template_blas_gemv("Transpose", &i__2, &i__3, &c_b6, &w_ref(1, iw + 1)
00247 , ldw, &a_ref(1, i__), &c__1, &c_b16, &w_ref(i__
00248 + 1, iw), &c__1);
00249 i__2 = i__ - 1;
00250 i__3 = *n - i__;
00251 template_blas_gemv("No transpose", &i__2, &i__3, &c_b5, &a_ref(1, i__
00252 + 1), lda, &w_ref(i__ + 1, iw), &c__1, &c_b6, &
00253 w_ref(1, iw), &c__1);
00254 i__2 = i__ - 1;
00255 i__3 = *n - i__;
00256 template_blas_gemv("Transpose", &i__2, &i__3, &c_b6, &a_ref(1, i__ +
00257 1), lda, &a_ref(1, i__), &c__1, &c_b16, &w_ref(
00258 i__ + 1, iw), &c__1);
00259 i__2 = i__ - 1;
00260 i__3 = *n - i__;
00261 template_blas_gemv("No transpose", &i__2, &i__3, &c_b5, &w_ref(1, iw
00262 + 1), ldw, &w_ref(i__ + 1, iw), &c__1, &c_b6, &
00263 w_ref(1, iw), &c__1);
00264 }
00265 i__2 = i__ - 1;
00266 template_blas_scal(&i__2, &tau[i__ - 1], &w_ref(1, iw), &c__1);
00267 i__2 = i__ - 1;
00268 alpha = tau[i__ - 1] * -.5 * template_blas_dot(&i__2, &w_ref(1, iw), &
00269 c__1, &a_ref(1, i__), &c__1);
00270 i__2 = i__ - 1;
00271 template_blas_axpy(&i__2, &alpha, &a_ref(1, i__), &c__1, &w_ref(1, iw), &
00272 c__1);
00273 }
00274
00275
00276 }
00277 } else {
00278
00279
00280
00281 i__1 = *nb;
00282 for (i__ = 1; i__ <= i__1; ++i__) {
00283
00284
00285
00286 i__2 = *n - i__ + 1;
00287 i__3 = i__ - 1;
00288 template_blas_gemv("No transpose", &i__2, &i__3, &c_b5, &a_ref(i__, 1), lda, &
00289 w_ref(i__, 1), ldw, &c_b6, &a_ref(i__, i__), &c__1);
00290 i__2 = *n - i__ + 1;
00291 i__3 = i__ - 1;
00292 template_blas_gemv("No transpose", &i__2, &i__3, &c_b5, &w_ref(i__, 1), ldw, &
00293 a_ref(i__, 1), lda, &c_b6, &a_ref(i__, i__), &c__1);
00294 if (i__ < *n) {
00295
00296
00297
00298
00299
00300 i__2 = i__ + 2;
00301 i__3 = *n - i__;
00302 template_lapack_larfg(&i__3, &a_ref(i__ + 1, i__), &a_ref(minMACRO(i__2,*n), i__)
00303 , &c__1, &tau[i__]);
00304 e[i__] = a_ref(i__ + 1, i__);
00305 a_ref(i__ + 1, i__) = 1.;
00306
00307
00308
00309 i__2 = *n - i__;
00310 template_blas_symv("Lower", &i__2, &c_b6, &a_ref(i__ + 1, i__ + 1), lda, &
00311 a_ref(i__ + 1, i__), &c__1, &c_b16, &w_ref(i__ + 1,
00312 i__), &c__1);
00313 i__2 = *n - i__;
00314 i__3 = i__ - 1;
00315 template_blas_gemv("Transpose", &i__2, &i__3, &c_b6, &w_ref(i__ + 1, 1),
00316 ldw, &a_ref(i__ + 1, i__), &c__1, &c_b16, &w_ref(1,
00317 i__), &c__1);
00318 i__2 = *n - i__;
00319 i__3 = i__ - 1;
00320 template_blas_gemv("No transpose", &i__2, &i__3, &c_b5, &a_ref(i__ + 1, 1)
00321 , lda, &w_ref(1, i__), &c__1, &c_b6, &w_ref(i__ + 1,
00322 i__), &c__1);
00323 i__2 = *n - i__;
00324 i__3 = i__ - 1;
00325 template_blas_gemv("Transpose", &i__2, &i__3, &c_b6, &a_ref(i__ + 1, 1),
00326 lda, &a_ref(i__ + 1, i__), &c__1, &c_b16, &w_ref(1,
00327 i__), &c__1);
00328 i__2 = *n - i__;
00329 i__3 = i__ - 1;
00330 template_blas_gemv("No transpose", &i__2, &i__3, &c_b5, &w_ref(i__ + 1, 1)
00331 , ldw, &w_ref(1, i__), &c__1, &c_b6, &w_ref(i__ + 1,
00332 i__), &c__1);
00333 i__2 = *n - i__;
00334 template_blas_scal(&i__2, &tau[i__], &w_ref(i__ + 1, i__), &c__1);
00335 i__2 = *n - i__;
00336 alpha = tau[i__] * -.5 * template_blas_dot(&i__2, &w_ref(i__ + 1, i__), &
00337 c__1, &a_ref(i__ + 1, i__), &c__1);
00338 i__2 = *n - i__;
00339 template_blas_axpy(&i__2, &alpha, &a_ref(i__ + 1, i__), &c__1, &w_ref(i__
00340 + 1, i__), &c__1);
00341 }
00342
00343
00344 }
00345 }
00346
00347 return 0;
00348
00349
00350
00351 }
00352
00353 #undef w_ref
00354 #undef a_ref
00355
00356
00357 #endif