! Mark 24 Release. NAG Copyright 2012. Module e05ucfe_mod ! E05UCF Example Program Module: ! Parameters and User-defined Routines ! .. Use Statements .. Use nag_library, Only: nag_wp ! .. Implicit None Statement .. Implicit None Contains Subroutine schwefel_obj(mode,n,x,objf,objgrd,nstate,iuser,ruser) ! .. 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(*) ! .. Intrinsic Procedures .. Intrinsic :: abs, cos, sin, sqrt, sum ! .. Executable Statements .. If (mode==0 .Or. mode==2) Then ! Evaluate the objective function. objf = sum(x(1:n)*sin(sqrt(abs(x(1:n))))) End If If (mode==1 .Or. mode==2) Then ! Calculate the gradient of the objective function. objgrd(1:n) = sin(sqrt(abs(x(1:n)))) + 0.5_nag_wp*sqrt(abs(x(1:n)))* & cos(sqrt(abs(x(1:n)))) End If Return End Subroutine schwefel_obj Subroutine schwefel_confun(mode,ncnln,n,ldcjsl,needc,x,c,cjsl,nstate, & iuser,ruser) ! .. 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,n), ruser(*) Real (Kind=nag_wp), Intent (In) :: x(n) Integer, Intent (Inout) :: iuser(*) Integer, Intent (In) :: needc(ncnln) ! .. Local Scalars .. Real (Kind=nag_wp) :: t1, t2 Integer :: k ! .. Intrinsic Procedures .. Intrinsic :: cos, sin ! .. Executable Statements .. If (mode==0 .Or. mode==2) Then ! Constraint values are required. ! Only those for which needc is non-zero need be set. Do k = 1, ncnln If (needc(k)>0) Then Select Case (k) Case (1) c(k) = x(1)**2 - x(2)**2 + 3.0_nag_wp*x(1)*x(2) Case (2) c(k) = cos((x(1)/200.0_nag_wp)**2+(x(2)/100.0_nag_wp)) End Select End If End Do End If If (mode==1 .Or. mode==2) Then ! Constraint derivatives are required. Do k = 1, ncnln Select Case (k) Case (1) cjsl(k,1) = 2.0_nag_wp*x(1) + 3.0_nag_wp*x(2) cjsl(k,2) = -2.0_nag_wp*x(2) + 3.0_nag_wp*x(1) Case (2) t1 = x(1)/200.0_nag_wp t2 = x(2)/100.0_nag_wp cjsl(k,1) = -sin(t1**2+t2)*2.0_nag_wp*t1/200.0_nag_wp cjsl(k,2) = -sin(t1**2+t2)/100.0_nag_wp End Select End Do End If Return End Subroutine schwefel_confun Subroutine mystart(npts,quas,n,repeat,bl,bu,iuser,ruser,mode) ! Sets an initial point. ! A typical user-defined start procedure (not used here) ! .. Scalar Arguments .. Integer, Intent (Inout) :: mode Integer, Intent (In) :: n, npts Logical, Intent (In) :: repeat ! .. 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(*) ! .. Executable Statements .. ! quas(1:n,1:npts) is pre-assigned to zero. quas(1,1) = 420.9687_nag_wp quas(2,2) = 1.0_nag_wp !etc Return End Subroutine mystart End Module e05ucfe_mod Program e05ucfe ! E05UCF Example Main Program ! .. Use Statements .. Use nag_library, Only: dgemv, e05ucf, e05ucz, e05zkf, nag_wp Use e05ucfe_mod, Only: schwefel_confun, schwefel_obj ! .. Implicit None Statement .. Implicit None ! .. Parameters .. Integer, Parameter :: n = 2, nclin = 1, ncnln = 2, & nin = 5, nout = 6 ! .. Local Scalars .. Integer :: i, ifail, j, k, l, lda, ldc, & ldcjac, ldclda, ldobjd, ldr, & ldx, liopts, listat, lopts, nb, & npts, sda, sdcjac, sdr Logical :: repeat ! .. Local Arrays .. Real (Kind=nag_wp), Allocatable :: a(:,:), bl(:), bu(:), c(:,:), & cjac(:,:,:), clamda(:,:), & objf(:), objgrd(:,:), opts(:), & r(:,:,:), work(:), x(:,:) Real (Kind=nag_wp) :: ruser(1) Integer, Allocatable :: info(:), iopts(:), istate(:,:), & iter(:) Integer :: iuser(1) ! .. Executable Statements .. Write (nout,*) 'E05UCF Example Program Results' Flush (nout) ! Skip heading in data file Read (nin,*) Read (nin,*) nb, npts Read (nin,*) repeat lda = nclin If (nclin>0) Then sda = n Else sda = 1 End If ldx = n ldobjd = n ldc = ncnln ldcjac = ncnln If (ncnln>0) Then sdcjac = n Else sdcjac = 0 End If ldr = n sdr = n ldclda = n + nclin + ncnln listat = n + nclin + ncnln liopts = 740 lopts = 485 Allocate (a(lda,sda),bl(n+nclin+ncnln),bu(n+nclin+ncnln),x(ldx,nb), & objf(nb),objgrd(ldobjd,nb),iter(nb),c(ldc,nb),cjac(ldcjac,sdcjac,nb), & r(ldr,sdr,nb),clamda(ldclda,nb),istate(listat,nb),iopts(liopts), & opts(lopts),info(nb),work(nclin)) bl(1:n+nclin+ncnln) = (/-500.0_nag_wp,-500.0_nag_wp,-10000.0_nag_wp, & -1.0_nag_wp,-0.9_nag_wp/) bu(1:n+nclin+ncnln) = (/500.0_nag_wp,500.0_nag_wp,10.0_nag_wp, & 500000.0_nag_wp,0.9_nag_wp/) a(1,1) = 3.0_nag_wp a(1,2) = -2.0_nag_wp ! Initialize the solver. ifail = 0 Call e05zkf('Initialize = E05UCF',iopts,liopts,opts,lopts,ifail) ! Solve the problem. ! USE mystart from e05ucfe_mod and pass this instead of e05ucz to e05ucf ! to set your own start points. ifail = -1 Call e05ucf(n,nclin,ncnln,a,lda,bl,bu,schwefel_confun,schwefel_obj,npts, & x,ldx,e05ucz,repeat,nb,objf,objgrd,ldobjd,iter,c,ldc,cjac,ldcjac, & sdcjac,r,ldr,sdr,clamda,ldclda,istate,listat,iopts,opts,iuser,ruser, & info,ifail) Select Case (ifail) Case (0) l = nb Case (8) l = info(nb) Write (nout,'(1X,I16,A)') iter(nb), ' starting points converged' Case Default Go To 100 End Select loop: Do i = 1, l Write (nout,99999) i Write (nout,99998) info(i) Write (nout,99997) 'Varbl' Do j = 1, n Write (nout,99996) 'V', j, istate(j,i), x(j,i), clamda(j,i) End Do If (nclin>0) Then Write (nout,99997) '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,99996) 'L', j, istate(k,i), work(j), clamda(k,i) End Do End If If (ncnln>0) Then Write (nout,99997) 'N Con' Do k = n + nclin + 1, n + nclin + ncnln j = k - n - nclin Write (nout,99996) 'N', j, istate(k,i), c(j,i), clamda(k,i) End Do End If Write (nout,99995) objf(i) Write (nout,99994) Write (nout,99993)(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,'Solution number',I16) 99998 Format (/1X,'Local minimization exited with code',I5) 99997 Format (/1X,A,2X,'Istate',3X,'Value',9X,'Lagr Mult'/) 99996 Format (1X,A,2(1X,I3),4X,F12.4,2X,F12.4) 99995 Format (/1X,'Final objective value = ',1X,F12.4) 99994 Format (/1X,'QP multipliers') 99993 Format (1X,F12.4) End Program e05ucfe