!   E04UCA Example Program Text
!   Mark 26.1 Release. NAG Copyright 2016.

    Module e04ucae_mod

!     E04UCA Example Program Module:
!            Parameters and User-defined Routines

!     .. Use Statements ..
      Use nag_library, Only: nag_wp
!     .. Implicit None Statement ..
      Implicit None
!     .. Accessibility Statements ..
      Private
      Public                           :: confun, objfun
!     .. Parameters ..
      Real (Kind=nag_wp), Parameter, Public :: one = 1.0_nag_wp
      Real (Kind=nag_wp), Parameter, Public :: zero = 0.0_nag_wp
      Integer, Parameter, Public       :: inc1 = 1, lcwsav = 1, liwsav = 610,  &
                                          llwsav = 120, lrwsav = 475, nin = 5, &
                                          nout = 6
    Contains
      Subroutine objfun(mode,n,x,objf,objgrd,nstate,iuser,ruser)
!       Routine to evaluate objective function and its 1st derivatives.

!       .. Scalar Arguments ..
        Real (Kind=nag_wp), Intent (Out) :: objf
        Integer, Intent (Inout)        :: mode
        Integer, Intent (In)           :: n, nstate
!       .. Array Arguments ..
        Real (Kind=nag_wp), Intent (Inout) :: objgrd(n), ruser(*)
        Real (Kind=nag_wp), Intent (In) :: x(n)
        Integer, Intent (Inout)        :: iuser(*)
!       .. Executable Statements ..
        If (mode==0 .Or. mode==2) Then
          objf = x(1)*x(4)*(x(1)+x(2)+x(3)) + x(3)
        End If

        If (mode==1 .Or. mode==2) Then
          objgrd(1) = x(4)*(x(1)+x(1)+x(2)+x(3))
          objgrd(2) = x(1)*x(4)
          objgrd(3) = x(1)*x(4) + one
          objgrd(4) = x(1)*(x(1)+x(2)+x(3))
        End If

        Return

      End Subroutine objfun
      Subroutine confun(mode,ncnln,n,ldcj,needc,x,c,cjac,nstate,iuser,ruser)
!       Routine to evaluate the nonlinear constraints and their 1st
!       derivatives.

!       .. Scalar Arguments ..
        Integer, Intent (In)           :: ldcj, n, ncnln, nstate
        Integer, Intent (Inout)        :: mode
!       .. Array Arguments ..
        Real (Kind=nag_wp), Intent (Out) :: c(ncnln)
        Real (Kind=nag_wp), Intent (Inout) :: cjac(ldcj,n), ruser(*)
        Real (Kind=nag_wp), Intent (In) :: x(n)
        Integer, Intent (Inout)        :: iuser(*)
        Integer, Intent (In)           :: needc(ncnln)
!       .. Executable Statements ..
        If (nstate==1) Then

!         First call to CONFUN.  Set all Jacobian elements to zero.
!         Note that this will only work when 'Derivative Level = 3'
!         (the default; see Section 11.2).

          cjac(1:ncnln,1:n) = zero
        End If

        If (needc(1)>0) Then

          If (mode==0 .Or. mode==2) Then
            c(1) = x(1)**2 + x(2)**2 + x(3)**2 + x(4)**2
          End If

          If (mode==1 .Or. mode==2) Then
            cjac(1,1) = x(1) + x(1)
            cjac(1,2) = x(2) + x(2)
            cjac(1,3) = x(3) + x(3)
            cjac(1,4) = x(4) + x(4)
          End If

        End If

        If (needc(2)>0) Then

          If (mode==0 .Or. mode==2) Then
            c(2) = x(1)*x(2)*x(3)*x(4)
          End If

          If (mode==1 .Or. mode==2) Then
            cjac(2,1) = x(2)*x(3)*x(4)
            cjac(2,2) = x(1)*x(3)*x(4)
            cjac(2,3) = x(1)*x(2)*x(4)
            cjac(2,4) = x(1)*x(2)*x(3)
          End If

        End If

        Return

      End Subroutine confun
    End Module e04ucae_mod
    Program e04ucae

!     E04UCA Example Main Program

!     .. Use Statements ..
      Use e04ucae_mod, Only: confun, inc1, lcwsav, liwsav, llwsav, lrwsav,     &
                             nin, nout, objfun, one, zero
      Use nag_library, Only: dgemv, e04uca, e04wbf, nag_wp
!     .. Implicit None Statement ..
      Implicit None
!     .. Local Scalars ..
      Real (Kind=nag_wp)               :: objf
      Integer                          :: i, ifail, iter, j, lda, ldcj, ldr,   &
                                          liwork, lwork, n, nclin, ncnln, sda, &
                                          sdcjac
!     .. Local Arrays ..
      Real (Kind=nag_wp), Allocatable  :: a(:,:), bl(:), bu(:), c(:),          &
                                          cjac(:,:), clamda(:), objgrd(:),     &
                                          r(:,:), work(:), x(:)
      Real (Kind=nag_wp)               :: ruser(1), rwsav(lrwsav)
      Integer, Allocatable             :: istate(:), iwork(:)
      Integer                          :: iuser(1), iwsav(liwsav)
      Logical                          :: lwsav(llwsav)
      Character (80)                   :: cwsav(lcwsav)
!     .. Intrinsic Procedures ..
      Intrinsic                        :: max
!     .. Executable Statements ..
      Write (nout,*) 'E04UCA Example Program Results'

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

      Read (nin,*) n, nclin, ncnln
      liwork = 3*n + nclin + 2*ncnln
      lda = max(1,nclin)

      If (nclin>0) Then
        sda = n
      Else
        sda = 1
      End If

      ldcj = max(1,ncnln)

      If (ncnln>0) Then
        sdcjac = n
      Else
        sdcjac = 1
      End If

      ldr = n

      If (ncnln==0 .And. nclin>0) Then
        lwork = 2*n**2 + 20*n + 11*nclin
      Else If (ncnln>0 .And. nclin>=0) Then
        lwork = 2*n**2 + n*nclin + 2*n*ncnln + 20*n + 11*nclin + 21*ncnln
      Else
        lwork = 20*n
      End If

      Allocate (istate(n+nclin+ncnln),iwork(liwork),a(lda,sda),                &
        bl(n+nclin+ncnln),bu(n+nclin+ncnln),c(max(1,                           &
        ncnln)),cjac(ldcj,sdcjac),clamda(n+nclin+ncnln),objgrd(n),r(ldr,n),    &
        x(n),work(lwork))

      If (nclin>0) Then
        Read (nin,*)(a(i,1:sda),i=1,nclin)
      End If

      Read (nin,*) bl(1:(n+nclin+ncnln))
      Read (nin,*) bu(1:(n+nclin+ncnln))
      Read (nin,*) x(1:n)

!     Initialise E04UCA

      ifail = 0
      Call e04wbf('E04UCA',cwsav,lcwsav,lwsav,llwsav,iwsav,liwsav,rwsav,       &
        lrwsav,ifail)

!     Solve the problem

      ifail = -1
      Call e04uca(n,nclin,ncnln,lda,ldcj,ldr,a,bl,bu,confun,objfun,iter,       &
        istate,c,cjac,clamda,objf,objgrd,r,x,iwork,liwork,work,lwork,iuser,    &
        ruser,lwsav,iwsav,rwsav,ifail)

      Select Case (ifail)
      Case (0:6,8)
        Write (nout,*)
        Write (nout,99999)
        Write (nout,*)

        Do i = 1, n
          Write (nout,99998) i, istate(i), x(i), clamda(i)
        End Do

        If (nclin>0) Then

!         A*x --> work.
!         The NAG name equivalent of dgemv is f06paf
          Call dgemv('N',nclin,n,one,a,lda,x,inc1,zero,work,inc1)

          Write (nout,*)
          Write (nout,*)
          Write (nout,99997)
          Write (nout,*)

          Do i = n + 1, n + nclin
            j = i - n
            Write (nout,99996) j, istate(i), work(j), clamda(i)
          End Do

        End If

        If (ncnln>0) Then
          Write (nout,*)
          Write (nout,*)
          Write (nout,99995)
          Write (nout,*)

          Do i = n + nclin + 1, n + nclin + ncnln
            j = i - n - nclin
            Write (nout,99994) j, istate(i), c(j), clamda(i)
          End Do

        End If

        Write (nout,*)
        Write (nout,*)
        Write (nout,99993) objf
      End Select

99999 Format (1X,'Varbl',2X,'Istate',3X,'Value',9X,'Lagr Mult')
99998 Format (1X,'V',2(1X,I3),4X,1P,G14.6,2X,1P,G12.4)
99997 Format (1X,'L Con',2X,'Istate',3X,'Value',9X,'Lagr Mult')
99996 Format (1X,'L',2(1X,I3),4X,1P,G14.6,2X,1P,G12.4)
99995 Format (1X,'N Con',2X,'Istate',3X,'Value',9X,'Lagr Mult')
99994 Format (1X,'N',2(1X,I3),4X,1P,G14.6,2X,1P,G12.4)
99993 Format (1X,'Final objective value = ',1P,G15.7)
    End Program e04ucae