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

PURPOSE

483 LAPACK auxiliary routine (versionNo3v.e2m)ber 2008                       DLAQR2(1)