root / src / lapack / double / dlasq6.f @ 10
Historique | Voir | Annoter | Télécharger (5,09 ko)
1 | 1 | pfleura2 | SUBROUTINE DLASQ6( I0, N0, Z, PP, DMIN, DMIN1, DMIN2, DN, |
---|---|---|---|
2 | 1 | pfleura2 | $ DNM1, DNM2 ) |
3 | 1 | pfleura2 | * |
4 | 1 | pfleura2 | * -- LAPACK routine (version 3.2) -- |
5 | 1 | pfleura2 | * |
6 | 1 | pfleura2 | * -- Contributed by Osni Marques of the Lawrence Berkeley National -- |
7 | 1 | pfleura2 | * -- Laboratory and Beresford Parlett of the Univ. of California at -- |
8 | 1 | pfleura2 | * -- Berkeley -- |
9 | 1 | pfleura2 | * -- November 2008 -- |
10 | 1 | pfleura2 | * |
11 | 1 | pfleura2 | * -- LAPACK is a software package provided by Univ. of Tennessee, -- |
12 | 1 | pfleura2 | * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- |
13 | 1 | pfleura2 | * |
14 | 1 | pfleura2 | * .. Scalar Arguments .. |
15 | 1 | pfleura2 | INTEGER I0, N0, PP |
16 | 1 | pfleura2 | DOUBLE PRECISION DMIN, DMIN1, DMIN2, DN, DNM1, DNM2 |
17 | 1 | pfleura2 | * .. |
18 | 1 | pfleura2 | * .. Array Arguments .. |
19 | 1 | pfleura2 | DOUBLE PRECISION Z( * ) |
20 | 1 | pfleura2 | * .. |
21 | 1 | pfleura2 | * |
22 | 1 | pfleura2 | * Purpose |
23 | 1 | pfleura2 | * ======= |
24 | 1 | pfleura2 | * |
25 | 1 | pfleura2 | * DLASQ6 computes one dqd (shift equal to zero) transform in |
26 | 1 | pfleura2 | * ping-pong form, with protection against underflow and overflow. |
27 | 1 | pfleura2 | * |
28 | 1 | pfleura2 | * Arguments |
29 | 1 | pfleura2 | * ========= |
30 | 1 | pfleura2 | * |
31 | 1 | pfleura2 | * I0 (input) INTEGER |
32 | 1 | pfleura2 | * First index. |
33 | 1 | pfleura2 | * |
34 | 1 | pfleura2 | * N0 (input) INTEGER |
35 | 1 | pfleura2 | * Last index. |
36 | 1 | pfleura2 | * |
37 | 1 | pfleura2 | * Z (input) DOUBLE PRECISION array, dimension ( 4*N ) |
38 | 1 | pfleura2 | * Z holds the qd array. EMIN is stored in Z(4*N0) to avoid |
39 | 1 | pfleura2 | * an extra argument. |
40 | 1 | pfleura2 | * |
41 | 1 | pfleura2 | * PP (input) INTEGER |
42 | 1 | pfleura2 | * PP=0 for ping, PP=1 for pong. |
43 | 1 | pfleura2 | * |
44 | 1 | pfleura2 | * DMIN (output) DOUBLE PRECISION |
45 | 1 | pfleura2 | * Minimum value of d. |
46 | 1 | pfleura2 | * |
47 | 1 | pfleura2 | * DMIN1 (output) DOUBLE PRECISION |
48 | 1 | pfleura2 | * Minimum value of d, excluding D( N0 ). |
49 | 1 | pfleura2 | * |
50 | 1 | pfleura2 | * DMIN2 (output) DOUBLE PRECISION |
51 | 1 | pfleura2 | * Minimum value of d, excluding D( N0 ) and D( N0-1 ). |
52 | 1 | pfleura2 | * |
53 | 1 | pfleura2 | * DN (output) DOUBLE PRECISION |
54 | 1 | pfleura2 | * d(N0), the last value of d. |
55 | 1 | pfleura2 | * |
56 | 1 | pfleura2 | * DNM1 (output) DOUBLE PRECISION |
57 | 1 | pfleura2 | * d(N0-1). |
58 | 1 | pfleura2 | * |
59 | 1 | pfleura2 | * DNM2 (output) DOUBLE PRECISION |
60 | 1 | pfleura2 | * d(N0-2). |
61 | 1 | pfleura2 | * |
62 | 1 | pfleura2 | * ===================================================================== |
63 | 1 | pfleura2 | * |
64 | 1 | pfleura2 | * .. Parameter .. |
65 | 1 | pfleura2 | DOUBLE PRECISION ZERO |
66 | 1 | pfleura2 | PARAMETER ( ZERO = 0.0D0 ) |
67 | 1 | pfleura2 | * .. |
68 | 1 | pfleura2 | * .. Local Scalars .. |
69 | 1 | pfleura2 | INTEGER J4, J4P2 |
70 | 1 | pfleura2 | DOUBLE PRECISION D, EMIN, SAFMIN, TEMP |
71 | 1 | pfleura2 | * .. |
72 | 1 | pfleura2 | * .. External Function .. |
73 | 1 | pfleura2 | DOUBLE PRECISION DLAMCH |
74 | 1 | pfleura2 | EXTERNAL DLAMCH |
75 | 1 | pfleura2 | * .. |
76 | 1 | pfleura2 | * .. Intrinsic Functions .. |
77 | 1 | pfleura2 | INTRINSIC MIN |
78 | 1 | pfleura2 | * .. |
79 | 1 | pfleura2 | * .. Executable Statements .. |
80 | 1 | pfleura2 | * |
81 | 1 | pfleura2 | IF( ( N0-I0-1 ).LE.0 ) |
82 | 1 | pfleura2 | $ RETURN |
83 | 1 | pfleura2 | * |
84 | 1 | pfleura2 | SAFMIN = DLAMCH( 'Safe minimum' ) |
85 | 1 | pfleura2 | J4 = 4*I0 + PP - 3 |
86 | 1 | pfleura2 | EMIN = Z( J4+4 ) |
87 | 1 | pfleura2 | D = Z( J4 ) |
88 | 1 | pfleura2 | DMIN = D |
89 | 1 | pfleura2 | * |
90 | 1 | pfleura2 | IF( PP.EQ.0 ) THEN |
91 | 1 | pfleura2 | DO 10 J4 = 4*I0, 4*( N0-3 ), 4 |
92 | 1 | pfleura2 | Z( J4-2 ) = D + Z( J4-1 ) |
93 | 1 | pfleura2 | IF( Z( J4-2 ).EQ.ZERO ) THEN |
94 | 1 | pfleura2 | Z( J4 ) = ZERO |
95 | 1 | pfleura2 | D = Z( J4+1 ) |
96 | 1 | pfleura2 | DMIN = D |
97 | 1 | pfleura2 | EMIN = ZERO |
98 | 1 | pfleura2 | ELSE IF( SAFMIN*Z( J4+1 ).LT.Z( J4-2 ) .AND. |
99 | 1 | pfleura2 | $ SAFMIN*Z( J4-2 ).LT.Z( J4+1 ) ) THEN |
100 | 1 | pfleura2 | TEMP = Z( J4+1 ) / Z( J4-2 ) |
101 | 1 | pfleura2 | Z( J4 ) = Z( J4-1 )*TEMP |
102 | 1 | pfleura2 | D = D*TEMP |
103 | 1 | pfleura2 | ELSE |
104 | 1 | pfleura2 | Z( J4 ) = Z( J4+1 )*( Z( J4-1 ) / Z( J4-2 ) ) |
105 | 1 | pfleura2 | D = Z( J4+1 )*( D / Z( J4-2 ) ) |
106 | 1 | pfleura2 | END IF |
107 | 1 | pfleura2 | DMIN = MIN( DMIN, D ) |
108 | 1 | pfleura2 | EMIN = MIN( EMIN, Z( J4 ) ) |
109 | 1 | pfleura2 | 10 CONTINUE |
110 | 1 | pfleura2 | ELSE |
111 | 1 | pfleura2 | DO 20 J4 = 4*I0, 4*( N0-3 ), 4 |
112 | 1 | pfleura2 | Z( J4-3 ) = D + Z( J4 ) |
113 | 1 | pfleura2 | IF( Z( J4-3 ).EQ.ZERO ) THEN |
114 | 1 | pfleura2 | Z( J4-1 ) = ZERO |
115 | 1 | pfleura2 | D = Z( J4+2 ) |
116 | 1 | pfleura2 | DMIN = D |
117 | 1 | pfleura2 | EMIN = ZERO |
118 | 1 | pfleura2 | ELSE IF( SAFMIN*Z( J4+2 ).LT.Z( J4-3 ) .AND. |
119 | 1 | pfleura2 | $ SAFMIN*Z( J4-3 ).LT.Z( J4+2 ) ) THEN |
120 | 1 | pfleura2 | TEMP = Z( J4+2 ) / Z( J4-3 ) |
121 | 1 | pfleura2 | Z( J4-1 ) = Z( J4 )*TEMP |
122 | 1 | pfleura2 | D = D*TEMP |
123 | 1 | pfleura2 | ELSE |
124 | 1 | pfleura2 | Z( J4-1 ) = Z( J4+2 )*( Z( J4 ) / Z( J4-3 ) ) |
125 | 1 | pfleura2 | D = Z( J4+2 )*( D / Z( J4-3 ) ) |
126 | 1 | pfleura2 | END IF |
127 | 1 | pfleura2 | DMIN = MIN( DMIN, D ) |
128 | 1 | pfleura2 | EMIN = MIN( EMIN, Z( J4-1 ) ) |
129 | 1 | pfleura2 | 20 CONTINUE |
130 | 1 | pfleura2 | END IF |
131 | 1 | pfleura2 | * |
132 | 1 | pfleura2 | * Unroll last two steps. |
133 | 1 | pfleura2 | * |
134 | 1 | pfleura2 | DNM2 = D |
135 | 1 | pfleura2 | DMIN2 = DMIN |
136 | 1 | pfleura2 | J4 = 4*( N0-2 ) - PP |
137 | 1 | pfleura2 | J4P2 = J4 + 2*PP - 1 |
138 | 1 | pfleura2 | Z( J4-2 ) = DNM2 + Z( J4P2 ) |
139 | 1 | pfleura2 | IF( Z( J4-2 ).EQ.ZERO ) THEN |
140 | 1 | pfleura2 | Z( J4 ) = ZERO |
141 | 1 | pfleura2 | DNM1 = Z( J4P2+2 ) |
142 | 1 | pfleura2 | DMIN = DNM1 |
143 | 1 | pfleura2 | EMIN = ZERO |
144 | 1 | pfleura2 | ELSE IF( SAFMIN*Z( J4P2+2 ).LT.Z( J4-2 ) .AND. |
145 | 1 | pfleura2 | $ SAFMIN*Z( J4-2 ).LT.Z( J4P2+2 ) ) THEN |
146 | 1 | pfleura2 | TEMP = Z( J4P2+2 ) / Z( J4-2 ) |
147 | 1 | pfleura2 | Z( J4 ) = Z( J4P2 )*TEMP |
148 | 1 | pfleura2 | DNM1 = DNM2*TEMP |
149 | 1 | pfleura2 | ELSE |
150 | 1 | pfleura2 | Z( J4 ) = Z( J4P2+2 )*( Z( J4P2 ) / Z( J4-2 ) ) |
151 | 1 | pfleura2 | DNM1 = Z( J4P2+2 )*( DNM2 / Z( J4-2 ) ) |
152 | 1 | pfleura2 | END IF |
153 | 1 | pfleura2 | DMIN = MIN( DMIN, DNM1 ) |
154 | 1 | pfleura2 | * |
155 | 1 | pfleura2 | DMIN1 = DMIN |
156 | 1 | pfleura2 | J4 = J4 + 4 |
157 | 1 | pfleura2 | J4P2 = J4 + 2*PP - 1 |
158 | 1 | pfleura2 | Z( J4-2 ) = DNM1 + Z( J4P2 ) |
159 | 1 | pfleura2 | IF( Z( J4-2 ).EQ.ZERO ) THEN |
160 | 1 | pfleura2 | Z( J4 ) = ZERO |
161 | 1 | pfleura2 | DN = Z( J4P2+2 ) |
162 | 1 | pfleura2 | DMIN = DN |
163 | 1 | pfleura2 | EMIN = ZERO |
164 | 1 | pfleura2 | ELSE IF( SAFMIN*Z( J4P2+2 ).LT.Z( J4-2 ) .AND. |
165 | 1 | pfleura2 | $ SAFMIN*Z( J4-2 ).LT.Z( J4P2+2 ) ) THEN |
166 | 1 | pfleura2 | TEMP = Z( J4P2+2 ) / Z( J4-2 ) |
167 | 1 | pfleura2 | Z( J4 ) = Z( J4P2 )*TEMP |
168 | 1 | pfleura2 | DN = DNM1*TEMP |
169 | 1 | pfleura2 | ELSE |
170 | 1 | pfleura2 | Z( J4 ) = Z( J4P2+2 )*( Z( J4P2 ) / Z( J4-2 ) ) |
171 | 1 | pfleura2 | DN = Z( J4P2+2 )*( DNM1 / Z( J4-2 ) ) |
172 | 1 | pfleura2 | END IF |
173 | 1 | pfleura2 | DMIN = MIN( DMIN, DN ) |
174 | 1 | pfleura2 | * |
175 | 1 | pfleura2 | Z( J4+2 ) = DN |
176 | 1 | pfleura2 | Z( 4*N0-PP ) = EMIN |
177 | 1 | pfleura2 | RETURN |
178 | 1 | pfleura2 | * |
179 | 1 | pfleura2 | * End of DLASQ6 |
180 | 1 | pfleura2 | * |
181 | 1 | pfleura2 | END |