!   E04UDA Example Program Text
!   Mark 23 Release. NAG Copyright 2011.

    MODULE e04udae_mod

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

!      .. Use Statements ..
       USE nag_library, ONLY : nag_wp
!      .. Implicit None Statement ..
       IMPLICIT NONE
!      .. Parameters ..
       REAL (KIND=nag_wp), PARAMETER       :: one = 1.0_nag_wp
       REAL (KIND=nag_wp), PARAMETER       :: zero = 0.0_nag_wp
       INTEGER, PARAMETER                  :: inc1 = 1, iset = 1, lcwsav = 1,  &
                                              liwsav = 610, llwsav = 120,      &
                                              lrwsav = 475, nin = 5,           &
                                              ninopt = 7, nout = 6
    CONTAINS
!dec$ attributes offload : mic :: objfun
       SUBROUTINE objfun(mode,n,x,objf,objgrd,nstate,iuser,ruser)
!         Routine to evaluate objective function and its 1st derivatives.

!         .. Implicit None Statement ..
          IMPLICIT NONE
!         .. 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

!dec$ attributes offload : mic :: objfun
       END SUBROUTINE objfun
!dec$ attributes offload : mic :: confun
       SUBROUTINE confun(mode,ncnln,n,ldcj,needc,x,c,cjac,nstate,iuser,ruser)
!         Routine to evaluate the nonlinear constraints and their 1st
!         derivatives.

!         .. Implicit None Statement ..
          IMPLICIT NONE
!         .. 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

!dec$ attributes offload : mic :: confun
       END SUBROUTINE confun
    END MODULE e04udae_mod
    PROGRAM e04udae

!      E04UDA Example Main Program

!      .. Use Statements ..
       USE nag_library, ONLY : dgemv, e04uca, e04uda, e04uea, e04wbf, nag_wp,  &
                               x04abf, x04acf, x04baf
       USE e04udae_mod, ONLY : confun, inc1, iset, lcwsav, liwsav, llwsav,     &
                               lrwsav, nin, ninopt, nout, objfun, one, zero
!      .. Implicit None Statement ..
       IMPLICIT NONE
!      .. Parameters ..
       CHARACTER (*), PARAMETER            :: fname = 'e04udae.opt'
!      .. Local Scalars ..
       REAL (KIND=nag_wp)                  :: objf
       INTEGER                             :: i, ifail, inform, iter, j, lda,  &
                                              ldcj, ldr, liwork, lwork, mode,  &
                                              n, nclin, ncnln, outchn, sda,    &
                                              sdcjac
       CHARACTER (80)                      :: rec
!      .. 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 Functions ..
       INTRINSIC                              max
!      .. Executable Statements ..
       WRITE (rec,99991) 'E04UDA Example Program Results'
       CALL x04baf(nout,rec)

!      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)

!      Set the unit number for advisory messages to OUTCHN

       outchn = nout
       CALL x04abf(iset,outchn)

!      Initialise E04UCA

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

!      Set two options using E04UEA

       CALL e04uea(' Infinite Bound Size = 1.0D+25 ',lwsav,iwsav,rwsav,inform)

       IF (inform==0) THEN

          CALL e04uea(' Verify Level = -1 ',lwsav,iwsav,rwsav,inform)

       END IF

       IF (inform/=0) THEN
          WRITE (rec,99999) 'E04UEA terminated with ' // 'INFORM = ', inform
          CALL x04baf(nout,rec)
          GO TO 20
       END IF

!      Open the options file for reading

       mode = 0

       ifail = 0
       CALL x04acf(ninopt,fname,mode,ifail)

!      Read the options file for the remaining options

       CALL e04uda(ninopt,lwsav,iwsav,rwsav,inform)

       IF (inform/=0) THEN
          WRITE (rec,99999) 'E04UDA terminated with ' // 'INFORM = ', inform
          CALL x04baf(nout,rec)
       END IF

!      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 (rec,'()')
          CALL x04baf(nout,rec)
          WRITE (rec,99998)
          CALL x04baf(nout,rec)
          WRITE (rec,'()')
          CALL x04baf(nout,rec)

          DO i = 1, n
             WRITE (rec,99997) i, istate(i), x(i), clamda(i)
             CALL x04baf(nout,rec)
          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 (rec,'()')
             CALL x04baf(nout,rec)
             WRITE (rec,'()')
             CALL x04baf(nout,rec)
             WRITE (rec,99996)
             CALL x04baf(nout,rec)
             WRITE (rec,'()')
             CALL x04baf(nout,rec)

             DO i = n + 1, n + nclin
                j = i - n
                WRITE (rec,99995) j, istate(i), work(j), clamda(i)
                CALL x04baf(nout,rec)
             END DO

          END IF

          IF (ncnln>0) THEN
             WRITE (rec,'()')
             CALL x04baf(nout,rec)
             WRITE (rec,'()')
             CALL x04baf(nout,rec)
             WRITE (rec,99994)
             CALL x04baf(nout,rec)
             WRITE (rec,'()')
             CALL x04baf(nout,rec)

             DO i = n + nclin + 1, n + nclin + ncnln
                j = i - n - nclin
                WRITE (rec,99993) j, istate(i), c(j), clamda(i)
                CALL x04baf(nout,rec)
             END DO

          END IF

          WRITE (rec,'()')
          CALL x04baf(nout,rec)
          WRITE (rec,'()')
          CALL x04baf(nout,rec)
          WRITE (rec,99992) objf
          CALL x04baf(nout,rec)
       END SELECT

20     CONTINUE

99999  FORMAT (1X,A,I5)
99998  FORMAT (1X,'Varbl',2X,'Istate',3X,'Value',9X,'Lagr Mult')
99997  FORMAT (1X,'V',2(1X,I3),4X,1P,G14.6,2X,1P,G12.4)
99996  FORMAT (1X,'L Con',2X,'Istate',3X,'Value',9X,'Lagr Mult')
99995  FORMAT (1X,'L',2(1X,I3),4X,1P,G14.6,2X,1P,G12.4)
99994  FORMAT (1X,'N Con',2X,'Istate',3X,'Value',9X,'Lagr Mult')
99993  FORMAT (1X,'N',2(1X,I3),4X,1P,G14.6,2X,1P,G12.4)
99992  FORMAT (1X,'Final objective value = ',1P,G15.7)
99991  FORMAT (1X,A)
    END PROGRAM e04udae
