Program f02jcfe

!     F02JCF Example Program Text

!     Mark 28.5 Release. NAG Copyright 2022.

!     .. Use Statements ..
      Use nag_library, Only: f02jcf, nag_wp, x02alf, x04caf
!     .. Implicit None Statement ..
      Implicit None
!     .. Parameters ..
      Real (Kind=nag_wp), Parameter    :: one = 1.0_nag_wp
      Real (Kind=nag_wp), Parameter    :: tol = 0.0E0_nag_wp
      Real (Kind=nag_wp), Parameter    :: zero = 0.0_nag_wp
      Integer, Parameter               :: nin = 5, nout = 6
!     .. Local Scalars ..
      Real (Kind=nag_wp)               :: inf, tmp
      Integer                          :: i, ifail, iwarn, j, lda, ldb, ldc,   &
                                          ldvl, ldvr, n, scal, sense, tdvl,    &
                                          tdvr
!     .. Local Arrays ..
      Real (Kind=nag_wp), Allocatable  :: a(:,:), alphai(:), alphar(:),        &
                                          b(:,:), beta(:), bevl(:), bevr(:),   &
                                          c(:,:), ei(:), er(:), s(:), vl(:,:), &
                                          vr(:,:)
!     .. Intrinsic Procedures ..
      Intrinsic                        :: abs
!     .. Executable Statements ..
      Write (nout,*) 'F02JCF Example Program Results'
!     Skip heading in data file
      Read (nin,*)
      Read (nin,*) n
      lda = n
      ldb = n
      ldc = n
      ldvl = n
      ldvr = n
      tdvl = 2*n
      tdvr = 2*n
      Allocate (a(lda,n),b(ldb,n),c(ldc,n),alphai(2*n),alphar(2*n),beta(2*n),  &
        ei(2*n),er(2*n),vl(ldvl,tdvl),vr(ldvr,tdvr),s(2*n),bevr(2*n),          &
        bevl(2*n))

!     Read in the matrices A, B and C
      Read (nin,*)(a(i,1:n),i=1,n)
      Read (nin,*)(b(i,1:n),i=1,n)
      Read (nin,*)(c(i,1:n),i=1,n)

!     Use default scaling and compute eigenvalue condition numbers and
!     backward errors for both left and right eigenpairs
      scal = 1
      sense = 7

!     Solve the quadratic eigenvalue problem

      ifail = -1
      Call f02jcf(scal,'V','V',sense,tol,n,a,lda,b,ldb,c,ldc,alphar,alphai,    &
        beta,vl,ldvl,vr,ldvr,s,bevl,bevr,iwarn,ifail)

      If (iwarn/=0) Then
        Write (nout,*)
        Write (nout,99999) 'Warning from f02jcf. IWARN =', iwarn
      End If

      Write (nout,*)
      If (ifail/=0) Then
        Write (nout,99999) 'Failure in f02jcf. IFAIL =', ifail
      Else
!       Infinity
        inf = x02alf()
        Do j = 1, 2*n
          If (beta(j)>=one) Then
            er(j) = alphar(j)/beta(j)
            ei(j) = alphai(j)/beta(j)
          Else
            tmp = inf*beta(j)
            If ((abs(alphar(j))<tmp) .And. (abs(alphai(j))<tmp)) Then
              er(j) = alphar(j)/beta(j)
              ei(j) = alphai(j)/beta(j)
            Else
              er(j) = inf
              ei(j) = zero
            End If
          End If
          If (er(j)<inf) Then
            Write (nout,99998) 'Eigenvalue(', j, ') = (', er(j), ', ', ei(j),  &
              ')'
          Else
            Write (nout,99997) 'Eigenvalue(', j, ') is infinite'
          End If
        End Do

        Write (nout,*)
        Flush (nout)
        ifail = 0
        Call x04caf('General',' ',n,2*n,vr,ldvr,                               &
          'Right eigenvectors (matrix VR)',ifail)

        Write (nout,*)
        Flush (nout)
        ifail = 0
        Call x04caf('General',' ',n,2*n,vl,ldvl,                               &
          'Left eigenvectors (matrix VL)',ifail)

        Write (nout,*)
        Write (nout,*) 'Eigenvalue Condition numbers'
        Do j = 1, 2*n
          Write (nout,99996) s(j)
        End Do

        Write (nout,*)
        Write (nout,*)                                                         &
          'Backward errors for eigenvalues and right eigenvectors'
        Do j = 1, 2*n
          Write (nout,99996) bevr(j)
        End Do

        Write (nout,*)
        Write (nout,*) 'Backward errors for eigenvalues and left eigenvectors'
        Do j = 1, 2*n
          Write (nout,99996) bevl(j)
        End Do
      End If

99999 Format (1X,A,I4)
99998 Format (1X,A,I3,A,1P,E11.4,A,1P,E11.4,A)
99997 Format (1X,A,I3,A)
99996 Format (1X,1P,E11.4)
    End Program f02jcfe