# Try out NAG Library functions

Explore NAG maths and stats routines with interactive demos
Function ID
E05USF
Name
nagf_glopt_nlp_multistart_sqp_lsq
Description
Global optimization of a sum of squares problem using multi-start, nonlinear constraints
Keywords
multi-start algorithm | sequential QP method | sum of squares
This example is based on Problem 57 in Hock and Schittkowski (1981) and involves the minimization of the sum of squares function
 $Fx=12∑i=144yi-fix2,$
where
 $fix=x1+0.49-x1e-x2ai-8$
and
 $iyiaiiyiai10.498230.412220.498240.402230.4810250.422440.4710260.402450.4810270.402460.4710280.412670.4612290.402680.4612300.412690.4512310.4128100.4312320.4028110.4514330.4030120.4314340.4030130.4314350.3830140.4416360.4132150.4316370.4032160.4316380.4034170.4618390.4136180.4518400.3836190.4220410.4038200.4220420.4038210.4320430.3940220.4122440.3942$
subject to the bounds
 $x1≥-0.4x2≥-4.0$
to the general linear constraint
 $x1+x2≥1.0$
and to the nonlinear constraint
 $0.49x2-x1x2≥0.09.$
!   E05USF Example Program Text
!   Mark 26.1 Release. NAG Copyright 2016.

Module e05usfe_mod

!     E05USF 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, start
Contains
Subroutine objfun(mode,m,n,ldfjsl,needfi,x,f,fjsl,nstate,iuser,ruser)
!       Evaluates the subfunctions and their 1st derivatives.

!       .. Parameters ..
Real (Kind=nag_wp), Parameter  :: a(44) = (/8._nag_wp,8._nag_wp,       &
10._nag_wp,10._nag_wp,10._nag_wp,    &
10._nag_wp,12._nag_wp,12._nag_wp,    &
12._nag_wp,12._nag_wp,14._nag_wp,    &
14._nag_wp,14._nag_wp,16._nag_wp,    &
16._nag_wp,16._nag_wp,18._nag_wp,    &
18._nag_wp,20._nag_wp,20._nag_wp,    &
20._nag_wp,22._nag_wp,22._nag_wp,    &
22._nag_wp,24._nag_wp,24._nag_wp,    &
24._nag_wp,26._nag_wp,26._nag_wp,    &
26._nag_wp,28._nag_wp,28._nag_wp,    &
30._nag_wp,30._nag_wp,30._nag_wp,    &
32._nag_wp,32._nag_wp,34._nag_wp,    &
36._nag_wp,36._nag_wp,38._nag_wp,    &
38._nag_wp,40._nag_wp,42._nag_wp/)
!       .. Scalar Arguments ..
Integer, Intent (In)           :: ldfjsl, m, n, needfi, nstate
Integer, Intent (Inout)        :: mode
!       .. Array Arguments ..
Real (Kind=nag_wp), Intent (Out) :: f(m)
Real (Kind=nag_wp), Intent (Inout) :: fjsl(ldfjsl,*), ruser(*)
Real (Kind=nag_wp), Intent (In) :: x(n)
Integer, Intent (Inout)        :: iuser(*)
!       .. Local Scalars ..
Real (Kind=nag_wp)             :: temp
Integer                        :: i
!       .. Intrinsic Procedures ..
Intrinsic                      :: exp
!       .. Executable Statements ..

!       This is a two-dimensional objective function.
!       As an example of using the mode mechanism,
!       terminate if any other problem size is supplied.

If (n/=2) Then
mode = -1
End If

If (nstate==1) Then
!         This is the first call.
!         Take any special action here if desired.
Continue
End If

If (mode==0 .And. needfi>0) Then
f(needfi) = x(1) + (0.49_nag_wp-x(1))*exp(-x(2)*(a(needfi)-          &
8.0_nag_wp))
Else
Do i = 1, m
temp = exp(-x(2)*(a(i)-8._nag_wp))

If (mode==0 .Or. mode==2) Then
f(i) = x(1) + (0.49_nag_wp-x(1))*temp
End If

If (mode==1 .Or. mode==2) Then
fjsl(i,1) = 1._nag_wp - temp
fjsl(i,2) = -(0.49_nag_wp-x(1))*(a(i)-8._nag_wp)*temp
End If

End Do
End If

Return
End Subroutine objfun
Subroutine confun(mode,ncnln,n,ldcjsl,needc,x,c,cjsl,nstate,iuser,ruser)
!       Evaluates the nonlinear constraints and their 1st derivatives.

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

!       This problem has only one constraint.
!       As an example of using the mode mechanism,
!       terminate if any other size is supplied.

If (ncnln/=1) Then
mode = -1
End If

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.1 of the E04USA document).

cjsl(1:ncnln,1:n) = 0._nag_wp
End If

If (needc(1)>0) Then

If (mode==0 .Or. mode==2) Then
c(1) = -0.09_nag_wp - x(1)*x(2) + 0.49_nag_wp*x(2)
End If

If (mode==1 .Or. mode==2) Then
cjsl(1,1) = -x(2)
cjsl(1,2) = -x(1) + 0.49_nag_wp
End If

End If

Return
End Subroutine confun
Subroutine start(npts,quas,n,repeat1,bl,bu,iuser,ruser,mode)

!       Sets the initial points.
!       A typical user-defined start procedure.

!       .. Use Statements ..
Use nag_library, Only: g05kgf, g05saf
!       .. Scalar Arguments ..
Integer, Intent (Inout)        :: mode
Integer, Intent (In)           :: n, npts
Logical, Intent (In)           :: repeat1
!       .. Array Arguments ..
Real (Kind=nag_wp), Intent (In) :: bl(n), bu(n)
Real (Kind=nag_wp), Intent (Inout) :: quas(n,npts), ruser(*)
Integer, Intent (Inout)        :: iuser(*)
!       .. Local Scalars ..
Integer                        :: genid, i, ifail, lstate, subid
!       .. Local Arrays ..
Integer, Allocatable           :: state(:)
!       .. Executable Statements ..
!       quas is pre-assigned to zero.
If (repeat1) Then
quas(1,1) = 0.4_nag_wp
quas(2,2) = 1._nag_wp
Else
!         Generate a non-repeatable spread of points between bl and bu.
genid = 2
subid = 53
lstate = -1
Allocate (state(lstate))
ifail = 0
Call g05kgf(genid,subid,state,lstate,ifail)
Deallocate (state)
Allocate (state(lstate))
ifail = 0
Call g05kgf(genid,subid,state,lstate,ifail)
Do i = 1, npts
ifail = 0
Call g05saf(n,state,quas(1,i),ifail)
quas(1:n,i) = bl(1:n) + (bu(1:n)-bl(1:n))*quas(1:n,i)
End Do
Deallocate (state)
End If
!       Set mode negative to terminate execution for any reason.
mode = 0
Return
End Subroutine start
End Module e05usfe_mod
Program e05usfe

!     E05USF Example Main Program

!     .. Use Statements ..
Use e05usfe_mod, Only: confun, objfun, start
Use nag_library, Only: dgemv, e05usf, e05zkf, nag_wp
!     .. Implicit None Statement ..
Implicit None
!     .. Parameters ..
Integer, Parameter               :: liopts = 740, lopts = 485, m = 44,   &
n = 2, nb = 1, nclin = 1, ncnln = 1, &
nin = 5, nout = 6, npts = 3
Integer, Parameter               :: sdfjac = n
Integer, Parameter               :: lda = nclin
Integer, Parameter               :: ldc = ncnln
Integer, Parameter               :: ldcjac = ncnln
Integer, Parameter               :: ldclda = n + nclin + ncnln
Integer, Parameter               :: ldfjac = m
Integer, Parameter               :: ldx = n
Integer, Parameter               :: listat = n + nclin + ncnln
Logical, Parameter               :: repeat1 = .True.
!     .. Local Scalars ..
Integer                          :: i, ifail, j, k, l, sda, sdcjac
!     .. Local Arrays ..
Real (Kind=nag_wp), Allocatable  :: a(:,:), bl(:), bu(:), c(:,:),        &
cjac(:,:,:), clamda(:,:), f(:,:),    &
fjac(:,:,:), work(:), x(:,:), y(:)
Real (Kind=nag_wp)               :: objf(nb), opts(lopts), ruser(1)
Integer                          :: info(nb), iopts(liopts), iter(nb),   &
iuser(1)
Integer, Allocatable             :: istate(:,:)
!     .. Intrinsic Procedures ..
Intrinsic                        :: max
!     .. Executable Statements ..
Write (nout,*) 'E05USF Example Program Results'
Flush (nout)

!     Skip heading in data file.

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

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

Allocate (a(lda,sda),bl(n+nclin+ncnln),bu(n+nclin+ncnln),y(m),c(ldc,nb), &
cjac(ldcjac,sdcjac,nb),f(m,nb),fjac(ldfjac,sdfjac,nb),                 &
clamda(ldclda,nb),istate(listat,nb),x(ldx,nb),work(max(1,nclin)))

If (nclin>0) Then
End If

!     Initialize the solver.

ifail = 0
Call e05zkf('Initialize = E05USF',iopts,liopts,opts,lopts,ifail)

!     Solve the problem.

ifail = -1
Call e05usf(m,n,nclin,ncnln,a,lda,bl,bu,y,confun,objfun,npts,x,ldx,      &
start,repeat1,nb,objf,f,fjac,ldfjac,sdfjac,iter,c,ldc,cjac,ldcjac,     &
sdcjac,clamda,ldclda,istate,listat,iopts,opts,iuser,ruser,info,ifail)

Select Case (ifail)
Case (0)
l = nb
Case (8)
l = info(nb)
Write (nout,99999) iter(nb)
Case Default
Go To 100
End Select

loop: Do i = 1, l
Write (nout,99998) i
Write (nout,99997) info(i)
Write (nout,99996) 'Varbl'
Do j = 1, n
Write (nout,99995) 'V', j, istate(j,i), x(j,i), clamda(j,i)
End Do
If (nclin>0) Then
Write (nout,99996) 'L Con'

!         Below is a call to the level 2 BLAS routine DGEMV.
!         This performs the matrix vector multiplication A*X
!         (linear constraint values) and puts the result in
!         the first NCLIN locations of WORK.

Call dgemv('N',nclin,n,1.0_nag_wp,a,lda,x(1,i),1,0.0_nag_wp,work,1)

Do k = n + 1, n + nclin
j = k - n
Write (nout,99995) 'L', j, istate(k,i), work(j), clamda(k,i)
End Do
End If
If (ncnln>0) Then
Write (nout,99996) 'N Con'
Do k = n + nclin + 1, n + nclin + ncnln
j = k - n - nclin
Write (nout,99995) 'N', j, istate(k,i), c(j,i), clamda(k,i)
End Do
End If
Write (nout,99994) objf(i)
Write (nout,99993)
Write (nout,99992)(clamda(k,i),k=1,n+nclin+ncnln)

If (l==1) Then
Exit loop
End If

Write (nout,*)

Write (nout,*)                                                         &
' ------------------------------------------------------ '

End Do loop

100   Continue

99999 Format (1X,I20,'starting points converged')
99998 Format (/,1X,'Solution number',I16)
99997 Format (/,1X,'Local minimization exited with code',I5)
99996 Format (/,1X,A,2X,'Istate',3X,'Value',9X,'Lagr Mult',/)
99995 Format (1X,A,2(1X,I3),4X,F12.4,2X,F12.4)
99994 Format (/,1X,'Final objective value = ',1X,F12.4)
99993 Format (/,1X,'QP multipliers')
99992 Format (1X,F12.4)
End Program e05usfe

The NAG Library
The world’s largest collection of robust, documented, tested and maintained numerical algorithms.