LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ strexc()

subroutine strexc ( character  COMPQ,
integer  N,
real, dimension( ldt, * )  T,
integer  LDT,
real, dimension( ldq, * )  Q,
integer  LDQ,
integer  IFST,
integer  ILST,
real, dimension( * )  WORK,
integer  INFO 
)

STREXC

Download STREXC + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 STREXC reorders the real Schur factorization of a real matrix
 A = Q*T*Q**T, so that the diagonal block of T with row index IFST is
 moved to row ILST.

 The real Schur form T is reordered by an orthogonal similarity
 transformation Z**T*T*Z, and optionally the matrix Q of Schur vectors
 is updated by postmultiplying it with Z.

 T must be in Schur canonical form (as returned by SHSEQR), that is,
 block upper triangular with 1-by-1 and 2-by-2 diagonal blocks; each
 2-by-2 diagonal block has its diagonal elements equal and its
 off-diagonal elements of opposite sign.
Parameters
[in]COMPQ
          COMPQ is CHARACTER*1
          = 'V':  update the matrix Q of Schur vectors;
          = 'N':  do not update Q.
[in]N
          N is INTEGER
          The order of the matrix T. N >= 0.
          If N == 0 arguments ILST and IFST may be any value.
[in,out]T
          T is REAL array, dimension (LDT,N)
          On entry, the upper quasi-triangular matrix T, in Schur
          Schur canonical form.
          On exit, the reordered upper quasi-triangular matrix, again
          in Schur canonical form.
[in]LDT
          LDT is INTEGER
          The leading dimension of the array T. LDT >= max(1,N).
[in,out]Q
          Q is REAL array, dimension (LDQ,N)
          On entry, if COMPQ = 'V', the matrix Q of Schur vectors.
          On exit, if COMPQ = 'V', Q has been postmultiplied by the
          orthogonal transformation matrix Z which reorders T.
          If COMPQ = 'N', Q is not referenced.
[in]LDQ
          LDQ is INTEGER
          The leading dimension of the array Q.  LDQ >= 1, and if
          COMPQ = 'V', LDQ >= max(1,N).
[in,out]IFST
          IFST is INTEGER
[in,out]ILST
          ILST is INTEGER

          Specify the reordering of the diagonal blocks of T.
          The block with row index IFST is moved to row ILST, by a
          sequence of transpositions between adjacent blocks.
          On exit, if IFST pointed on entry to the second row of a
          2-by-2 block, it is changed to point to the first row; ILST
          always points to the first row of the block in its final
          position (which may differ from its input value by +1 or -1).
          1 <= IFST <= N; 1 <= ILST <= N.
[out]WORK
          WORK is REAL array, dimension (N)
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          = 1:  two adjacent blocks were too close to swap (the problem
                is very ill-conditioned); T may have been partially
                reordered, and ILST points to the first row of the
                current position of the block being moved.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 146 of file strexc.f.

148 *
149 * -- LAPACK computational routine --
150 * -- LAPACK is a software package provided by Univ. of Tennessee, --
151 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
152 *
153 * .. Scalar Arguments ..
154  CHARACTER COMPQ
155  INTEGER IFST, ILST, INFO, LDQ, LDT, N
156 * ..
157 * .. Array Arguments ..
158  REAL Q( LDQ, * ), T( LDT, * ), WORK( * )
159 * ..
160 *
161 * =====================================================================
162 *
163 * .. Parameters ..
164  REAL ZERO
165  parameter( zero = 0.0e+0 )
166 * ..
167 * .. Local Scalars ..
168  LOGICAL WANTQ
169  INTEGER HERE, NBF, NBL, NBNEXT
170 * ..
171 * .. External Functions ..
172  LOGICAL LSAME
173  EXTERNAL lsame
174 * ..
175 * .. External Subroutines ..
176  EXTERNAL slaexc, xerbla
177 * ..
178 * .. Intrinsic Functions ..
179  INTRINSIC max
180 * ..
181 * .. Executable Statements ..
182 *
183 * Decode and test the input arguments.
184 *
185  info = 0
186  wantq = lsame( compq, 'V' )
187  IF( .NOT.wantq .AND. .NOT.lsame( compq, 'N' ) ) THEN
188  info = -1
189  ELSE IF( n.LT.0 ) THEN
190  info = -2
191  ELSE IF( ldt.LT.max( 1, n ) ) THEN
192  info = -4
193  ELSE IF( ldq.LT.1 .OR. ( wantq .AND. ldq.LT.max( 1, n ) ) ) THEN
194  info = -6
195  ELSE IF(( ifst.LT.1 .OR. ifst.GT.n ).AND.( n.GT.0 )) THEN
196  info = -7
197  ELSE IF(( ilst.LT.1 .OR. ilst.GT.n ).AND.( n.GT.0 )) THEN
198  info = -8
199  END IF
200  IF( info.NE.0 ) THEN
201  CALL xerbla( 'STREXC', -info )
202  RETURN
203  END IF
204 *
205 * Quick return if possible
206 *
207  IF( n.LE.1 )
208  $ RETURN
209 *
210 * Determine the first row of specified block
211 * and find out it is 1 by 1 or 2 by 2.
212 *
213  IF( ifst.GT.1 ) THEN
214  IF( t( ifst, ifst-1 ).NE.zero )
215  $ ifst = ifst - 1
216  END IF
217  nbf = 1
218  IF( ifst.LT.n ) THEN
219  IF( t( ifst+1, ifst ).NE.zero )
220  $ nbf = 2
221  END IF
222 *
223 * Determine the first row of the final block
224 * and find out it is 1 by 1 or 2 by 2.
225 *
226  IF( ilst.GT.1 ) THEN
227  IF( t( ilst, ilst-1 ).NE.zero )
228  $ ilst = ilst - 1
229  END IF
230  nbl = 1
231  IF( ilst.LT.n ) THEN
232  IF( t( ilst+1, ilst ).NE.zero )
233  $ nbl = 2
234  END IF
235 *
236  IF( ifst.EQ.ilst )
237  $ RETURN
238 *
239  IF( ifst.LT.ilst ) THEN
240 *
241 * Update ILST
242 *
243  IF( nbf.EQ.2 .AND. nbl.EQ.1 )
244  $ ilst = ilst - 1
245  IF( nbf.EQ.1 .AND. nbl.EQ.2 )
246  $ ilst = ilst + 1
247 *
248  here = ifst
249 *
250  10 CONTINUE
251 *
252 * Swap block with next one below
253 *
254  IF( nbf.EQ.1 .OR. nbf.EQ.2 ) THEN
255 *
256 * Current block either 1 by 1 or 2 by 2
257 *
258  nbnext = 1
259  IF( here+nbf+1.LE.n ) THEN
260  IF( t( here+nbf+1, here+nbf ).NE.zero )
261  $ nbnext = 2
262  END IF
263  CALL slaexc( wantq, n, t, ldt, q, ldq, here, nbf, nbnext,
264  $ work, info )
265  IF( info.NE.0 ) THEN
266  ilst = here
267  RETURN
268  END IF
269  here = here + nbnext
270 *
271 * Test if 2 by 2 block breaks into two 1 by 1 blocks
272 *
273  IF( nbf.EQ.2 ) THEN
274  IF( t( here+1, here ).EQ.zero )
275  $ nbf = 3
276  END IF
277 *
278  ELSE
279 *
280 * Current block consists of two 1 by 1 blocks each of which
281 * must be swapped individually
282 *
283  nbnext = 1
284  IF( here+3.LE.n ) THEN
285  IF( t( here+3, here+2 ).NE.zero )
286  $ nbnext = 2
287  END IF
288  CALL slaexc( wantq, n, t, ldt, q, ldq, here+1, 1, nbnext,
289  $ work, info )
290  IF( info.NE.0 ) THEN
291  ilst = here
292  RETURN
293  END IF
294  IF( nbnext.EQ.1 ) THEN
295 *
296 * Swap two 1 by 1 blocks, no problems possible
297 *
298  CALL slaexc( wantq, n, t, ldt, q, ldq, here, 1, nbnext,
299  $ work, info )
300  here = here + 1
301  ELSE
302 *
303 * Recompute NBNEXT in case 2 by 2 split
304 *
305  IF( t( here+2, here+1 ).EQ.zero )
306  $ nbnext = 1
307  IF( nbnext.EQ.2 ) THEN
308 *
309 * 2 by 2 Block did not split
310 *
311  CALL slaexc( wantq, n, t, ldt, q, ldq, here, 1,
312  $ nbnext, work, info )
313  IF( info.NE.0 ) THEN
314  ilst = here
315  RETURN
316  END IF
317  here = here + 2
318  ELSE
319 *
320 * 2 by 2 Block did split
321 *
322  CALL slaexc( wantq, n, t, ldt, q, ldq, here, 1, 1,
323  $ work, info )
324  CALL slaexc( wantq, n, t, ldt, q, ldq, here+1, 1, 1,
325  $ work, info )
326  here = here + 2
327  END IF
328  END IF
329  END IF
330  IF( here.LT.ilst )
331  $ GO TO 10
332 *
333  ELSE
334 *
335  here = ifst
336  20 CONTINUE
337 *
338 * Swap block with next one above
339 *
340  IF( nbf.EQ.1 .OR. nbf.EQ.2 ) THEN
341 *
342 * Current block either 1 by 1 or 2 by 2
343 *
344  nbnext = 1
345  IF( here.GE.3 ) THEN
346  IF( t( here-1, here-2 ).NE.zero )
347  $ nbnext = 2
348  END IF
349  CALL slaexc( wantq, n, t, ldt, q, ldq, here-nbnext, nbnext,
350  $ nbf, work, info )
351  IF( info.NE.0 ) THEN
352  ilst = here
353  RETURN
354  END IF
355  here = here - nbnext
356 *
357 * Test if 2 by 2 block breaks into two 1 by 1 blocks
358 *
359  IF( nbf.EQ.2 ) THEN
360  IF( t( here+1, here ).EQ.zero )
361  $ nbf = 3
362  END IF
363 *
364  ELSE
365 *
366 * Current block consists of two 1 by 1 blocks each of which
367 * must be swapped individually
368 *
369  nbnext = 1
370  IF( here.GE.3 ) THEN
371  IF( t( here-1, here-2 ).NE.zero )
372  $ nbnext = 2
373  END IF
374  CALL slaexc( wantq, n, t, ldt, q, ldq, here-nbnext, nbnext,
375  $ 1, work, info )
376  IF( info.NE.0 ) THEN
377  ilst = here
378  RETURN
379  END IF
380  IF( nbnext.EQ.1 ) THEN
381 *
382 * Swap two 1 by 1 blocks, no problems possible
383 *
384  CALL slaexc( wantq, n, t, ldt, q, ldq, here, nbnext, 1,
385  $ work, info )
386  here = here - 1
387  ELSE
388 *
389 * Recompute NBNEXT in case 2 by 2 split
390 *
391  IF( t( here, here-1 ).EQ.zero )
392  $ nbnext = 1
393  IF( nbnext.EQ.2 ) THEN
394 *
395 * 2 by 2 Block did not split
396 *
397  CALL slaexc( wantq, n, t, ldt, q, ldq, here-1, 2, 1,
398  $ work, info )
399  IF( info.NE.0 ) THEN
400  ilst = here
401  RETURN
402  END IF
403  here = here - 2
404  ELSE
405 *
406 * 2 by 2 Block did split
407 *
408  CALL slaexc( wantq, n, t, ldt, q, ldq, here, 1, 1,
409  $ work, info )
410  CALL slaexc( wantq, n, t, ldt, q, ldq, here-1, 1, 1,
411  $ work, info )
412  here = here - 2
413  END IF
414  END IF
415  END IF
416  IF( here.GT.ilst )
417  $ GO TO 20
418  END IF
419  ilst = here
420 *
421  RETURN
422 *
423 * End of STREXC
424 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
subroutine slaexc(WANTQ, N, T, LDT, Q, LDQ, J1, N1, N2, WORK, INFO)
SLAEXC swaps adjacent diagonal blocks of a real upper quasi-triangular matrix in Schur canonical form...
Definition: slaexc.f:138
Here is the call graph for this function:
Here is the caller graph for this function: