複素線形方程式: 混合精度演算を用いた実連立一次方程式の解

ZCGESV Example Program Data

   4    1                                                 :Value of N, R

 (-1.34, 2.55) ( 0.28, 3.17) (-6.39,-2.20) ( 0.72,-0.92)
 (-0.17,-1.41) ( 3.31,-0.15) (-0.15, 1.34) ( 1.29, 1.38)
 (-3.29,-2.39) (-1.91, 4.42) (-0.14,-1.35) ( 1.72, 1.35)
 ( 2.41, 0.39) (-0.56, 1.47) (-0.83,-0.69) (-1.96, 0.67)  :End of matrix A

 (26.26,51.78) ( 6.43,-8.68) (-5.75,25.31) ( 1.16, 2.57)  :End of vector b

 ZCGESV Example Program Results

 Solution
    ( 1.0000, 1.0000) ( 2.0000,-3.0000) (-4.0000,-5.0000) ( 0.0000, 6.0000)

 Pivot indices
           3          2          3          4

    Program zcgesv_example

!     ZCGESV Example Program Text

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

!     .. Use Statements ..
      Use lapack_interfaces, Only: zcgesv
      Use lapack_precision, Only: sp, dp
!     .. Implicit None Statement ..
      Implicit None
!     .. Parameters ..
      Integer, Parameter :: nin = 5, nout = 6
!     .. Local Scalars ..
      Integer :: i, info, iter, j, lda, ldb, ldx, n, r
!     .. Local Arrays ..
      Complex (Kind=dp), Allocatable :: a(:, :), b(:, :), work(:), x(:, :)
      Complex (Kind=sp), Allocatable :: swork(:)
      Real (Kind=dp), Allocatable :: rwork(:)
      Integer, Allocatable :: ipiv(:)
!     .. Executable Statements ..
      Write (nout, *) 'ZCGESV Example Program Results'
      Write (nout, *)
!     Skip heading in data file
      Read (nin, *)
      Read (nin, *) n, r
      lda = n
      ldb = n
      ldx = n
      Allocate (a(lda,n), b(ldb,r), work(n*r), x(ldx,r), swork(n*(n+ &
        r)), ipiv(n), rwork(n))

!     Read A and B from data file

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

!     Solve the equations Ax = b for x

      Call zcgesv(n, r, a, lda, ipiv, b, ldb, x, ldx, work, swork, rwork, &
        iter, info)

      If (info==0) Then

!       Print solution

        Write (nout, *) 'Solution'
        Write (nout, 100)((x(i,r),j=1,r), i=1, n)

!       Print pivot indices

        Write (nout, *)
        Write (nout, *) 'Pivot indices'
        Write (nout, 110) ipiv(1:n)

      Else
        Write (nout, 120) 'The (', info, ',', info, ')', &
          ' element of the factor U is zero'
      End If

100   Format ((3X,4(' (',F7.4,',',F7.4,')',:)))
110   Format (1X, 7I11)
120   Format (1X, A, I3, A, I3, A, A)
    End Program