LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
dlasq6.f
Go to the documentation of this file.
1 *> \brief \b DLASQ6 computes one dqd transform in ping-pong form. Used by sbdsqr and sstegr.
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 *> \htmlonly
9 *> Download DLASQ6 + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlasq6.f">
11 *> [TGZ]</a>
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlasq6.f">
13 *> [ZIP]</a>
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlasq6.f">
15 *> [TXT]</a>
16 *> \endhtmlonly
17 *
18 * Definition:
19 * ===========
20 *
21 * SUBROUTINE DLASQ6( I0, N0, Z, PP, DMIN, DMIN1, DMIN2, DN,
22 * DNM1, DNM2 )
23 *
24 * .. Scalar Arguments ..
25 * INTEGER I0, N0, PP
26 * DOUBLE PRECISION DMIN, DMIN1, DMIN2, DN, DNM1, DNM2
27 * ..
28 * .. Array Arguments ..
29 * DOUBLE PRECISION Z( * )
30 * ..
31 *
32 *
33 *> \par Purpose:
34 * =============
35 *>
36 *> \verbatim
37 *>
38 *> DLASQ6 computes one dqd (shift equal to zero) transform in
39 *> ping-pong form, with protection against underflow and overflow.
40 *> \endverbatim
41 *
42 * Arguments:
43 * ==========
44 *
45 *> \param[in] I0
46 *> \verbatim
47 *> I0 is INTEGER
48 *> First index.
49 *> \endverbatim
50 *>
51 *> \param[in] N0
52 *> \verbatim
53 *> N0 is INTEGER
54 *> Last index.
55 *> \endverbatim
56 *>
57 *> \param[in] Z
58 *> \verbatim
59 *> Z is DOUBLE PRECISION array, dimension ( 4*N )
60 *> Z holds the qd array. EMIN is stored in Z(4*N0) to avoid
61 *> an extra argument.
62 *> \endverbatim
63 *>
64 *> \param[in] PP
65 *> \verbatim
66 *> PP is INTEGER
67 *> PP=0 for ping, PP=1 for pong.
68 *> \endverbatim
69 *>
70 *> \param[out] DMIN
71 *> \verbatim
72 *> DMIN is DOUBLE PRECISION
73 *> Minimum value of d.
74 *> \endverbatim
75 *>
76 *> \param[out] DMIN1
77 *> \verbatim
78 *> DMIN1 is DOUBLE PRECISION
79 *> Minimum value of d, excluding D( N0 ).
80 *> \endverbatim
81 *>
82 *> \param[out] DMIN2
83 *> \verbatim
84 *> DMIN2 is DOUBLE PRECISION
85 *> Minimum value of d, excluding D( N0 ) and D( N0-1 ).
86 *> \endverbatim
87 *>
88 *> \param[out] DN
89 *> \verbatim
90 *> DN is DOUBLE PRECISION
91 *> d(N0), the last value of d.
92 *> \endverbatim
93 *>
94 *> \param[out] DNM1
95 *> \verbatim
96 *> DNM1 is DOUBLE PRECISION
97 *> d(N0-1).
98 *> \endverbatim
99 *>
100 *> \param[out] DNM2
101 *> \verbatim
102 *> DNM2 is DOUBLE PRECISION
103 *> d(N0-2).
104 *> \endverbatim
105 *
106 * Authors:
107 * ========
108 *
109 *> \author Univ. of Tennessee
110 *> \author Univ. of California Berkeley
111 *> \author Univ. of Colorado Denver
112 *> \author NAG Ltd.
113 *
114 *> \date September 2012
115 *
116 *> \ingroup auxOTHERcomputational
117 *
118 * =====================================================================
119  SUBROUTINE dlasq6( I0, N0, Z, PP, DMIN, DMIN1, DMIN2, DN,
120  $ dnm1, dnm2 )
121 *
122 * -- LAPACK computational routine (version 3.4.2) --
123 * -- LAPACK is a software package provided by Univ. of Tennessee, --
124 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
125 * September 2012
126 *
127 * .. Scalar Arguments ..
128  INTEGER I0, N0, PP
129  DOUBLE PRECISION DMIN, DMIN1, DMIN2, DN, DNM1, DNM2
130 * ..
131 * .. Array Arguments ..
132  DOUBLE PRECISION Z( * )
133 * ..
134 *
135 * =====================================================================
136 *
137 * .. Parameter ..
138  DOUBLE PRECISION ZERO
139  parameter ( zero = 0.0d0 )
140 * ..
141 * .. Local Scalars ..
142  INTEGER J4, J4P2
143  DOUBLE PRECISION D, EMIN, SAFMIN, TEMP
144 * ..
145 * .. External Function ..
146  DOUBLE PRECISION DLAMCH
147  EXTERNAL dlamch
148 * ..
149 * .. Intrinsic Functions ..
150  INTRINSIC min
151 * ..
152 * .. Executable Statements ..
153 *
154  IF( ( n0-i0-1 ).LE.0 )
155  $ RETURN
156 *
157  safmin = dlamch( 'Safe minimum' )
158  j4 = 4*i0 + pp - 3
159  emin = z( j4+4 )
160  d = z( j4 )
161  dmin = d
162 *
163  IF( pp.EQ.0 ) THEN
164  DO 10 j4 = 4*i0, 4*( n0-3 ), 4
165  z( j4-2 ) = d + z( j4-1 )
166  IF( z( j4-2 ).EQ.zero ) THEN
167  z( j4 ) = zero
168  d = z( j4+1 )
169  dmin = d
170  emin = zero
171  ELSE IF( safmin*z( j4+1 ).LT.z( j4-2 ) .AND.
172  $ safmin*z( j4-2 ).LT.z( j4+1 ) ) THEN
173  temp = z( j4+1 ) / z( j4-2 )
174  z( j4 ) = z( j4-1 )*temp
175  d = d*temp
176  ELSE
177  z( j4 ) = z( j4+1 )*( z( j4-1 ) / z( j4-2 ) )
178  d = z( j4+1 )*( d / z( j4-2 ) )
179  END IF
180  dmin = min( dmin, d )
181  emin = min( emin, z( j4 ) )
182  10 CONTINUE
183  ELSE
184  DO 20 j4 = 4*i0, 4*( n0-3 ), 4
185  z( j4-3 ) = d + z( j4 )
186  IF( z( j4-3 ).EQ.zero ) THEN
187  z( j4-1 ) = zero
188  d = z( j4+2 )
189  dmin = d
190  emin = zero
191  ELSE IF( safmin*z( j4+2 ).LT.z( j4-3 ) .AND.
192  $ safmin*z( j4-3 ).LT.z( j4+2 ) ) THEN
193  temp = z( j4+2 ) / z( j4-3 )
194  z( j4-1 ) = z( j4 )*temp
195  d = d*temp
196  ELSE
197  z( j4-1 ) = z( j4+2 )*( z( j4 ) / z( j4-3 ) )
198  d = z( j4+2 )*( d / z( j4-3 ) )
199  END IF
200  dmin = min( dmin, d )
201  emin = min( emin, z( j4-1 ) )
202  20 CONTINUE
203  END IF
204 *
205 * Unroll last two steps.
206 *
207  dnm2 = d
208  dmin2 = dmin
209  j4 = 4*( n0-2 ) - pp
210  j4p2 = j4 + 2*pp - 1
211  z( j4-2 ) = dnm2 + z( j4p2 )
212  IF( z( j4-2 ).EQ.zero ) THEN
213  z( j4 ) = zero
214  dnm1 = z( j4p2+2 )
215  dmin = dnm1
216  emin = zero
217  ELSE IF( safmin*z( j4p2+2 ).LT.z( j4-2 ) .AND.
218  $ safmin*z( j4-2 ).LT.z( j4p2+2 ) ) THEN
219  temp = z( j4p2+2 ) / z( j4-2 )
220  z( j4 ) = z( j4p2 )*temp
221  dnm1 = dnm2*temp
222  ELSE
223  z( j4 ) = z( j4p2+2 )*( z( j4p2 ) / z( j4-2 ) )
224  dnm1 = z( j4p2+2 )*( dnm2 / z( j4-2 ) )
225  END IF
226  dmin = min( dmin, dnm1 )
227 *
228  dmin1 = dmin
229  j4 = j4 + 4
230  j4p2 = j4 + 2*pp - 1
231  z( j4-2 ) = dnm1 + z( j4p2 )
232  IF( z( j4-2 ).EQ.zero ) THEN
233  z( j4 ) = zero
234  dn = z( j4p2+2 )
235  dmin = dn
236  emin = zero
237  ELSE IF( safmin*z( j4p2+2 ).LT.z( j4-2 ) .AND.
238  $ safmin*z( j4-2 ).LT.z( j4p2+2 ) ) THEN
239  temp = z( j4p2+2 ) / z( j4-2 )
240  z( j4 ) = z( j4p2 )*temp
241  dn = dnm1*temp
242  ELSE
243  z( j4 ) = z( j4p2+2 )*( z( j4p2 ) / z( j4-2 ) )
244  dn = z( j4p2+2 )*( dnm1 / z( j4-2 ) )
245  END IF
246  dmin = min( dmin, dn )
247 *
248  z( j4+2 ) = dn
249  z( 4*n0-pp ) = emin
250  RETURN
251 *
252 * End of DLASQ6
253 *
254  END
subroutine dlasq6(I0, N0, Z, PP, DMIN, DMIN1, DMIN2, DN, DNM1, DNM2)
DLASQ6 computes one dqd transform in ping-pong form. Used by sbdsqr and sstegr.
Definition: dlasq6.f:121