LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ slagv2()

subroutine slagv2 ( real, dimension( lda, * )  A,
integer  LDA,
real, dimension( ldb, * )  B,
integer  LDB,
real, dimension( 2 )  ALPHAR,
real, dimension( 2 )  ALPHAI,
real, dimension( 2 )  BETA,
real  CSL,
real  SNL,
real  CSR,
real  SNR 
)

SLAGV2 computes the Generalized Schur factorization of a real 2-by-2 matrix pencil (A,B) where B is upper triangular.

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

Purpose:
 SLAGV2 computes the Generalized Schur factorization of a real 2-by-2
 matrix pencil (A,B) where B is upper triangular. This routine
 computes orthogonal (rotation) matrices given by CSL, SNL and CSR,
 SNR such that

 1) if the pencil (A,B) has two real eigenvalues (include 0/0 or 1/0
    types), then

    [ a11 a12 ] := [  CSL  SNL ] [ a11 a12 ] [  CSR -SNR ]
    [  0  a22 ]    [ -SNL  CSL ] [ a21 a22 ] [  SNR  CSR ]

    [ b11 b12 ] := [  CSL  SNL ] [ b11 b12 ] [  CSR -SNR ]
    [  0  b22 ]    [ -SNL  CSL ] [  0  b22 ] [  SNR  CSR ],

 2) if the pencil (A,B) has a pair of complex conjugate eigenvalues,
    then

    [ a11 a12 ] := [  CSL  SNL ] [ a11 a12 ] [  CSR -SNR ]
    [ a21 a22 ]    [ -SNL  CSL ] [ a21 a22 ] [  SNR  CSR ]

    [ b11  0  ] := [  CSL  SNL ] [ b11 b12 ] [  CSR -SNR ]
    [  0  b22 ]    [ -SNL  CSL ] [  0  b22 ] [  SNR  CSR ]

    where b11 >= b22 > 0.
Parameters
[in,out]A
          A is REAL array, dimension (LDA, 2)
          On entry, the 2 x 2 matrix A.
          On exit, A is overwritten by the ``A-part'' of the
          generalized Schur form.
[in]LDA
          LDA is INTEGER
          THe leading dimension of the array A.  LDA >= 2.
[in,out]B
          B is REAL array, dimension (LDB, 2)
          On entry, the upper triangular 2 x 2 matrix B.
          On exit, B is overwritten by the ``B-part'' of the
          generalized Schur form.
[in]LDB
          LDB is INTEGER
          THe leading dimension of the array B.  LDB >= 2.
[out]ALPHAR
          ALPHAR is REAL array, dimension (2)
[out]ALPHAI
          ALPHAI is REAL array, dimension (2)
[out]BETA
          BETA is REAL array, dimension (2)
          (ALPHAR(k)+i*ALPHAI(k))/BETA(k) are the eigenvalues of the
          pencil (A,B), k=1,2, i = sqrt(-1).  Note that BETA(k) may
          be zero.
[out]CSL
          CSL is REAL
          The cosine of the left rotation matrix.
[out]SNL
          SNL is REAL
          The sine of the left rotation matrix.
[out]CSR
          CSR is REAL
          The cosine of the right rotation matrix.
[out]SNR
          SNR is REAL
          The sine of the right rotation matrix.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Contributors:
Mark Fahey, Department of Mathematics, Univ. of Kentucky, USA

Definition at line 155 of file slagv2.f.

157 *
158 * -- LAPACK auxiliary routine --
159 * -- LAPACK is a software package provided by Univ. of Tennessee, --
160 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
161 *
162 * .. Scalar Arguments ..
163  INTEGER LDA, LDB
164  REAL CSL, CSR, SNL, SNR
165 * ..
166 * .. Array Arguments ..
167  REAL A( LDA, * ), ALPHAI( 2 ), ALPHAR( 2 ),
168  $ B( LDB, * ), BETA( 2 )
169 * ..
170 *
171 * =====================================================================
172 *
173 * .. Parameters ..
174  REAL ZERO, ONE
175  parameter( zero = 0.0e+0, one = 1.0e+0 )
176 * ..
177 * .. Local Scalars ..
178  REAL ANORM, ASCALE, BNORM, BSCALE, H1, H2, H3, QQ,
179  $ R, RR, SAFMIN, SCALE1, SCALE2, T, ULP, WI, WR1,
180  $ WR2
181 * ..
182 * .. External Subroutines ..
183  EXTERNAL slag2, slartg, slasv2, srot
184 * ..
185 * .. External Functions ..
186  REAL SLAMCH, SLAPY2
187  EXTERNAL slamch, slapy2
188 * ..
189 * .. Intrinsic Functions ..
190  INTRINSIC abs, max
191 * ..
192 * .. Executable Statements ..
193 *
194  safmin = slamch( 'S' )
195  ulp = slamch( 'P' )
196 *
197 * Scale A
198 *
199  anorm = max( abs( a( 1, 1 ) )+abs( a( 2, 1 ) ),
200  $ abs( a( 1, 2 ) )+abs( a( 2, 2 ) ), safmin )
201  ascale = one / anorm
202  a( 1, 1 ) = ascale*a( 1, 1 )
203  a( 1, 2 ) = ascale*a( 1, 2 )
204  a( 2, 1 ) = ascale*a( 2, 1 )
205  a( 2, 2 ) = ascale*a( 2, 2 )
206 *
207 * Scale B
208 *
209  bnorm = max( abs( b( 1, 1 ) ), abs( b( 1, 2 ) )+abs( b( 2, 2 ) ),
210  $ safmin )
211  bscale = one / bnorm
212  b( 1, 1 ) = bscale*b( 1, 1 )
213  b( 1, 2 ) = bscale*b( 1, 2 )
214  b( 2, 2 ) = bscale*b( 2, 2 )
215 *
216 * Check if A can be deflated
217 *
218  IF( abs( a( 2, 1 ) ).LE.ulp ) THEN
219  csl = one
220  snl = zero
221  csr = one
222  snr = zero
223  a( 2, 1 ) = zero
224  b( 2, 1 ) = zero
225  wi = zero
226 *
227 * Check if B is singular
228 *
229  ELSE IF( abs( b( 1, 1 ) ).LE.ulp ) THEN
230  CALL slartg( a( 1, 1 ), a( 2, 1 ), csl, snl, r )
231  csr = one
232  snr = zero
233  CALL srot( 2, a( 1, 1 ), lda, a( 2, 1 ), lda, csl, snl )
234  CALL srot( 2, b( 1, 1 ), ldb, b( 2, 1 ), ldb, csl, snl )
235  a( 2, 1 ) = zero
236  b( 1, 1 ) = zero
237  b( 2, 1 ) = zero
238  wi = zero
239 *
240  ELSE IF( abs( b( 2, 2 ) ).LE.ulp ) THEN
241  CALL slartg( a( 2, 2 ), a( 2, 1 ), csr, snr, t )
242  snr = -snr
243  CALL srot( 2, a( 1, 1 ), 1, a( 1, 2 ), 1, csr, snr )
244  CALL srot( 2, b( 1, 1 ), 1, b( 1, 2 ), 1, csr, snr )
245  csl = one
246  snl = zero
247  a( 2, 1 ) = zero
248  b( 2, 1 ) = zero
249  b( 2, 2 ) = zero
250  wi = zero
251 *
252  ELSE
253 *
254 * B is nonsingular, first compute the eigenvalues of (A,B)
255 *
256  CALL slag2( a, lda, b, ldb, safmin, scale1, scale2, wr1, wr2,
257  $ wi )
258 *
259  IF( wi.EQ.zero ) THEN
260 *
261 * two real eigenvalues, compute s*A-w*B
262 *
263  h1 = scale1*a( 1, 1 ) - wr1*b( 1, 1 )
264  h2 = scale1*a( 1, 2 ) - wr1*b( 1, 2 )
265  h3 = scale1*a( 2, 2 ) - wr1*b( 2, 2 )
266 *
267  rr = slapy2( h1, h2 )
268  qq = slapy2( scale1*a( 2, 1 ), h3 )
269 *
270  IF( rr.GT.qq ) THEN
271 *
272 * find right rotation matrix to zero 1,1 element of
273 * (sA - wB)
274 *
275  CALL slartg( h2, h1, csr, snr, t )
276 *
277  ELSE
278 *
279 * find right rotation matrix to zero 2,1 element of
280 * (sA - wB)
281 *
282  CALL slartg( h3, scale1*a( 2, 1 ), csr, snr, t )
283 *
284  END IF
285 *
286  snr = -snr
287  CALL srot( 2, a( 1, 1 ), 1, a( 1, 2 ), 1, csr, snr )
288  CALL srot( 2, b( 1, 1 ), 1, b( 1, 2 ), 1, csr, snr )
289 *
290 * compute inf norms of A and B
291 *
292  h1 = max( abs( a( 1, 1 ) )+abs( a( 1, 2 ) ),
293  $ abs( a( 2, 1 ) )+abs( a( 2, 2 ) ) )
294  h2 = max( abs( b( 1, 1 ) )+abs( b( 1, 2 ) ),
295  $ abs( b( 2, 1 ) )+abs( b( 2, 2 ) ) )
296 *
297  IF( ( scale1*h1 ).GE.abs( wr1 )*h2 ) THEN
298 *
299 * find left rotation matrix Q to zero out B(2,1)
300 *
301  CALL slartg( b( 1, 1 ), b( 2, 1 ), csl, snl, r )
302 *
303  ELSE
304 *
305 * find left rotation matrix Q to zero out A(2,1)
306 *
307  CALL slartg( a( 1, 1 ), a( 2, 1 ), csl, snl, r )
308 *
309  END IF
310 *
311  CALL srot( 2, a( 1, 1 ), lda, a( 2, 1 ), lda, csl, snl )
312  CALL srot( 2, b( 1, 1 ), ldb, b( 2, 1 ), ldb, csl, snl )
313 *
314  a( 2, 1 ) = zero
315  b( 2, 1 ) = zero
316 *
317  ELSE
318 *
319 * a pair of complex conjugate eigenvalues
320 * first compute the SVD of the matrix B
321 *
322  CALL slasv2( b( 1, 1 ), b( 1, 2 ), b( 2, 2 ), r, t, snr,
323  $ csr, snl, csl )
324 *
325 * Form (A,B) := Q(A,B)Z**T where Q is left rotation matrix and
326 * Z is right rotation matrix computed from SLASV2
327 *
328  CALL srot( 2, a( 1, 1 ), lda, a( 2, 1 ), lda, csl, snl )
329  CALL srot( 2, b( 1, 1 ), ldb, b( 2, 1 ), ldb, csl, snl )
330  CALL srot( 2, a( 1, 1 ), 1, a( 1, 2 ), 1, csr, snr )
331  CALL srot( 2, b( 1, 1 ), 1, b( 1, 2 ), 1, csr, snr )
332 *
333  b( 2, 1 ) = zero
334  b( 1, 2 ) = zero
335 *
336  END IF
337 *
338  END IF
339 *
340 * Unscaling
341 *
342  a( 1, 1 ) = anorm*a( 1, 1 )
343  a( 2, 1 ) = anorm*a( 2, 1 )
344  a( 1, 2 ) = anorm*a( 1, 2 )
345  a( 2, 2 ) = anorm*a( 2, 2 )
346  b( 1, 1 ) = bnorm*b( 1, 1 )
347  b( 2, 1 ) = bnorm*b( 2, 1 )
348  b( 1, 2 ) = bnorm*b( 1, 2 )
349  b( 2, 2 ) = bnorm*b( 2, 2 )
350 *
351  IF( wi.EQ.zero ) THEN
352  alphar( 1 ) = a( 1, 1 )
353  alphar( 2 ) = a( 2, 2 )
354  alphai( 1 ) = zero
355  alphai( 2 ) = zero
356  beta( 1 ) = b( 1, 1 )
357  beta( 2 ) = b( 2, 2 )
358  ELSE
359  alphar( 1 ) = anorm*wr1 / scale1 / bnorm
360  alphai( 1 ) = anorm*wi / scale1 / bnorm
361  alphar( 2 ) = alphar( 1 )
362  alphai( 2 ) = -alphai( 1 )
363  beta( 1 ) = one
364  beta( 2 ) = one
365  END IF
366 *
367  RETURN
368 *
369 * End of SLAGV2
370 *
subroutine slasv2(F, G, H, SSMIN, SSMAX, SNR, CSR, SNL, CSL)
SLASV2 computes the singular value decomposition of a 2-by-2 triangular matrix.
Definition: slasv2.f:138
subroutine slartg(f, g, c, s, r)
SLARTG generates a plane rotation with real cosine and real sine.
Definition: slartg.f90:113
real function slapy2(X, Y)
SLAPY2 returns sqrt(x2+y2).
Definition: slapy2.f:63
subroutine slag2(A, LDA, B, LDB, SAFMIN, SCALE1, SCALE2, WR1, WR2, WI)
SLAG2 computes the eigenvalues of a 2-by-2 generalized eigenvalue problem, with scaling as necessary ...
Definition: slag2.f:156
subroutine srot(N, SX, INCX, SY, INCY, C, S)
SROT
Definition: srot.f:92
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:68
Here is the call graph for this function:
Here is the caller graph for this function: