Program f03bn_t1w_fe

!     F03BN_T1W_F Example Program Text

!     Mark 29.2 Release. NAG Copyright 2023.

!     .. Use Statements ..
      Use iso_c_binding, Only: c_ptr
      Use nagad_library, Only: f03bn_t1w_f, f07ar_t1w_f,                       &
                               nagad_t1w_set_derivative, nagad_t1w_w_ctype,    &
                               x10aa_t1w_f, x10ab_t1w_f, Assignment (=)
      Use nag_library, Only: nag_wp, x04daf
!     .. Implicit None Statement ..
      Implicit None
!     .. Parameters ..
      Real (Kind=nag_wp), Parameter    :: two = 2.0_nag_wp
      Integer, Parameter               :: nin = 5, nout = 6
!     .. Local Scalars ..
      Type (c_ptr)                     :: ad_handle
      Type (nagad_t1w_w_ctype)         :: d
      Integer                          :: i, ifail, j, k, lda, n
!     .. Local Arrays ..
      Type (nagad_t1w_w_ctype), Allocatable :: a(:,:), a_in(:,:)
      Complex (Kind=nag_wp), Allocatable :: ac(:,:,:)
      Integer                          :: id(2)
      Integer, Allocatable             :: ipiv(:)
!     .. Intrinsic Procedures ..
      Intrinsic                        :: cmplx
!     .. Executable Statements ..
      Write (nout,*) 'F03BN_T1W_F Example Program Results'

!     Skip heading in data file
      Read (nin,*)

      Read (nin,*) n
      lda = n
      Allocate (a(lda,n),ac(lda,n,2),a_in(lda,n),ipiv(n))

      Read (nin,*)(ac(i,1:n,1),i=1,n)

      a_in(1:n,1:n) = ac(1:n,1:n,1)

!     Create AD configuration data object
      ifail = 0
      Call x10aa_t1w_f(ad_handle,ifail)

      Do i = 1, n
        Do j = 1, n
          Do k = 1, 2
            If (k==1) Then
              Call nagad_t1w_set_derivative(a_in(i,j)%vreal,1.0_nag_wp)
            Else
              Call nagad_t1w_set_derivative(a_in(i,j)%vimag,1.0_nag_wp)
            End If

!           LU Factorize A
!           The NAG name equivalent of zgetrf is f07arf
            a = a_in
            ifail = 0
            Call f07ar_t1w_f(ad_handle,n,n,a,lda,ipiv,ifail)
            If (i==1 .And. j==1 .And. k==1) Then
              Write (nout,*)
              Flush (nout)

              ac(1:n,1:n,1) = a(1:n,1:n)
              ifail = 0
              Call x04daf('g','n',n,n,ac,lda,'Array A after factorization',    &
                ifail)

              Write (nout,*)
              Write (nout,*) 'Pivots'
              Write (nout,99999) ipiv(1:n)
              Write (nout,*)
              Flush (nout)
99999         Format (1X,8(1X,I13))
            End If

            ifail = 0
            Call f03bn_t1w_f(ad_handle,n,a,lda,ipiv,d,id,ifail)

            If (i==1 .And. j==1 .And. k==1) Then
              Write (nout,99998) d%vreal%value, d%vimag%value, id
              Write (nout,*)
              Write (nout,99997) two**id(1)*d%vreal%value,                     &
                two**id(2)*d%vimag%value

99998         Format (1X,'D = (',F13.5,',',F13.5,'), ID = (',I0,', ',I0,')')
99997         Format (1X,'Value of determinant = (',E13.5,',',E13.5,')')
            End If
            If (k==1) Then
              ac(i,j,1) = d%vreal%tangent
              ac(i,j,2) = d%vimag%tangent
              a_in(i,j)%vreal%tangent = 0.0_nag_wp
            Else
              ac(i,j,1) = ac(i,j,1) + cmplx(0.0_nag_wp,d%vreal%tangent,kind=   &
                nag_wp)
              ac(i,j,2) = ac(i,j,2) + cmplx(0.0_nag_wp,d%vimag%tangent,kind=   &
                nag_wp)
              a_in(i,j)%vimag%tangent = 0.0_nag_wp
            End If
          End Do
        End Do
      End Do
      Write (nout,*)
      Write (nout,*) ' Derivatives calculated: First order tangents'
      Write (nout,*) ' Computational mode    : algorithmic'

      Write (nout,*)
      Write (nout,*) ' Derivatives of scaled determinant of A w.r.t. A'
      Write (nout,*)

      Call x04daf('General',' ',n,n,ac,n,'   d(dreal)/da',ifail)

      Write (nout,*)
      Call x04daf('General',' ',n,n,ac(1,1,2),n,'   d(dimag)/da',ifail)

!     Remove computational data object
      Call x10ab_t1w_f(ad_handle,ifail)

    End Program f03bn_t1w_fe