C ================================================================= C File: algencanma.f C ================================================================= C ================================================================= C Module: Main program C ================================================================= C Last update of any of the component of this module: C C February 14, 2006. C Users are encouraged to download periodically updated versions of C this code at the TANGO home page: C C www.ime.usp.br/~egbirgin/tango/ C ***************************************************************** C ***************************************************************** C ALGENCAN solves problems of the form C ------------------------------------ C C min f(x) C C subject to C C c_j(x) = 0, j in E, C c_j(x) <= 0, j in I, C l <= x <= u, C C where E is the set of indices of the equality constraints, I is C the set of indices of the inequality constraints, and there are C n variables and m constraints. C C ALGENCAN is an Augmented Lagrangian method that uses GENCAN to C solve the bound-constrained problems. C C ALGENCAN is part of the TANGO Project. C C Visit the TANGO home page in the web: C C www.ime.usp.br/~egbirgin/tango/ C ***************************************************************** C TANGO LICENSE: C -------------- C C The codes of TANGO are free for non-commercial use. They may be C distributed within a "Research Group" but not distributed to C third parties without the authorization of E. G. Birgin or J. M. C Mart�ez. The code must not be included in any commercial package C without previous agreement. C C Every published work that uses ALGENCAN must cite: C C R. Andreani, E. G. Birgin, J. M. Mart�ez and M. L. Schuverdt, C "On Augmented Lagrangian methods with general lower-level C constraints", Technical Report MCDO-050304, Department of Applied C Mathematics, UNICAMP, Brazil, 2005 C C and C C R. Andreani, E. G. Birgin, J. M. Mart�ez and M. L. Schuverdt, C "Augmented Lagrangian methods under the Constant Positive Linear C Dependence constraint qualification", to appear in Mathematical C Programming, 2005. C C Every published work that uses GENCAN must cite: C C E. G. Birgin and J. M. Mart�ez, "Large-scale active-set C box-constrained optimization method with spectral projected C gradients", Computational Optimization and Applications 23, pp. C 101-125, 2002. C C (See other related works at the TANGO home page.) C ***************************************************************** C HOW TO USE ALGENCAN TO SOLVE YOUR OWN PROBLEM C -------------------------------------------- C C You will find below 10 subroutines and that YOU SHOULD MODIFY to C solve your own problem. The modifications that you must do are: C C 1) In the SUBROUTINE INIP you must set the number of variables C of your problem (n), the initial point (x), the lower and upper C bounds on the variables (l and u), the number of constraints (m) C (set m equal to zero if your problems has no constraints), a C vector (equatn) that indicates, for each constraint, whether it C is and equality constraint or an inequality (least-or-equal), C a vector (linear) that indicates, for each constraint, whether it C is a linear constraint or not, the initial estimation of the C Lagrange multipliers (lambda) (you can set it equal to zero if C you do not have a better estimative) and the initial penalty C parameters (rho). C C 2) SUBROUTINES EVALF and EVALC MUST also be MODIFIED to evaluate C the objective function and the constraints, respectively. You C will find the way in which these subroutines must be modified C below. C C 3) SUBROUTINES EVALG and EVALJAC, to compute the gradient vector C of the objective function and the gradients of the constraints, C respectively, are NOT MANDATORY but HIGHLY RECOMMENDED. They C must be provided by the user if he/she set gtype = 0. If the C user set gtype = 1 then finite differences approximations will C be used. (See the description of gtype parameter in subroutine C param.) You will find the way in which these subroutines must C be modified below. C C 4) OPTIONALY, you can modify SUBROUTINES EVALH and EVALHC to C evaluate the Hessian matrices of the objective function and the C constraints. Both subroutines are just optional subroutines and C are required when hptype = 0. (See the description of hptype C parameter in subroutine param.) C C 5) OPTIONALY, another way of using second-order information is C to provide, not individual subroutines to compute the Hessians C of the objective function and the constraints, but to provide C a unique SUBROUTINE EVALHLP to compute the product of an C arbitrary vector times the Lagrangian Hessian matrix. There is C not a clear advantage in coding EVALHLP instead of coding C EVALH and EVALHC. This alternative was incorporated to be C primarily used in the AMPL and CUTEr interfaces. This C subroutine is optional and is required when hptype = 1. (See C the description of hptype parameter in subroutine param.) C C 6) OPTIONALLY, in the SUBROUTINE PARAM you may modify tolerances C and maximum of iterations related to the stopping criteria. There C are also two arguments related to the output in the screen (iprint C and ncomp), a parameter related to the derivatives calculation C (gtype), a parameter related to the Hessians (or matrix-vector C product) approximation (hptype) and a parameter related to the C Conjugates Gradients preconditioning (precond). Finally, there is C also a parameter that allows to test the analytic derivatives C coded by you using finite differences. If you prefer, leave all C these parameters with the suggested values. C C 7) Alternatively, if you are interested in modifying any of the C parameters described in (3) but you prefer not to modify subroutine C PARAM, you can use the algencan.dat SPECIFICATION FILE. The C advantage of using the specification file is that you can test C the method with different parameters without compiling it every C time. This file IS NOT MANDATORY. So, if you will not modify the C default parameters or if you fill more comfortable modifying C subroutine PARAM, you do not need to worry about algencan.dat C specification file. C C The specification file can have any number of lines, with at most C 80 characters per line. Blank lines are admitted and lines starting C with '*' are considered comment lines and will be ignored. Each C line must have just a "keyword" or a keyword followed by an integer C or a real value when required. As the interpreter is not case- C sensitive, you can write the keywords in lower case, upper case or C any combination of them. You can also insert blanks in any place C before or after the keywords. Moreover, just the first 10 characters C of the keyword are relevant and the rest will be ignored. Available C keywords are: C C ANALYTIC-GRADIENT C FINITE-DIFFERENCES-GRADIENT C HESSIANS-PROVIDED C LAGRHESS-PRODUCT-PROVIDED C INCREMENTAL-QUOTIENTS C BFGS-QN-APPROXIMATION C ADAPTIVE-HESSIAN C UNPRECONDITIONED-CG C BFGS-QN-PRECONDITIONER C CHECK-DERIVATIVES C FEASIBILITY-TOLERANCE C OPTIMALITY-TOLERANCE C MAX-OUTER-ITERATIONS C MAX-INNER-ITERATIONS C MAX-FUNCTION-EVALUATIONS C OUTPUT-DETAIL C NCOMP-ARRAY C C By default, ALGENCAN uses: C C ANALYTIC-GRADIENT C ADAPTIVE-HESSIAN C BFGS-QN-PRECONDITIONER C FEASIBILITY-TOLERANCE 1.0d-04 C OPTIMALITY-TOLERANCE 1.0d-04 C MAX-OUTER-ITERATIONS 50 C MAX-INNER-ITERATIONS 1000000 C MAX-FUNCTION-EVALUATIONS 5000000 C OUTPUT-DETAIL 1 C NCOMP-ARRAY 5 C C and derivatives are not tested. C C 8) Finally, SUBROUTINE ENDP can be modified by the user to write, C save or compute any type of extra information related to the C solution of the problem. Subroutine ENDP will be called just once C after algencan found a solution of the problem. You can modify C this subroutine for, for example, save the solution in a file C chosed by you and in a special format to be processed later by C another software. C C Besides the modifications that you must do to solve your own C problem using ALGENCAN, ALGENCAN has many other arguments to C control a large variety of issues of the algorithmic behaviour. C To see how to control these other arguments, see their detailed C description at the begining of subroutine algencan (for arguments C related to Augmented Lagrangians) and at the begining of C subroutine gencan (for arguments related to the bound-constraints C internal solver). C ****************************************************************** C ****************************************************************** program algencanma implicit none C ****************************************************************** C HERE WE SET THE MAXIMUM NUMBER OF VARIABLES (NMAX) AND THE MAXIMUM C NUMBER OF CONSTRAINTS (MMAX) THAT THE PROBLEM IS ALLOWED TO HAVE. C IF YOUR PROBLEM HAS MORE THAN THE DECLARED VALUES FOR NMAX AND C MMAX THEN YOU MUST MODIFY THESE VALUES IN SUCH A WAY THAT THE C NUMBER OF VARIABLES (N) OF YOUR PROBLEM BE SMALLER THAN OR EQUAL C TO NMAX AND THE NUMBER OF CONSTRAINTS (M) OF YOUR PROBLEM BE C SMALLER THAN OR EQUAL TO MMAX. C C ATTENTION: PARAMETERS NMAX AND MMAX ARE DECLARED IN FIVE PLACES, C THE MAIN PROGRAM AND SUBROUTINES EVALNAL, EVALHD, COMPP AND COMPH. C IF YOU CHANGE THEIR VALUES HERE THEN YOU MUST ALSO DO THE SAME C MODIFICATION IN THE OTHER PLACES TOO. THE VALUES OF NMAX AND MMAX C MUST BE THE SAME IN ALL PLACES! C ****************************************************************** C PARAMETERS integer mmax,nmax,samplemax,maxbadsample,MaxIter double precision INFINITO,ZERO,deltaX parameter ( MaxIter = 2000 ) parameter ( mmax = 1000000 ) parameter ( nmax = 1000000 ) parameter ( INFINITO = 1.0E+16 ) parameter ( ZERO = 1.0E-8 ) parameter ( deltaX = 1.0E-4 ) parameter ( samplemax = 10000000 ) parameter ( maxbadsample = 1000 ) C LOCAL SCALARS logical checkder,indsplpt,indreg,indRandomLinkage, + indPiecewiseLinearThreshould,indLissajous character * 6 precond character * 10 method integer gtype,hptype,inform,iprint,m,maxoutit,maxtotfc,maxtotit,n, + ncomp,outiter,totcgcnt,totfcnt,totgcnt,totiter, + ktrial,nsamples,err,fat,i,isample,itsample, + totmin,btrial,maxtrysample,trysample,badsample, + consbadsample,ltrial,totltrial,bltrial,btunnit,tunnit, + totevalf,totevalg,bevalf,bevalg,bnsamples,badltrial, + badtunnit,bbadltrial,bbadtunnit,maxbadchoice, + bfound double precision epsfeas,epsopt,f,nalpsupn,snorm,drand,seed,rk, + cpi,barea,sigmaf,gammavle,fmin,psearch, + fNormaVetor,varphi_k,vdummy real time,tottime,btime,etime,ctime C LOCAL ARRAYS real tarray(2) integer wi1(nmax) logical equatn(mmax),linear(mmax) double precision l(nmax),lambda(mmax),rho(mmax),u(nmax),wd1(mmax), + wd2(mmax),wd3(nmax),wd4(mmax),wd5(mmax),wd6(nmax), + wd7(nmax),wd8(nmax),wd9(nmax),wd10(nmax),wd11(nmax), + wd12(nmax),wd13(nmax),wd14(nmax),wd15(nmax),wd16(nmax), + wd17(nmax),wd18(nmax),wd19(nmax),x(nmax), + vetx(samplemax),vetf(samplemax),vAux(nmax), + xmin(samplemax) + C EXTERNAL SUBROUTINES external solver C SET UP PROBLEM DATA call inip(n,x,l,u,m,lambda,rho,equatn,linear) C SET METHOD CHOICE C set indRandomLinkage to .true. to enable Random Linkage method indRandomLinkage = .false. C if random linkage is enabled set the sigma vle sigmaf = 2.0d0 C set indLissajous to .true. to enable Tunnel with Lis method indLissajous = .true. C experimental only - not used on tests indPiecewiseLinearThreshould = .false. if( indPiecewiseLinearThreshould ) then method = 'LTRLINKAGE' ELSE IF( indRandomLinkage ) then method = 'RNDLINKAGE' else if( indLissajous ) then method = 'LIS TUNNEL' else method = 'MULTISTART' END IF sigmaf = 0.0d0 end if END IF write(*,*) " " write(*,*) " " write(*,0020) write(*,*) " " write(*,0019) method,sigmaf write(*,0021) write(*,0022) C SET LOCAL SCALARS maxtrysample = 50 nsamples = 0 fmin = INFINITO seed = 0.0d0 totmin = 0 badsample = 0 consbadsample = 0 totltrial = 0 tunnit = 0 badltrial = 0 badtunnit = 0 totevalf = 0 totevalg = 0 maxbadchoice = 5 bfound = 0 C START TIME COUNTER tottime = etime(tarray) ctime = tottime C SET PI VALUE cpi = 2.0d0 * acos(0.0d0) C SET BOX AREA barea = u(1) - l(1) do i = 2,n barea = barea * ( u(i) - l(i) ) end do C SET THE GAMMA FUNCTION VALUE FOR: 1 + n/2 gammavle = 1.0d0 C FOR n EVEN: if ( mod(n,2) .eq. 0 ) then fat = n / 2 do i = 2,fat gammavle = gammavle * i end do C FOR n ODD: else fat = ( n + 1 ) / 2 do i = fat,n gammavle = gammavle * i end do gammavle = gammavle * sqrt( cpi ) / 2.0d0 ** n end if C SET SOME SOLVER ARGUMENTS call param(gtype,hptype,precond,checkder,epsfeas,epsopt,maxoutit, +maxtotit,maxtotfc,iprint,ncomp) C PREPARE SOLUTION OUTPUT FILE C open(11,file='fbest.txt') C ITERATE ktrial = 1 f = INFINITO do while( ( ktrial .le. MaxIter ) .and. + ( consbadsample .le. maxbadsample ) ) C SAMPLE A NEW POINT C RANDOM LINKAGE METHOD itsample = 1 if( indRandomLinkage ) then indsplpt = .true. trysample = 0 do while( indsplpt ) C CHECK THE INITIAL PTS VECTOR SIZE if( ( nsamples * n ) .ge. ( samplemax ) ) then write(*,*) 'END OF INITIAL PTS VECTOR!!!' go to 650 end if seed = nsamples + 1 vdummy = drand(seed) do i = 1,n x(i) = l(i) + ( u(i) - l(i) ) * drand(seed) vetx( nsamples * n + i ) = x(i) end do nsamples = nsamples + 1 trysample = trysample + 1 C SET FUNCTION VALUE call evalf(n,x,f,err) vetf(nsamples) = f totevalf = totevalf + 1 C SET THE THRESHOLD rk rk = ( 1.0d0 / sqrt( cpi ) ) * + ( gammavle * barea * + sigmaf * log10( dble( nsamples ) ) / + dble( nsamples ) ) ** ( 1.0d0 / n ) C VERIFY THE LOCAL SEARCH START CONDITION indreg = .true. isample = 1 do while( indreg ) C SET THE PROBABILISTIC THRESHOLD do while( ( isample .lt. nsamples ) .and. + ( f .le. vetf(isample) ) ) isample = isample + 1 end do if( isample .lt. nsamples ) then do i = 1,n vAux(i) = x(i) - + vetx( (isample - 1) * n + i) end do C SET THE NORM BETWEEN x AND X_j psearch = fNormaVetor(n,vAux) C SET varphi_k if( indPiecewiseLinearThreshould ) then if( ( psearch ) .le. + ( 0.5d0 * rk ) ) then varphi_k = 0.0d0 else if( ( psearch ) .ge. + ( 1.5d0 * rk ) ) then varphi_k = 1.0d0 else varphi_k = (psearch - 0.5d0 + * rk) / rk end if end if seed = ( nsamples + isample ) / 2.0d0 + 1 vdummy = drand(seed) vdummy = drand(seed) else if( psearch .le. rk ) then varphi_k = 0.0d0 else varphi_k = 1.0d0 end if vdummy = 0.5d0 end if if( varphi_k .le. vdummy ) then indreg = .false. itsample = itsample + 1 else isample = isample + 1 end if else indreg = .false. indsplpt = .false. consbadsample = 0 end if end do if( ( indsplpt ) .and. + ( trysample .ge. maxtrysample ) ) then indsplpt = .false. badsample = badsample + 1 consbadsample = consbadsample + 1 C CHECK IF THE CRITIC DISTANCE NEEDS AN ADJUST IN SIGMAF if( indPiecewiseLinearThreshould ) then if( mod(badsample,maxbadchoice) .eq. 0 ) then sigmaf = sigmaf / 2.0d0 end if end if end if end do else C MULTI-START METHOD seed = nsamples + 1 vdummy = drand(seed) do i = 1,n x(i) = l(i) + ( u(i) - l(i) ) * drand(seed) end do nsamples = nsamples + 1 end if 550 continue C CALL OPTIMIZATION SOLVER call solver(n,x,l,u,m,lambda,rho,equatn,linear,gtype,hptype, + precond,checkder,epsfeas,epsopt,maxoutit,maxtotit,maxtotfc, + iprint,ncomp,f,snorm,nalpsupn,outiter,totiter,totfcnt,totgcnt, + totcgcnt,time,inform,wi1,wd1,wd2,wd3,wd4,wd5,wd6,wd7,wd8,wd9, + wd10,wd11,wd12,wd13,wd14,wd15,wd16,wd17,wd18,wd19) totevalf = totevalf + totfcnt totevalg = totevalg + totgcnt C TAKES ITERATION TIME tottime = etime(tarray) C VERIFY IF ITS THE BEST SOLUTION if( fmin .gt. f + 0.0001d0 * abs(f) + 1.0E-6 ) then btime = tottime btrial = ktrial btunnit = tunnit bltrial = totltrial bevalf = totevalf bevalg = totevalg bnsamples = nsamples bbadltrial = badltrial bbadtunnit = badtunnit fmin = f bfound = bfound + 1 do i = 1,n xmin(i) = x(i) end do write(*,0023) bfound,btrial,fmin,btime,bnsamples,btunnit end if C CHECK TIME SPENT if( tottime - ctime .gt. 60.0d0 ) then ctime = tottime write(*,0010) ctime,ktrial end if C ITERATE ktrial = ktrial + 1 C CHECKS THE TUNNELING METHOD if( indLissajous ) then call lissajous(n,x,l,u,ltrial,err) totevalf = totevalf + ltrial if( err .eq. 1 ) then tunnit = tunnit + 1 totltrial = totltrial + ltrial go to 550 else badtunnit = badtunnit + 1 badltrial = badltrial + ltrial end if end if end do 650 continue C FINISHES TIME COUNT tottime = etime(tarray) C WRITES THE BEST SOLUTION FOUND write(*,*) " " write(*,*) " " write(*,*) " " write(*,0030) write(*,0031) fmin,fmin write(*,0035) write(*,0032) btime,tottime write(*,0033) btrial,bnsamples write(*,0034) bevalf,bevalg write(*,0037) write(*,*) " " write(*,*) '* Original experiments executed with G77. This ', + 'result may reach a diferent soluction since G77 and ', + 'gfortran are not fully compatible. *' C WRITE ADDITIONAL OUTPUT INFORMATION CODED BY THE USER call endp(n,x,l,u,m,lambda,rho,equatn,linear) stop C NON-EXECUTABLE STATEMENTS 0010 format("----- TIME : ",F7.2," -- Iter : ",I4, + " ----------------------") 0019 format("-------------< ",A10," (sigma = ",F6.3,")>----", + "------------") 0020 format("-----------------< Best solution found ", + ">-------------------") 0021 format(" # Iter f* t*", + " #x0 iter tun") 0022 format("--- ----- ------------------ --------", + " ------- --------") 0023 format(I3,2X,I5,2X,F18.10,2X,F7.2,2X,I7,2x,I7) 0024 format(" # Iter f* t*", + " #x0") 0025 format(" # Iter f* t*", + " Iter TUNN proc") 0030 format("******************** FINAL SOLUTION ", + "********************") 0031 format("f(x*) : ",D25.16," ( ",F18.10," )" ) 0032 format("time(x*) : ",F12.5," | total time : ",F12.5 ) 0033 format("iter(x*) : ",I4,8X," | #x0 : ",I7 ) 0034 format("eval_f : ",I7,5X," | eval_g : ",I7 ) 0035 format("------------------------------------", + "--------------------") 0036 format("x( ",I4," ) = ",D25.16) 0037 format("************************************", + "********************") end C ****************************************************************** C UNLESS YOU ARE AN EXPERIENT USER, YOUR MODIFICATIONS FINISH HERE. C ****************************************************************** C ****************************************************************** C ALL THE INFORMATION THAT ALGENCAN PRINTS IN THE SCREEN WILL ALSO C BE SAVED IN TEXT FILE CALLED algencan.out. THE SOLUTION, AS WELL C AS THE FINAL LAGRANGE MULTIPLIERS APPROXIMATION AND THE FINAL C PENALTY PARAMETERS WILL BE SAVED IN A TEXT FILE CALLED C solution.txt. C ****************************************************************** C ****************************************************************** C ****************************************************************** subroutine param(gtype,hptype,precond,checkder,epsfeas,epsopt, +maxoutit,maxtotit,maxtotfc,iprint,ncomp) implicit none C Parameters of the subroutine: C ============================= C C On Entry: C ========= C C This subroutine has no input parameters. C C On Return: C ========== C C GTYPE integer C ------------- C C Type of first derivatives calculation according to the following C convention: C C 0 means true first derivatives. In this case, subroutines evalg C and evaljac must be modified by the user to compute the C gradient of the objective function and the gradients of the C constraints, respectively. C C 1 means that a finite difference approximation will be used. In C this case, subroutines evalg and evaljac may have an empty body C but must be present. It is also recommended that those C empty-body subroutines set flag = - 1. Last but not least, the C option gtype = 1 is not cheap neither safe. C C HPTYPE integer C -------------- C C Type of Hessian-vector product according to the following C convention: C C We will first describe the meaning if the choices in the Augmented C Lagrangian framework. The way in which the product of the Hessian C matrix of the Augmented Lagrangian by a vector will be done C depends on the value of the parameter hptype in the following way: C C 9 means that an incremental quotients approximation without any C extra consideration will be used. This option requires the C evaluation of an extra gradient at each Conjugate Gradient C iteration. If gtype = 0 then this gradient evaluation will be C done using the user supplied subroutines evalg and evaljac C (evalc will also be used). On the other hand, if gtype = 1, the C gradient calculation will be done using just calls to the user C provided subroutines evalf and evalc. nind calls will be done, C where nind is the dimension of the current face of the C active-set method. This option is not cheap neither safe. C C If you did not code subroutines evalg and evaljac, to compute C the gradient of the objective function and the Jacobian of the C constraints then your options finished here. C C 0 means that subroutines to compute the Hessian of the objective C function (evalh) and the Hessians of the constraints (evalhc) C were provided by the user. So, the product of the Hessian of the C Augmented Lagrangian times a vector will be computed using the C Hessians provided by these subroutines and then adding the first C order term (for the first order term the user-supplied C subroutine to compute the Jacobian of the constraints (evaljac) C is also used). C C 1 means that, instead of providing individual subroutines to C compute the Hessians of the objective function and the C constraints, the user provided a subroutine to compute the C product of an arbitrary vector times the Hessian of the C Lagrangian. C C 2 means that incremental quotients will be used. The difference C between hptype = 9 and hptype = 2 is that, in the latter case, C the non-differentiability of the Hessian of the Augmented C Lagrangian will be taken into account. In particular, when C computing the gradient of the Augmented Lagrangian at C (x + step p), the constraints that will be considered will be C the same constraints that contributed to the computation of the C gradient of the Augmented Lagrangian at x. C C If GENCAN is been used to solve a bound-constrained problem C which is not the subproblem of an Augmented Lagrangian method, C but an isolated bound-constrained problem, then there is no C difference between this option and hptype = 9. C C This option also requires the evaluation of an extra gradient at C each Conjugate Gradient iteration. C C 3 is similar to hptype = 2. The difference is that the C contribution of the linear constraints in the Hessian matrix C will be computed explicitly. If the problem has not linear C constraints then this option is identical to hptype = 2. C Moreover, if the problem has no constraints then this option is C equal to hptype = 9. C C This option also requires the evaluation of an extra gradient at C each Conjugate Gradient iteration. C C 4 means that the Hessian matrix will be approximated and then the C product of the Hessian approximation by the vector will be C computed exactly. In particular, the Hessian matrix will be C approximated doing a BFGS correction to the Gauss-Newton C approximation of the Hessian. Before the BFGS correction, a C structured spectral correction is done to force the Gauss-Newton C approximation to be positive definite. C C If the problem has not constraints then the approximation C reduces to a BFGS approximation of the Hessian (without memory) C and using the spectral approximation (instead of the identity) C as initial approximation. C C Numerical experiments suggested that this option is convenient C just for constrained problems. This motivated the introduction C of the next option. C C This option does NOT require an extra gradient evaluation per C iteration and, in this sense, each CG iteration is C computationally cheaper than a CG iteration of the previous C choices. However, the approximation of the Hessian matrix C requires some information (mainly the Jacobian of the C constraints) that must be saved during the gradient evaluation. C To save this information requires an amount of memory C proportional to the number of non-null elements of the Jacobian C matrix. C C Quadratic subproblems are convex with this choice. C C 5 is an adaptive strategy that choose, at every iteration, between C 2 and 4. When the gradient of the Augmented Lagrangian is C computed, it is verified if at least a constraint contributes to C the calculation. If this is the case, 4 is used. Otherwise, 2 is C used. C C For problems with equality constraints (that always contributes C to the Augmented Lagrangian function) this option is identical C to 4. C C For problems without constraints this option is identical to 2. C C 6 is identical to 5 but the choice is made between 3 and 4 instead C of between 2 and 4. C C For problems with equality constraints (that always contributes C to the Augmented Lagrangian function) this option is identical C to 4. C C For problems without constraints this option is identical to 3. C C We will now describe the meaning if the choices for unconstrained C and bound-constrained problems. In this context the way in which C the product of the Hessian matrix by a vector will be done depends C on the value of the parameter hptype in the following way: C C 0 means that the subroutine to compute the Hessian of the C objective function (evalh) was provided by the user. So, the C product of the Hessian times a vector will be computed using the C Hessian provided by this subroutine. C C 1 means that a subroutine (evalhlp) to compute the product of the C Hessian of the objective function times an arbitrary vector is C being provided by the user. C C 9 means that an incremental quotients approximation will be used. C This option requires the evaluation of an extra gradient at each C Conjugate Gradient iteration. If gtype = 0 then this gradient C evaluation will be done using the user supplied subroutine C evalg. On the other hand, if gtype = 1, the gradient calculation C will be done using just calls to the user provided subroutine C evalf. nind calls will be done, where nind is the dimension of C the current face of the active-set method. C C If you did not code subroutine evalg to compute the gradient of C the objective function then your options finished here. C C 4 means that the Hessian matrix will be approximated and then the C product of the Hessian approximation by the vector will be C computed exactly. The approximation is a BFGS approximation of C the Hessian (without memory) and using the spectral C approximation (instead of the identity) as initial C approximation. C C Numerical experiments suggested that this option is not C convenient for unconstrained or just bound-constrained problems. C (Note that this option was developed to be used in the Augmented C Lagrangian framework.) C C This option does NOT require an extra gradient evaluation per C iteration and, in this sense, each CG iteration is C computationally cheaper than a CG iteration of the previous C choices. C C Quadratic subproblems are convex with this choice. C C In the bound-constrained context, options hptype = 2,3,5 and 6 C are all identical to hptype = 9. C C PRECOND character * 6 C --------------------- C C Indicates the type of preconditioning that will be used for C Conjugates Gradients according to the following convention: C C 'NONE' means no preconditioner at all. C C 'QNAGNC' means Quasi-Newton Correction of the Gauss-Newton C approximation of the Hessian. The exact form is this C preconditioner is described in: C C E. G. Birgin and J. M. Mart�ez, "Structured minimal- C memory inexact quasi-Newton method and secant C preconditioners for Augmented Lagrangian Optimization", C submitted, 2005. C C CHECKDER logical C ---------------- C C If you are using finite differences aproximations for the C derivatives (gtype = 1) then checkder must be set to FALSE. C On the other hand, if you are using your own coded derivatives C you may would like to test them against a finite differences C approximation to verify if they are o.k. In this case, set C checkder = TRUE. C C EPSFEAS double precision C ------------------------ C C Feasibility tolerance for the sup-norm of the constraints. C (Ignored in the unconstrained and bound-constrained cases.) C C EPSOPT double precision C ----------------------- C C Optimality tolerance for the sup-norm of the projected gradient of C the Augmented Lagrangian in the constrained case and the sup-norm C of the projected gradient of the objective function in the C unconstrained and the bound-constrained cases. C C MAXOUTIT integer C ---------------- C C Maximum number of Augmented Lagrangian (outer) iterations. C (Ignored in the unconstrained and bound-constrained cases.) C C MAXTOTIT integer C ---------------- C C Maximum total number of inner iterations in the Augmented C Lagrangian context (total means summing up the inner iterations of C each outer iteration). In the unconstrained and bound-constrained C cases it means just the maximum number of iterations. C C MAXTOTFC integer C ---------------- C C Idem MAXTOTIT but for number of functional evaluations. C C IPRINT integer C -------------- C C Controls the ammount of information of the output according to the C following convention: C C 0 means no output at all. C C 1 means information at each outer iteration but without any C information of how the subproblems are being solved. C C 2 means the same as 1 plus information of each inner iteration. C C 3 means the same as 2 plus information of the line searches and C the calculation of the truncated Newton direction (using CG) of C each inner iteration. C C In all cases, an output file named solution.txt with the final C point, Lagrange mutipliers and penalty parameters will be C generated. Moreover, the same output of the screen will be saved C in a file named algencan.out. C C NCOMP integer C ------------- C C Every time a vector is printed, just its first ncomp component C will be displayed. C SCALAR ARGUMENTS logical checkder character * 6 precond integer gtype,hptype,iprint,maxoutit,maxtotfc,maxtotit,ncomp double precision epsfeas,epsopt C ****************************************************************** C FROM HERE ON YOU MUST MODIFY THE SUBROUTINE TO SET SOME ALGENCAN C ARGUMENTS RELATED TO DERIVATIVES, STOPPING CRITERIA AND OUTPUT: C ****************************************************************** gtype = 0 hptype = 6 C precond = 'QNCGNA' precond = 'NONE' checkder = .false. epsfeas = 1.0d-04 epsopt = 1.0d-04 maxoutit = 50 maxtotit = 1000000 maxtotfc = 5 * maxtotit iprint = 0 ncomp = 5 C ****************************************************************** C STOP HERE YOUR MODIFICATIONS OF SUBROUTINE PARAM. C ****************************************************************** end C ****************************************************************** C STOP HERE ALL YOUR MODIFICATIONS !!!!!! C ******************************************************************