2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2012,2013, by the GROMACS development team, led by
5 * David van der Spoel, Berk Hess, Erik Lindahl, and including many
6 * others, as listed in the AUTHORS file in the top-level source
7 * directory and at http://www.gromacs.org.
9 * GROMACS is free software; you can redistribute it and/or
10 * modify it under the terms of the GNU Lesser General Public License
11 * as published by the Free Software Foundation; either version 2.1
12 * of the License, or (at your option) any later version.
14 * GROMACS is distributed in the hope that it will be useful,
15 * but WITHOUT ANY WARRANTY; without even the implied warranty of
16 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17 * Lesser General Public License for more details.
19 * You should have received a copy of the GNU Lesser General Public
20 * License along with GROMACS; if not, see
21 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
22 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
24 * If you want to redistribute modifications to GROMACS, please
25 * consider that scientific software is very special. Version
26 * control is crucial - bugs must be traceable. We will be happy to
27 * consider code for inclusion in the official distribution, but
28 * derived work must not be called official GROMACS. Details are found
29 * in the README & COPYING files - if they are missing, get the
30 * official version at http://www.gromacs.org.
32 * To help us fund GROMACS development, we humbly ask that you cite
33 * the research papers on the package. Check out http://www.gromacs.org.
36 #include <types/simple.h>
38 #include "gmx_lapack.h"
39 #include "lapack_limits.h"
42 F77_FUNC(dlasq3,DLASQ3)(int *i0,
88 nn = (*n0 << 2) + *pp;
93 if (z__[nn - 5] > tol2 * (*sigma + z__[nn - 3]) && z__[nn - (*pp << 1) -
94 4] > tol2 * z__[nn - 7]) {
100 z__[(*n0 << 2) - 3] = z__[(*n0 << 2) + *pp - 3] + *sigma;
106 if (z__[nn - 9] > tol2 * *sigma && z__[nn - (*pp << 1) - 8] > tol2 * z__[
113 if (z__[nn - 3] > z__[nn - 7]) {
115 z__[nn - 3] = z__[nn - 7];
118 if (z__[nn - 5] > z__[nn - 3] * tol2) {
119 t = (z__[nn - 7] - z__[nn - 3] + z__[nn - 5]) * .5;
120 s = z__[nn - 3] * (z__[nn - 5] / t);
122 s = z__[nn - 3] * (z__[nn - 5] / (t * (sqrt(s / t + 1.) + 1.)));
124 s = z__[nn - 3] * (z__[nn - 5] / (t + sqrt(t) * sqrt(t + s)));
126 t = z__[nn - 7] + (s + z__[nn - 5]);
127 z__[nn - 3] *= z__[nn - 7] / t;
130 z__[(*n0 << 2) - 7] = z__[nn - 7] + *sigma;
131 z__[(*n0 << 2) - 3] = z__[nn - 3] + *sigma;
140 if (*dmin__ <= 0. || *n0 < n0in) {
141 if (z__[(*i0 << 2) + *pp - 3] * 1.5 < z__[(*n0 << 2) + *pp - 3]) {
142 ipn4 = 4*(*i0 + *n0);
143 i__1 = 2*(*i0 + *n0 - 1);
144 for (j4 = *i0 << 2; j4 <= i__1; j4 += 4) {
146 z__[j4 - 3] = z__[ipn4 - j4 - 3];
147 z__[ipn4 - j4 - 3] = temp;
149 z__[j4 - 2] = z__[ipn4 - j4 - 2];
150 z__[ipn4 - j4 - 2] = temp;
152 z__[j4 - 1] = z__[ipn4 - j4 - 5];
153 z__[ipn4 - j4 - 5] = temp;
155 z__[j4] = z__[ipn4 - j4 - 4];
156 z__[ipn4 - j4 - 4] = temp;
158 if (*n0 - *i0 <= 4) {
159 z__[(*n0 << 2) + *pp - 1] = z__[(*i0 << 2) + *pp - 1];
160 z__[(*n0 << 2) - *pp] = z__[(*i0 << 2) - *pp];
162 d__1 = dmin2, d__2 = z__[(*n0 << 2) + *pp - 1];
163 dmin2 = ((d__1<d__2) ? d__1 : d__2);
164 d__1 = z__[(*n0 << 2) + *pp - 1], d__2 = z__[(*i0 << 2) + *pp - 1]
165 , d__1 = ((d__1<d__2) ? d__1 : d__2), d__2 = z__[(*i0 << 2) + *pp + 3];
166 z__[(*n0 << 2) + *pp - 1] = ((d__1<d__2) ? d__1 : d__2);
167 d__1 = z__[(*n0 << 2) - *pp], d__2 = z__[(*i0 << 2) - *pp], d__1 =
168 ((d__1<d__2) ? d__1 : d__2), d__2 = z__[(*i0 << 2) - *pp + 4];
169 z__[(*n0 << 2) - *pp] = ((d__1<d__2) ? d__1 : d__2);
171 d__2 = z__[(*i0 << 2) + *pp - 3];
172 d__1 = (d__1>d__2) ? d__1 : d__2;
173 d__2 = z__[(*i0 << 2) + *pp + 1];
174 *qmax = ((d__1>d__2) ? d__1 : d__2);
180 F77_FUNC(dlasq4,DLASQ4)(i0, n0, &z__[1], pp, &n0in, dmin__, &dmin1, &dmin2, &dn, &dn1, &
185 F77_FUNC(dlasq5,DLASQ5)(i0, n0, &z__[1], pp, &tau, dmin__, &dmin1, &dmin2, &dn, &dn1, &
188 *ndiv += *n0 - *i0 + 2;
191 if (*dmin__ >= 0. && dmin1 > 0.) {
195 } else if (*dmin__ < 0. && dmin1 > 0. && z__[4*(*n0 - 1) - *pp] < tol *
196 (*sigma + dn1) && fabs(dn) < tol * *sigma) {
198 z__[4*(*n0 - 1) - *pp + 2] = 0.;
201 } else if (*dmin__ < 0.) {
207 } else if (dmin1 > 0.) {
209 tau = (tau + *dmin__) * (1. - eps * 2.);
224 F77_FUNC(dlasq6,DLASQ6)(i0, n0, &z__[1], pp, dmin__, &dmin1, &dmin2, &dn, &dn1, &dn2);
225 *ndiv += *n0 - *i0 + 2;
233 *desig -= t - *sigma;
236 *desig = *sigma - (t - tau) + *desig;