1      SUBROUTINE DLASQ5( I0, N0, Z, PP, TAU, DMIN, DMIN1, DMIN2, DN,
2     $                   DNM1, DNM2, IEEE )
3*
4*  -- LAPACK auxiliary routine (version 3.0) --
5*     Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd.,
6*     Courant Institute, Argonne National Lab, and Rice University
7*     May 17, 2000
8*
9*     .. Scalar Arguments ..
10      LOGICAL            IEEE
11      INTEGER            I0, N0, PP
12      DOUBLE PRECISION   DMIN, DMIN1, DMIN2, DN, DNM1, DNM2, TAU
13*     ..
14*     .. Array Arguments ..
15      DOUBLE PRECISION   Z( * )
16*     ..
17*
18*  Purpose
19*  =======
20*
21*  DLASQ5 computes one dqds transform in ping-pong form, one
22*  version for IEEE machines another for non IEEE machines.
23*
24*  Arguments
25*  =========
26*
27*  I0    (input) INTEGER
28*        First index.
29*
30*  N0    (input) INTEGER
31*        Last index.
32*
33*  Z     (input) DOUBLE PRECISION array, dimension ( 4*N )
34*        Z holds the qd array. EMIN is stored in Z(4*N0) to avoid
35*        an extra argument.
36*
37*  PP    (input) INTEGER
38*        PP=0 for ping, PP=1 for pong.
39*
40*  TAU   (input) DOUBLE PRECISION
41*        This is the shift.
42*
43*  DMIN  (output) DOUBLE PRECISION
44*        Minimum value of d.
45*
46*  DMIN1 (output) DOUBLE PRECISION
47*        Minimum value of d, excluding D( N0 ).
48*
49*  DMIN2 (output) DOUBLE PRECISION
50*        Minimum value of d, excluding D( N0 ) and D( N0-1 ).
51*
52*  DN    (output) DOUBLE PRECISION
53*        d(N0), the last value of d.
54*
55*  DNM1  (output) DOUBLE PRECISION
56*        d(N0-1).
57*
58*  DNM2  (output) DOUBLE PRECISION
59*        d(N0-2).
60*
61*  IEEE  (input) LOGICAL
62*        Flag for IEEE or non IEEE arithmetic.
63*
64*  =====================================================================
65*
66*     .. Parameter ..
67      DOUBLE PRECISION   ZERO
68      PARAMETER          ( ZERO = 0.0D0 )
69*     ..
70*     .. Local Scalars ..
71      INTEGER            J4, J4P2
72      DOUBLE PRECISION   D, EMIN, TEMP
73*     ..
74*     .. Intrinsic Functions ..
75      INTRINSIC          MIN
76*     ..
77*     .. Executable Statements ..
78*
79      IF( ( N0-I0-1 ).LE.0 )
80     $   RETURN
81*
82      J4 = 4*I0 + PP - 3
83      EMIN = Z( J4+4 )
84      D = Z( J4 ) - TAU
85      DMIN = D
86      DMIN1 = -Z( J4 )
87*
88      IF( IEEE ) THEN
89*
90*        Code for IEEE arithmetic.
91*
92         IF( PP.EQ.0 ) THEN
93            DO 10 J4 = 4*I0, 4*( N0-3 ), 4
94               Z( J4-2 ) = D + Z( J4-1 )
95               TEMP = Z( J4+1 ) / Z( J4-2 )
96               D = D*TEMP - TAU
97               DMIN = MIN( DMIN, D )
98               Z( J4 ) = Z( J4-1 )*TEMP
99               EMIN = MIN( Z( J4 ), EMIN )
100   10       CONTINUE
101         ELSE
102            DO 20 J4 = 4*I0, 4*( N0-3 ), 4
103               Z( J4-3 ) = D + Z( J4 )
104               TEMP = Z( J4+2 ) / Z( J4-3 )
105               D = D*TEMP - TAU
106               DMIN = MIN( DMIN, D )
107               Z( J4-1 ) = Z( J4 )*TEMP
108               EMIN = MIN( Z( J4-1 ), EMIN )
109   20       CONTINUE
110         END IF
111*
112*        Unroll last two steps.
113*
114         DNM2 = D
115         DMIN2 = DMIN
116         J4 = 4*( N0-2 ) - PP
117         J4P2 = J4 + 2*PP - 1
118         Z( J4-2 ) = DNM2 + Z( J4P2 )
119         Z( J4 ) = Z( J4P2+2 )*( Z( J4P2 ) / Z( J4-2 ) )
120         DNM1 = Z( J4P2+2 )*( DNM2 / Z( J4-2 ) ) - TAU
121         DMIN = MIN( DMIN, DNM1 )
122*
123         DMIN1 = DMIN
124         J4 = J4 + 4
125         J4P2 = J4 + 2*PP - 1
126         Z( J4-2 ) = DNM1 + Z( J4P2 )
127         Z( J4 ) = Z( J4P2+2 )*( Z( J4P2 ) / Z( J4-2 ) )
128         DN = Z( J4P2+2 )*( DNM1 / Z( J4-2 ) ) - TAU
129         DMIN = MIN( DMIN, DN )
130*
131      ELSE
132*
133*        Code for non IEEE arithmetic.
134*
135         IF( PP.EQ.0 ) THEN
136            DO 30 J4 = 4*I0, 4*( N0-3 ), 4
137               Z( J4-2 ) = D + Z( J4-1 )
138               IF( D.LT.ZERO ) THEN
139                  RETURN
140               ELSE
141                  Z( J4 ) = Z( J4+1 )*( Z( J4-1 ) / Z( J4-2 ) )
142                  D = Z( J4+1 )*( D / Z( J4-2 ) ) - TAU
143               END IF
144               DMIN = MIN( DMIN, D )
145               EMIN = MIN( EMIN, Z( J4 ) )
146   30       CONTINUE
147         ELSE
148            DO 40 J4 = 4*I0, 4*( N0-3 ), 4
149               Z( J4-3 ) = D + Z( J4 )
150               IF( D.LT.ZERO ) THEN
151                  RETURN
152               ELSE
153                  Z( J4-1 ) = Z( J4+2 )*( Z( J4 ) / Z( J4-3 ) )
154                  D = Z( J4+2 )*( D / Z( J4-3 ) ) - TAU
155               END IF
156               DMIN = MIN( DMIN, D )
157               EMIN = MIN( EMIN, Z( J4-1 ) )
158   40       CONTINUE
159         END IF
160*
161*        Unroll last two steps.
162*
163         DNM2 = D
164         DMIN2 = DMIN
165         J4 = 4*( N0-2 ) - PP
166         J4P2 = J4 + 2*PP - 1
167         Z( J4-2 ) = DNM2 + Z( J4P2 )
168         IF( DNM2.LT.ZERO ) THEN
169            RETURN
170         ELSE
171            Z( J4 ) = Z( J4P2+2 )*( Z( J4P2 ) / Z( J4-2 ) )
172            DNM1 = Z( J4P2+2 )*( DNM2 / Z( J4-2 ) ) - TAU
173         END IF
174         DMIN = MIN( DMIN, DNM1 )
175*
176         DMIN1 = DMIN
177         J4 = J4 + 4
178         J4P2 = J4 + 2*PP - 1
179         Z( J4-2 ) = DNM1 + Z( J4P2 )
180         IF( DNM1.LT.ZERO ) THEN
181            RETURN
182         ELSE
183            Z( J4 ) = Z( J4P2+2 )*( Z( J4P2 ) / Z( J4-2 ) )
184            DN = Z( J4P2+2 )*( DNM1 / Z( J4-2 ) ) - TAU
185         END IF
186         DMIN = MIN( DMIN, DN )
187*
188      END IF
189*
190      Z( J4+2 ) = DN
191      Z( 4*N0-PP ) = EMIN
192      RETURN
193*
194*     End of DLASQ5
195*
196      END
197