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_LASQ5_HEADER
00036 #define TEMPLATE_LAPACK_LASQ5_HEADER
00037
00038 template<class Treal>
00039 int template_lapack_lasq5(integer *i0, integer *n0, Treal *z__,
00040 integer *pp, Treal *tau, Treal *dmin__, Treal *dmin1,
00041 Treal *dmin2, Treal *dn, Treal *dnm1, Treal *dnm2,
00042 logical *ieee)
00043 {
00044
00045 integer i__1;
00046 Treal d__1, d__2;
00047
00048
00049 Treal d__;
00050 integer j4, j4p2;
00051 Treal emin, temp;
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 --z__;
00127
00128
00129 if (*n0 - *i0 - 1 <= 0) {
00130 return 0;
00131 }
00132
00133 j4 = (*i0 << 2) + *pp - 3;
00134 emin = z__[j4 + 4];
00135 d__ = z__[j4] - *tau;
00136 *dmin__ = d__;
00137 *dmin1 = -z__[j4];
00138
00139 if (*ieee) {
00140
00141
00142
00143 if (*pp == 0) {
00144 i__1 = ( *n0 - 3 ) << 2;
00145 for (j4 = *i0 << 2; j4 <= i__1; j4 += 4) {
00146 z__[j4 - 2] = d__ + z__[j4 - 1];
00147 temp = z__[j4 + 1] / z__[j4 - 2];
00148 d__ = d__ * temp - *tau;
00149 *dmin__ = minMACRO(*dmin__,d__);
00150 z__[j4] = z__[j4 - 1] * temp;
00151
00152 d__1 = z__[j4];
00153 emin = minMACRO(d__1,emin);
00154
00155 }
00156 } else {
00157 i__1 = ( *n0 - 3 ) << 2;
00158 for (j4 = *i0 << 2; j4 <= i__1; j4 += 4) {
00159 z__[j4 - 3] = d__ + z__[j4];
00160 temp = z__[j4 + 2] / z__[j4 - 3];
00161 d__ = d__ * temp - *tau;
00162 *dmin__ = minMACRO(*dmin__,d__);
00163 z__[j4 - 1] = z__[j4] * temp;
00164
00165 d__1 = z__[j4 - 1];
00166 emin = minMACRO(d__1,emin);
00167
00168 }
00169 }
00170
00171
00172
00173 *dnm2 = d__;
00174 *dmin2 = *dmin__;
00175 j4 = ( ( *n0 - 2 ) << 2) - *pp;
00176 j4p2 = j4 + (*pp << 1) - 1;
00177 z__[j4 - 2] = *dnm2 + z__[j4p2];
00178 z__[j4] = z__[j4p2 + 2] * (z__[j4p2] / z__[j4 - 2]);
00179 *dnm1 = z__[j4p2 + 2] * (*dnm2 / z__[j4 - 2]) - *tau;
00180 *dmin__ = minMACRO(*dmin__,*dnm1);
00181
00182 *dmin1 = *dmin__;
00183 j4 += 4;
00184 j4p2 = j4 + (*pp << 1) - 1;
00185 z__[j4 - 2] = *dnm1 + z__[j4p2];
00186 z__[j4] = z__[j4p2 + 2] * (z__[j4p2] / z__[j4 - 2]);
00187 *dn = z__[j4p2 + 2] * (*dnm1 / z__[j4 - 2]) - *tau;
00188 *dmin__ = minMACRO(*dmin__,*dn);
00189
00190 } else {
00191
00192
00193
00194 if (*pp == 0) {
00195 i__1 = ( *n0 - 3 ) << 2;
00196 for (j4 = *i0 << 2; j4 <= i__1; j4 += 4) {
00197 z__[j4 - 2] = d__ + z__[j4 - 1];
00198 if (d__ < 0.) {
00199 return 0;
00200 } else {
00201 z__[j4] = z__[j4 + 1] * (z__[j4 - 1] / z__[j4 - 2]);
00202 d__ = z__[j4 + 1] * (d__ / z__[j4 - 2]) - *tau;
00203 }
00204 *dmin__ = minMACRO(*dmin__,d__);
00205
00206 d__1 = emin, d__2 = z__[j4];
00207 emin = minMACRO(d__1,d__2);
00208
00209 }
00210 } else {
00211 i__1 = ( *n0 - 3 ) << 2;
00212 for (j4 = *i0 << 2; j4 <= i__1; j4 += 4) {
00213 z__[j4 - 3] = d__ + z__[j4];
00214 if (d__ < 0.) {
00215 return 0;
00216 } else {
00217 z__[j4 - 1] = z__[j4 + 2] * (z__[j4] / z__[j4 - 3]);
00218 d__ = z__[j4 + 2] * (d__ / z__[j4 - 3]) - *tau;
00219 }
00220 *dmin__ = minMACRO(*dmin__,d__);
00221
00222 d__1 = emin, d__2 = z__[j4 - 1];
00223 emin = minMACRO(d__1,d__2);
00224
00225 }
00226 }
00227
00228
00229
00230 *dnm2 = d__;
00231 *dmin2 = *dmin__;
00232 j4 = ( ( *n0 - 2 ) << 2) - *pp;
00233 j4p2 = j4 + (*pp << 1) - 1;
00234 z__[j4 - 2] = *dnm2 + z__[j4p2];
00235 if (*dnm2 < 0.) {
00236 return 0;
00237 } else {
00238 z__[j4] = z__[j4p2 + 2] * (z__[j4p2] / z__[j4 - 2]);
00239 *dnm1 = z__[j4p2 + 2] * (*dnm2 / z__[j4 - 2]) - *tau;
00240 }
00241 *dmin__ = minMACRO(*dmin__,*dnm1);
00242
00243 *dmin1 = *dmin__;
00244 j4 += 4;
00245 j4p2 = j4 + (*pp << 1) - 1;
00246 z__[j4 - 2] = *dnm1 + z__[j4p2];
00247 if (*dnm1 < 0.) {
00248 return 0;
00249 } else {
00250 z__[j4] = z__[j4p2 + 2] * (z__[j4p2] / z__[j4 - 2]);
00251 *dn = z__[j4p2 + 2] * (*dnm1 / z__[j4 - 2]) - *tau;
00252 }
00253 *dmin__ = minMACRO(*dmin__,*dn);
00254
00255 }
00256
00257 z__[j4 + 2] = *dn;
00258 z__[(*n0 << 2) - *pp] = emin;
00259 return 0;
00260
00261
00262
00263 }
00264
00265 #endif