37 #ifndef TEMPLATE_LAPACK_LAG2_HEADER
38 #define TEMPLATE_LAPACK_LAG2_HEADER
43 const integer *ldb,
const Treal *safmin, Treal *scale1, Treal *
44 scale2, Treal *wr1, Treal *wr2, Treal *wi)
134 integer a_dim1, a_offset, b_dim1, b_offset;
135 Treal d__1, d__2, d__3, d__4, d__5, d__6;
137 Treal diff, bmin, wbig, wabs, wdet, r__, binv11, binv22,
138 discr, anorm, bnorm, bsize, shift, c1, c2, c3, c4, c5, rtmin,
139 rtmax, wsize, s1, s2, a11, a12, a21, a22, b11, b12, b22, ascale,
140 bscale, pp, qq, ss, wscale, safmax, wsmall, as11, as12, as22, sum,
142 #define a_ref(a_1,a_2) a[(a_2)*a_dim1 + a_1]
143 #define b_ref(a_1,a_2) b[(a_2)*b_dim1 + a_1]
146 a_offset = 1 + a_dim1 * 1;
149 b_offset = 1 + b_dim1 * 1;
155 safmax = 1. / *safmin;
165 a11 = ascale *
a_ref(1, 1);
166 a21 = ascale *
a_ref(2, 1);
167 a12 = ascale *
a_ref(1, 2);
168 a22 = ascale *
a_ref(2, 2);
178 bmin = rtmin *
maxMACRO(d__1,rtmin);
208 as12 = a12 - s1 * b12;
209 as22 = a22 - s1 * b22;
210 ss = a21 * (binv11 * binv22);
211 abi22 = as22 * binv22 - ss * b12;
215 as12 = a12 - s2 * b12;
216 as11 = a11 - s2 * b11;
217 ss = a21 * (binv11 * binv22);
219 pp = (as11 * binv11 + abi22) * .5;
223 if ((d__1 = pp * rtmin,
absMACRO(d__1)) >= 1.) {
226 discr = d__1 * d__1 + qq * *safmin;
231 if (d__1 * d__1 +
absMACRO(qq) <= *safmin) {
234 discr = d__1 * d__1 + qq * safmax;
239 discr = d__1 * d__1 + qq;
250 if (discr >= 0. || r__ == 0.) {
257 wsmall = shift + diff;
261 wdet = (a11 * a22 - a12 * a21) * (binv11 * binv22);
262 wsmall = wdet / wbig;
297 c1 = bsize * (*safmin *
maxMACRO(1.,ascale));
299 c3 = bsize * *safmin;
300 if (ascale <= 1. && bsize <= 1.) {
302 d__1 = 1., d__2 = ascale / *safmin * bsize;
307 if (ascale <= 1. || bsize <= 1.) {
309 d__1 = 1., d__2 = ascale * bsize;
320 d__3 = c4, d__4 =
maxMACRO(wabs,c5) * .5;
321 d__1 =
maxMACRO(*safmin,c1), d__2 = (wabs * c2 + c3) * 1.0000100000000001,
338 *scale1 = ascale * bsize;
349 d__3 = c4, d__4 =
maxMACRO(d__5,c5) * .5;
363 *scale2 = ascale * bsize;