計算ルーチン: 実行列の特異値分解 : (オプションで左および／または右特異ベクトル（高速ヤコビ）)

DGESVJ Example Program Data

   6      4                   :Values of M and N

   2.27  -1.54   1.15  -1.94
   0.28  -1.67   0.94  -0.78
  -0.48  -3.09   0.99  -0.21
   1.07   1.22   0.79   0.63
  -2.35   2.93  -1.45   2.30
   0.62  -7.39   1.03  -2.57  :End of matrix A

 DGESVJ Example Program Results

 Singular values
     9.9966  3.6831  1.3569  0.5000
 Left singular vectors
          1       2       3       4
 1  -0.2774  0.6003 -0.1277  0.1323
 2  -0.2020  0.0301  0.2805  0.7034
 3  -0.2918 -0.3348  0.6453  0.1906
 4   0.0938  0.3699  0.6781 -0.5399
 5   0.4213 -0.5266  0.0413 -0.0575
 6  -0.7816 -0.3353 -0.1645 -0.3957

 Right singular vectors
          1       2       3       4
 1  -0.1921  0.8030  0.0041 -0.5642
 2   0.8794  0.3926 -0.0752  0.2587
 3  -0.2140  0.2980  0.7827  0.5027
 4   0.3795 -0.3351  0.6178 -0.6017

 Error estimate for the singular values
        2.2E-15

 Error estimates for left singular vectors
        3.5E-16    9.5E-16    2.6E-15    4.4E-15

 Error estimates for right singular vectors
        3.5E-16    9.5E-16    2.6E-15    2.6E-15

    Program dgesvj_example

!     DGESVJ Example Program Text

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

!     .. Use Statements ..
      Use lapack_example_aux, Only: nagf_file_print_matrix_real_gen
      Use lapack_interfaces, Only: ddisna, dgesvj
      Use lapack_precision, Only: dp
!     .. Implicit None Statement ..
      Implicit None
!     .. Parameters ..
      Integer, Parameter :: nin = 5, nout = 6
!     .. Local Scalars ..
      Real (Kind=dp) :: eps, serrbd
      Integer :: i, ifail, info, j, lda, ldv, lwork, m, n
!     .. Local Arrays ..
      Real (Kind=dp), Allocatable :: a(:, :), rcondu(:), rcondv(:), s(:), &
        v(:, :), work(:)
!     .. Intrinsic Procedures ..
      Intrinsic :: abs, epsilon
!     .. Executable Statements ..
      Write (nout, *) 'DGESVJ Example Program Results'
      Write (nout, *)
      Flush (nout)
!     Skip heading in data file
      Read (nin, *)
      Read (nin, *) m, n
      lda = m
      ldv = n
      lwork = n + m
      Allocate (a(lda,n), rcondu(m), rcondv(m), s(n), v(ldv,n), work(lwork))

!     Read the m by n matrix A from data file

      Read (nin, *)((a(i,j),j=1,n), i=1, m)

!     Compute the singular values and left and right singular vectors
!     of A (A = U*S*V, m.ge.n)

      Call dgesvj('G', 'U', 'V', m, n, a, lda, s, 0, v, ldv, work, lwork, &
        info)

      If (info==0) Then

!       Compute the approximate error bound for the computed singular values
!       using the 2-norm, s(1) = norm(A), and machine precision, eps.
        eps = epsilon(1.0E0_dp)
        serrbd = eps*s(1)

!       Print solution
        Write (nout, *) 'Singular values'
        Write (nout, 100)(s(j), j=1, n)

        If (abs(work(1)-1.0_dp)>eps) Then
          Write (nout, 130) 'Values need scaling by factor = ', work(1)
        End If
        Flush (nout)

!       ifail: behaviour on error exit
!              =0 for hard exit, =1 for quiet-soft, =-1 for noisy-soft
        ifail = 0
        Call nagf_file_print_matrix_real_gen('General', ' ', m, n, a, lda, &
          'Left singular vectors', ifail)

        Write (nout, *)
        Flush (nout)

        ifail = 0
        Call nagf_file_print_matrix_real_gen('General', ' ', n, n, v, ldv, &
          'Right singular vectors', ifail)

!       Call DDISNA to estimate reciprocal condition
!       numbers for the singular vectors
        Call ddisna('Left', m, n, s, rcondu, info)
        Call ddisna('Right', m, n, s, rcondv, info)

!       Print the approximate error bounds for the singular values
!       and vectors
        Write (nout, '(/1X,A)') 'Error estimate for the singular values'
        Write (nout, 110) serrbd
        Write (nout, '(/1X,A)') 'Error estimates for left singular vectors'
        Write (nout, 110)(serrbd/rcondu(i), i=1, n)
        Write (nout, '(/1X,A)') 'Error estimates for right singular vectors'
        Write (nout, 110)(serrbd/rcondv(i), i=1, n)
      Else
        Write (nout, 120) 'Failure in DGESVJ. INFO =', info
      End If

100   Format (3X, (8F8.4))
110   Format (4X, 1P, 6E11.1)
120   Format (1X, A, I4)
130   Format (/, 1X, A, 1P, E13.5)
    End Program