```!   E04UDA Example Program Text
!   Mark 26.1 Release. NAG Copyright 2016.

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
!     .. 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, iset = 1, lcwsav = 1,      &
liwsav = 610, llwsav = 120,          &
lrwsav = 475, nin = 5, ninopt = 7,   &
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 e04udae_mod
Program e04udae

!     E04UDA Example Main Program

!     .. Use Statements ..
Use e04udae_mod, Only: confun, inc1, iset, lcwsav, liwsav, llwsav,       &
lrwsav, nin, ninopt, nout, objfun, one, zero
Use nag_library, Only: dgemv, e04uca, e04uda, e04uea, e04wbf, nag_wp,    &
x04abf, x04acf, x04baf
!     .. 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 Procedures ..
Intrinsic                        :: max
!     .. Executable Statements ..
Write (rec,99991) 'E04UDA Example Program Results'
Call x04baf(nout,rec)

!     Skip heading in data file

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
End If

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

100   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
```