2 #include "gmx_lapack.h"
5 F77_FUNC(dlasd5,DLASD5)(int *i__,
13 double b, c__, w, del, tau, delsq;
20 del = d__[2] - d__[1];
21 delsq = del * (d__[2] + d__[1]);
23 w = *rho * 4. * (z__[2] * z__[2] / (d__[1] + d__[2] * 3.) - z__[1] *
24 z__[1] / (d__[1] * 3. + d__[2])) / del + 1.;
26 b = delsq + *rho * (z__[1] * z__[1] + z__[2] * z__[2]);
27 c__ = *rho * z__[1] * z__[1] * delsq;
29 tau = c__ * 2. / (b + sqrt(fabs(b * b - c__ * 4.)));
31 tau /= d__[1] + sqrt(d__[1] * d__[1] + tau);
32 *dsigma = d__[1] + tau;
35 work[1] = d__[1] * 2. + tau;
36 work[2] = d__[1] + tau + d__[2];
38 b = -delsq + *rho * (z__[1] * z__[1] + z__[2] * z__[2]);
39 c__ = *rho * z__[2] * z__[2] * delsq;
42 tau = c__ * -2. / (b + sqrt(b * b + c__ * 4.));
44 tau = (b - sqrt(b * b + c__ * 4.)) / 2.;
47 tau /= d__[2] + sqrt(fabs(d__[2] * d__[2] + tau));
48 *dsigma = d__[2] + tau;
49 delta[1] = -(del + tau);
51 work[1] = d__[1] + tau + d__[2];
52 work[2] = d__[2] * 2. + tau;
56 b = -delsq + *rho * (z__[1] * z__[1] + z__[2] * z__[2]);
57 c__ = *rho * z__[2] * z__[2] * delsq;
60 tau = (b + sqrt(b * b + c__ * 4.)) / 2.;
62 tau = c__ * 2. / (-b + sqrt(b * b + c__ * 4.));
64 tau /= d__[2] + sqrt(d__[2] * d__[2] + tau);
65 *dsigma = d__[2] + tau;
66 delta[1] = -(del + tau);
68 work[1] = d__[1] + tau + d__[2];
69 work[2] = d__[2] * 2. + tau;