計算ルーチン: 一般化実シュール正規の形式の実行列ペアの指定固有値・固有ベクトルの相反条件数の推定

DTGSNA Example Program Data
  4                    :Value of N
  4.0  1.0  1.0  2.0
  0.0  3.0 -1.0  1.0
  0.0  1.0  3.0  1.0
  0.0  0.0  0.0  6.0   :End of matrix A
  2.0  1.0  1.0  3.0
  0.0  1.0  0.0  1.0
  0.0  0.0  1.0  1.0
  0.0  0.0  0.0  2.0   :End of matrix B

Warning: Floating underflow occurred
 DTGSNA Example Program Results

 S
       1.6E+00    1.7E+00    1.7E+00    1.4E+00

 DIF
       5.4E-01    1.5E-01    1.5E-01    1.2E-01

 Approximate error estimates for eigenvalues of (A,B)
       1.7E-15    1.6E-15    1.6E-15    2.0E-15

 Approximate error estimates for right eigenvectors of (A,B)
       5.0E-15    1.8E-14    1.8E-14    2.2E-14

    Program dtgsna_example

!     DTGSNA Example Program Text

!     Copyright 2017, Numerical Algorithms Group Ltd. http://www.nag.com

!     .. Use Statements ..
      Use lapack_example_aux, Only: nagf_blas_dpyth
      Use lapack_interfaces, Only: dlange, dtgevc, dtgsna
      Use lapack_precision, Only: dp
!     .. Implicit None Statement ..
      Implicit None
!     .. Parameters ..
      Integer, Parameter :: nin = 5, nout = 6
!     .. Local Scalars ..
      Real (Kind=dp) :: eps, snorm, stnrm, tnorm
      Integer :: i, info, lda, ldb, ldvl, ldvr, lwork, m, n
!     .. Local Arrays ..
      Real (Kind=dp), Allocatable :: a(:, :), b(:, :), dif(:), s(:), vl(:, :), &
        vr(:, :), work(:)
      Integer, Allocatable :: iwork(:)
      Logical :: select(1)
!     .. Intrinsic Procedures ..
      Intrinsic :: epsilon
!     .. Executable Statements ..
      Write (nout, *) 'DTGSNA Example Program Results'
      Write (nout, *)
!     Skip heading in data file
      Read (nin, *)
      Read (nin, *) n
      lda = n
      ldb = n
      ldvl = n
      ldvr = n
      lwork = 2*n*(n+2) + 16
      Allocate (a(lda,n), b(ldb,n), dif(n), s(n), vl(ldvl,n), vr(ldvr,n), &
        work(lwork), iwork(n+6))

!     Read A and B from data file

      Read (nin, *)(a(i,1:n), i=1, n)
      Read (nin, *)(b(i,1:n), i=1, n)

!     Calculate the left and right generalized eigenvectors of the
!     pair (A,B). Note that DTGEVC requires WORK to be of dimension
!     at least 6*n.

      Call dtgevc('Both', 'All', select, n, a, lda, b, ldb, vl, ldvl, vr, &
        ldvr, n, m, work, info)

      If (info>0) Then
        Write (nout, 100) info, info + 1
      Else

!       Estimate condition numbers for all the generalized eigenvalues
!       and right eigenvectors of the pair (A,B)

        Call dtgsna('Both', 'All', select, n, a, lda, b, ldb, vl, ldvl, vr, &
          ldvr, s, dif, n, m, work, lwork, iwork, info)

!       Print condition numbers of eigenvalues and right eigenvectors

        Write (nout, *) 'S'
        Write (nout, 110) s(1:m)
        Write (nout, *)
        Write (nout, *) 'DIF'
        Write (nout, 110) dif(1:m)

!       Calculate approximate error estimates

!       Compute the 1-norms of A and B and then compute
!       SQRT(snorm**2 + tnorm**2)

        eps = epsilon(1.0E0_dp)
        snorm = dlange('1-norm', n, n, a, lda, work)
        tnorm = dlange('1-norm', n, n, b, ldb, work)
        stnrm = nagf_blas_dpyth(snorm, tnorm)
        Write (nout, *)
        Write (nout, *) 'Approximate error estimates for eigenvalues of (A,B)'
        Write (nout, 110)(eps*stnrm/s(i), i=1, m)
        Write (nout, *)
        Write (nout, *) 'Approximate error estimates for right ', &
          'eigenvectors of (A,B)'
        Write (nout, 110)(eps*stnrm/dif(i), i=1, m)
      End If

100   Format (' The 2-by-2 block (', I5, ':', I5, ') does not have a co', &
        'mplex eigenvalue')
110   Format ((3X,1P,7E11.1))
    End Program