1DLAQR2(1)           LAPACK auxiliary routine (version 3.1)           DLAQR2(1)
2
3
4

NAME

SYNOPSIS

7       SUBROUTINE DLAQR2( WANTT, WANTZ, N, KTOP, KBOT, NW, H, LDH, ILOZ, IHIZ,
8                          Z, LDZ, NS, ND, SR, SI, V, LDV, NH, T, LDT, NV,  WV,
9                          LDWV, WORK, LWORK )
10
11           INTEGER        IHIZ,  ILOZ,  KBOT,  KTOP, LDH, LDT, LDV, LDWV, LDZ,
12                          LWORK, N, ND, NH, NS, NV, NW
13
14           LOGICAL        WANTT, WANTZ
15
16           DOUBLE         PRECISION H( LDH, * ), SI( * ), SR( * ), T(  LDT,  *
17                          ),  V( LDV, * ), WORK( * ), WV( LDWV, * ), Z( LDZ, *
18                          )
19
20           DOUBLE         PRECISION ZERO, ONE
21
22           PARAMETER      ( ZERO = 0.0d0, ONE = 1.0d0 )
23
24           DOUBLE         PRECISION AA, BB, BETA, CC, CS, DD, EVI,  EVK,  FOO,
25                          S, SAFMAX, SAFMIN, SMLNUM, SN, TAU, ULP
26
27           INTEGER        I,  IFST,  ILST,  INFO, INFQR, J, JW, K, KCOL, KEND,
28                          KLN, KROW, KWTOP, LTOP, LWK1, LWK2, LWKOPT
29
30           LOGICAL        BULGE, SORTED
31
32           DOUBLE         PRECISION DLAMCH
33
34           EXTERNAL       DLAMCH
35
36           EXTERNAL       DCOPY,  DGEHRD,  DGEMM,  DLABAD,   DLACPY,   DLAHQR,
37                          DLANV2, DLARF, DLARFG, DLASET, DORGHR, DTREXC
38
39           INTRINSIC      ABS, DBLE, INT, MAX, MIN, SQRT
40
41           JW             = MIN( NW, KBOT-KTOP+1 )
42
43           IF(            JW.LE.2 ) THEN
44
45           LWKOPT         = 1
46
47           ELSE
48
49           CALL           DGEHRD( JW, 1, JW-1, T, LDT, WORK, WORK, -1, INFO )
50
51           LWK1           = INT( WORK( 1 ) )
52
53           CALL           DORGHR( JW, 1, JW-1, T, LDT, WORK, WORK, -1, INFO )
54
55           LWK2           = INT( WORK( 1 ) )
56
57           LWKOPT         = JW + MAX( LWK1, LWK2 )
58
59           END            IF
60
61           IF(            LWORK.EQ.-1 ) THEN
62
63           WORK(          1 ) = DBLE( LWKOPT )
64
65           RETURN
66
67           END            IF
68
69           NS             = 0
70
71           ND             = 0
72
73           IF(            KTOP.GT.KBOT ) RETURN
74
75           IF(            NW.LT.1 ) RETURN
76
77           SAFMIN         = DLAMCH( 'SAFE MINIMUM' )
78
79           SAFMAX         = ONE / SAFMIN
80
81           CALL           DLABAD( SAFMIN, SAFMAX )
82
83           ULP            = DLAMCH( 'PRECISION' )
84
85           SMLNUM         = SAFMIN*( DBLE( N ) / ULP )
86
87           JW             = MIN( NW, KBOT-KTOP+1 )
88
89           KWTOP          = KBOT - JW + 1
90
91           IF(            KWTOP.EQ.KTOP ) THEN
92
93           S              = ZERO
94
95           ELSE
96
97           S              = H( KWTOP, KWTOP-1 )
98
99           END            IF
100
101           IF(            KBOT.EQ.KWTOP ) THEN
102
103           SR(            KWTOP ) = H( KWTOP, KWTOP )
104
105           SI(            KWTOP ) = ZERO
106
107           NS             = 1
108
109           ND             = 0
110
111           IF(            ABS(  S ).LE.MAX( SMLNUM, ULP*ABS( H( KWTOP, KWTOP )
112                          ) ) ) THEN
113
114           NS             = 0
115
116           ND             = 1
117
118           IF(            KWTOP.GT.KTOP ) H( KWTOP, KWTOP-1 ) = ZERO
119
120           END            IF
121
122           RETURN
123
124           END            IF
125
126           CALL           DLACPY( 'U', JW, JW, H( KWTOP, KWTOP ), LDH, T,  LDT
127                          )
128
129           CALL           DCOPY(  JW-1, H( KWTOP+1, KWTOP ), LDH+1, T( 2, 1 ),
130                          LDT+1 )
131
132           CALL           DLASET( 'A', JW, JW, ZERO, ONE, V, LDV )
133
134           CALL           DLAHQR( .true., .true., JW, 1, JW, T, LDT, SR( KWTOP
135                          ), SI( KWTOP ), 1, JW, V, LDV, INFQR )
136
137           DO             10 J = 1, JW - 3
138
139           T(             J+2, J ) = ZERO
140
141           T(             J+3, J ) = ZERO
142
143           10             CONTINUE
144
145           IF(            JW.GT.2 ) T( JW, JW-2 ) = ZERO
146
147           NS             = JW
148
149           ILST           = INFQR + 1
150
151           20             CONTINUE
152
153           IF(            ILST.LE.NS ) THEN
154
155           IF(            NS.EQ.1 ) THEN
156
157           BULGE          = .FALSE.
158
159           ELSE
160
161           BULGE          = T( NS, NS-1 ).NE.ZERO
162
163           END            IF
164
165           IF(
166
167           FOO            = ABS( T( NS, NS ) )
168
169           IF(            FOO.EQ.ZERO ) FOO = ABS( S )
170
171           IF(            ABS( S*V( 1, NS ) ).LE.MAX( SMLNUM, ULP*FOO ) ) THEN
172
173           NS             = NS - 1
174
175           ELSE
176
177           IFST           = NS
178
179           CALL           DTREXC(  'V',  JW, T, LDT, V, LDV, IFST, ILST, WORK,
180                          INFO )
181
182           ILST           = ILST + 1
183
184           END            IF
185
186           ELSE
187
188           FOO            = ABS( T( NS, NS ) ) + SQRT( ABS( T( NS, NS-1 ) ) )*
189                          SQRT( ABS( T( NS-1, NS ) ) )
190
191           IF(            FOO.EQ.ZERO ) FOO = ABS( S )
192
193           IF(            MAX(  ABS(  S*V(  1,  NS  ) ), ABS( S*V( 1, NS-1 ) )
194                          ).LE.  MAX( SMLNUM, ULP*FOO ) ) THEN
195
196           NS             = NS - 2
197
198           ELSE
199
200           IFST           = NS
201
202           CALL           DTREXC( 'V', JW, T, LDT, V, LDV, IFST,  ILST,  WORK,
203                          INFO )
204
205           ILST           = ILST + 2
206
207           END            IF
208
209           END            IF
210
211           GO             TO 20
212
213           END            IF
214
215           IF(            NS.EQ.0 ) S = ZERO
216
217           IF(            NS.LT.JW ) THEN
218
219           SORTED         = .false.
220
221           I              = NS + 1
222
223           30             CONTINUE
224
225           IF(            SORTED ) GO TO 50
226
227           SORTED         = .true.
228
229           KEND           = I - 1
230
231           I              = INFQR + 1
232
233           IF(            I.EQ.NS ) THEN
234
235           K              = I + 1
236
237           ELSE           IF( T( I+1, I ).EQ.ZERO ) THEN
238
239           K              = I + 1
240
241           ELSE
242
243           K              = I + 2
244
245           END            IF
246
247           40             CONTINUE
248
249           IF(            K.LE.KEND ) THEN
250
251           IF(            K.EQ.I+1 ) THEN
252
253           EVI            = ABS( T( I, I ) )
254
255           ELSE
256
257           EVI            =  ABS(  T(  I,  I ) ) + SQRT( ABS( T( I+1, I ) ) )*
258                          SQRT( ABS( T( I, I+1 ) ) )
259
260           END            IF
261
262           IF(            K.EQ.KEND ) THEN
263
264           EVK            = ABS( T( K, K ) )
265
266           ELSE           IF( T( K+1, K ).EQ.ZERO ) THEN
267
268           EVK            = ABS( T( K, K ) )
269
270           ELSE
271
272           EVK            = ABS( T( K, K ) ) + SQRT( ABS( T( K+1,  K  )  )  )*
273                          SQRT( ABS( T( K, K+1 ) ) )
274
275           END            IF
276
277           IF(            EVI.GE.EVK ) THEN
278
279           I              = K
280
281           ELSE
282
283           SORTED         = .false.
284
285           IFST           = I
286
287           ILST           = K
288
289           CALL           DTREXC(  'V',  JW, T, LDT, V, LDV, IFST, ILST, WORK,
290                          INFO )
291
292           IF(            INFO.EQ.0 ) THEN
293
294           I              = ILST
295
296           ELSE
297
298           I              = K
299
300           END            IF
301
302           END            IF
303
304           IF(            I.EQ.KEND ) THEN
305
306           K              = I + 1
307
308           ELSE           IF( T( I+1, I ).EQ.ZERO ) THEN
309
310           K              = I + 1
311
312           ELSE
313
314           K              = I + 2
315
316           END            IF
317
318           GO             TO 40
319
320           END            IF
321
322           GO             TO 30
323
324           50             CONTINUE
325
326           END            IF
327
328           I              = JW
329
330           60             CONTINUE
331
332           IF(            I.GE.INFQR+1 ) THEN
333
334           IF(            I.EQ.INFQR+1 ) THEN
335
336           SR(            KWTOP+I-1 ) = T( I, I )
337
338           SI(            KWTOP+I-1 ) = ZERO
339
340           I              = I - 1
341
342           ELSE           IF( T( I, I-1 ).EQ.ZERO ) THEN
343
344           SR(            KWTOP+I-1 ) = T( I, I )
345
346           SI(            KWTOP+I-1 ) = ZERO
347
348           I              = I - 1
349
350           ELSE
351
352           AA             = T( I-1, I-1 )
353
354           CC             = T( I, I-1 )
355
356           BB             = T( I-1, I )
357
358           DD             = T( I, I )
359
360           CALL           DLANV2(  AA,  BB,  CC,  DD,  SR(  KWTOP+I-2  ),  SI(
361                          KWTOP+I-2  ),  SR( KWTOP+I-1 ), SI( KWTOP+I-1 ), CS,
362                          SN )
363
364           I              = I - 2
365
366           END            IF
367
368           GO             TO 60
369
370           END            IF
371
372           IF(            NS.LT.JW .OR. S.EQ.ZERO ) THEN
373
374           IF(            NS.GT.1 .AND. S.NE.ZERO ) THEN
375
376           CALL           DCOPY( NS, V, LDV, WORK, 1 )
377
378           BETA           = WORK( 1 )
379
380           CALL           DLARFG( NS, BETA, WORK( 2 ), 1, TAU )
381
382           WORK(          1 ) = ONE
383
384           CALL           DLASET( 'L', JW-2, JW-2, ZERO, ZERO, T( 3, 1 ),  LDT
385                          )
386
387           CALL           DLARF( 'L', NS, JW, WORK, 1, TAU, T, LDT, WORK( JW+1
388                          ) )
389
390           CALL           DLARF( 'R', NS, NS, WORK, 1, TAU, T, LDT, WORK( JW+1
391                          ) )
392
393           CALL           DLARF( 'R', JW, NS, WORK, 1, TAU, V, LDV, WORK( JW+1
394                          ) )
395
396           CALL           DGEHRD( JW, 1, NS,  T,  LDT,  WORK,  WORK(  JW+1  ),
397                          LWORK-JW, INFO )
398
399           END            IF
400
401           IF(            KWTOP.GT.1 ) H( KWTOP, KWTOP-1 ) = S*V( 1, 1 )
402
403           CALL           DLACPY(  'U', JW, JW, T, LDT, H( KWTOP, KWTOP ), LDH
404                          )
405
406           CALL           DCOPY( JW-1, T( 2, 1 ), LDT+1, H( KWTOP+1, KWTOP  ),
407                          LDH+1 )
408
409           IF(            NS.GT.1 .AND. S.NE.ZERO ) THEN
410
411           CALL           DORGHR(  JW,  1,  NS,  T,  LDT,  WORK, WORK( JW+1 ),
412                          LWORK-JW, INFO )
413
414           CALL           DGEMM( 'N', 'N', JW, NS, NS, ONE, V,  LDV,  T,  LDT,
415                          ZERO, WV, LDWV )
416
417           CALL           DLACPY( 'A', JW, NS, WV, LDWV, V, LDV )
418
419           END            IF
420
421           IF(            WANTT ) THEN
422
423           LTOP           = 1
424
425           ELSE
426
427           LTOP           = KTOP
428
429           END            IF
430
431           DO             70 KROW = LTOP, KWTOP - 1, NV
432
433           KLN            = MIN( NV, KWTOP-KROW )
434
435           CALL           DGEMM( 'N', 'N', KLN, JW, JW, ONE, H( KROW, KWTOP ),
436                          LDH, V, LDV, ZERO, WV, LDWV )
437
438           CALL           DLACPY( 'A', KLN, JW, WV, LDWV, H(  KROW,  KWTOP  ),
439                          LDH )
440
441           70             CONTINUE
442
443           IF(            WANTT ) THEN
444
445           DO             80 KCOL = KBOT + 1, N, NH
446
447           KLN            = MIN( NH, N-KCOL+1 )
448
449           CALL           DGEMM( 'C', 'N', JW, KLN, JW, ONE, V, LDV, H( KWTOP,
450                          KCOL ), LDH, ZERO, T, LDT )
451
452           CALL           DLACPY( 'A', JW, KLN, T, LDT, H( KWTOP, KCOL ),  LDH
453                          )
454
455           80             CONTINUE
456
457           END            IF
458
459           IF(            WANTZ ) THEN
460
461           DO             90 KROW = ILOZ, IHIZ, NV
462
463           KLN            = MIN( NV, IHIZ-KROW+1 )
464
465           CALL           DGEMM( 'N', 'N', KLN, JW, JW, ONE, Z( KROW, KWTOP ),
466                          LDZ, V, LDV, ZERO, WV, LDWV )
467
468           CALL           DLACPY( 'A', KLN, JW, WV, LDWV, Z(  KROW,  KWTOP  ),
469                          LDZ )
470
471           90             CONTINUE
472
473           END            IF
474
475           END            IF
476
477           ND             = JW - NS
478
479           NS             = NS - INFQR
480
481           WORK(          1 ) = DBLE( LWKOPT )
482
483           END
484

PURPOSE

486 LAPACK auxiliary routine (versionNo3v.e1m)ber 2006                       DLAQR2(1)