************************************************************************ * SUBROUTINE PNECU ALL SYSTEMS 97/01/22 * PURPOSE : * EASY TO USE SUBROUTINE FOR UNCONSTRAINED MINIMIZATION OF * FUNCTIONS WITH LARGE-SCALE SPARSE HESSIAN MATRICES * * PARAMETERS : * II NF NUMBER OF VARIABLES. * IU MH NUMBER OF NONZERO ELEMENTS IN THE UPPER PART OF THE HESSIAN * MATRIX. * RI X(NF) VECTOR OF VARIABLES. * II IH(NF+1) POINTERS OF DIAGONAL ELEMENTS OF THE MATRIX H. * IU JH(MH) INDICES OF NONZERO ELEMENTS OF THE MATRIX H * TOGETHER WITH AN ADDITIONAL SPACE USED FOR THE NUMERICAL * DIFFERENTIATION. * II IPAR(7) INTEGER PAREMETERS: * IPAR(1) MAXIMUM NUMBER OF ITERATIONS. * IPAR(2) MAXIMUM NUMBER OF FUNCTION EVALUATIONS. * IPAR(3) MAXIMUM NUMBER OF GRADIENT EVALUATIONS. * IPAR(4) ESTIMATION INDICATOR. IPAR(4)=0-MINIMUM IS NOT * ESTIMATED. IPAR(4)=1-MINIMUM IS ESTIMATED BY THE VALUE * RPAR(6). * IPAR(5) METHOD FOR COMPUTING A TRUST REGION STEP. IPAR(5)=1-THE * STEIHAUG-TOINT METHOD. IPAR(5)=2-THE SHIFTED STEIHAUG-TOINT * METHOD WITH FIVE LANCZOS STEPS. IPAR(5)>2-THE SHIFTED * STEIHAUG-TOINT METHOD WITH IPAR(5) LANCZOS STEPS. * IPAR(6) TYPE OF PRECONDITIONING. IPAR(6)=1-PRECONDITIONING IS * NOT USED. IPAR(6)=2-PRECONDITIONING BY THE INCOMPLETE * GILL-MURRAY DECOMPOSITION. IPAR(6)=3-PRECONDITIONING BY THE * INCOMPLETE GILL-MURRAY DECOMPOSITION WITH A PRELIMINARY * SOLUTION OF THE PRECONDITIONED SYSTEM WHICH IS USED IF IT * SATISFIES THE TERMINATION CRITERION. * IPAR(7) NUMBER DEFINING THE SPACE FOR FILL-IN (THE SIZE OF THIS * SPACE IS IFIL TIMES THE STANDARD SIZE). THE DEFAULT VALUE IS * IFIL=1. THE DEFAULT VALUE HAS TO BE INCREASED IF ITERM IS * LESS OR EQUAL TO -40. * RI RPAR(9) REAL PARAMETERS: * RPAR(1) MAXIMUM STEPSIZE. * RPAR(2) TOLERANCE FOR THE CHANGE OF VARIABLES. * RPAR(3) TOLERANCE FOR THE CHANGE OF FUNCTION VALUES. * RPAR(4) TOLERANCE FOR THE FUNCTION FALUE. * RPAR(5) TOLERANCE FOR THE GRADIENT NORM. * RPAR(6) ESTIMATION OF THE MINIMUM FUNCTION VALUE. * RPAR(7) INITIAL TRUST-REGION RADIUS. * RPAR(8) THIS PARAMETER IS NOT USED IN THE SUBROUTINE PNEC. * RPAR(9) THIS PARAMETER IS NOT USED IN THE SUBROUTINE PNEC. * RO F VALUE OF THE OBJECTIVE FUNCTION. * RO GMAX MAXIMUM PARTIAL DERIVATIVE. * II ISPAS INPUT SPARSE STRUCTURE. ISPAS=1-STANDARD COORDINATE * FORM. ISPAS=2-SPARSE STRUCTURE COMPRESSED BY ROWS. * II IPRNT PRINT SPECIFICATION. IPRNT=0-NO PRINT. * ABS(IPRNT)=1-PRINT OF FINAL RESULTS. * ABS(IPRNT)=2-PRINT OF FINAL RESULTS AND ITERATIONS. * IPRNT>0-BASIC FINAL RESULTS. IPRNT<0-EXTENDED FINAL * RESULTS. * IO ITERM VARIABLE THAT INDICATES THE CAUSE OF TERMINATION. * ITERM=1-IF ABS(X-XO) WAS LESS THAN OR EQUAL TO TOLX IN * MTESX (USUALLY TWO) SUBSEQUENT ITERATIONS. * ITERM=2-IF ABS(F-FO) WAS LESS THAN OR EQUAL TO TOLF IN * MTESF (USUALLY TWO) SUBSEQUENT ITERATIONS. * ITERM=3-IF F IS LESS THAN OR EQUAL TO TOLB. * ITERM=4-IF GMAX IS LESS THAN OR EQUAL TO TOLG. * ITERM=6-IF THE TERMINATION CRITERION WAS NOT SATISFIED, * BUT THE SOLUTION OBTAINED IS PROBABLY ACCEPTABLE. * ITERM=11-IF NIT EXCEEDED MIT. ITERM=12-IF NFV EXCEEDED MFV. * ITERM=13-IF NFG EXCEEDED MFG. ITERM<0-IF THE METHOD FAILED. * VALUES ITERM<=-40 DETECT A LACK OF SPACE. IN THIS CASE, * PARAMETER IPAR(7) HAS TO BE INCREASED. * * VARIABLES IN COMMON /STAT/ (STATISTICS) : * IO NRES NUMBER OF RESTARTS. * IO NDEC NUMBER OF MATRIX DECOMPOSITIONS. * IO NIN NUMBER OF INNER ITERATIONS. * IO NIT NUMBER OF ITERATIONS. * IO NFV NUMBER OF FUNCTION EVALUATIONS. * IO NFG NUMBER OF GRADIENT EVALUATIONS. * IO NFH NUMBER OF HESSIAN EVALUATIONS. * * SUBPROGRAMS USED : * S PNEC DISCRETE NEWTON METHOD WITH ITERATIVE SOLUTION OF THE * TRUST-REGION SUBPROBLEM. * S PFSED3 COMPRESSED SPARSE STRUCTURE OF THE HESSIAN MATRIX IS * COMPUTED FROM THE COORDINATE FORM. * S MXVICP COPYING OF AN INTEGER VECTOR. * * EXTERNAL SUBROUTINES : * SE OBJ COMPUTATION OF THE VALUE OF THE OBJECTIVE FUNCTION. * CALLING SEQUENCE: CALL OBJ(NF,X,FF) WHERE NF IS THE NUMBER * OF VARIABLES, X(NF) IS THE VECTOR OF VARIABLES AND FF IS * THE VALUE OF THE OBJECTIVE FUNCTION. * SE DOBJ COMPUTATION OF THE GRADIENT OF THE OBJECTIVE FUNCTION. * CALLING SEQUENCE: CALL DOBJ(NF,X,GF) WHERE NF IS THE NUMBER * OF VARIABLES, X(NF) IS THE VECTOR OF VARIABLES AND GF(NF) * IS THE GRADIENT OF THE OBJECTIVE FUNCTION. * SUBROUTINE PNECU(NF,MH,X,IH,JH,IPAR,RPAR,F,GMAX,ISPAS,IPRNT, & ITERM) INTEGER NF,MH,IH(*),JH(*),IPAR(7),ISPAS,IPRNT,ITERM DOUBLE PRECISION X(*),RPAR(9),F,GMAX INTEGER NB,LGF,LHF,LS,LXO,LGO,LXS,LGS,LCOL,LWN11,LWN12,LIW,LJH, & IFIL,IER INTEGER NRES,NDEC,NIN,NIT,NFV,NFG,NFH COMMON /STAT/ NRES,NDEC,NIN,NIT,NFV,NFG,NFH INTEGER IA(:) DOUBLE PRECISION RA(:) ALLOCATABLE IA,RA IFIL=IPAR(7) IF (IFIL.LE.0) IFIL=1 IF (ISPAS.LE.1) THEN CALL PFSED3(NF,MH,IH,JH,IER) IF (IER.NE.0) THEN WRITE (6,'(''INPUT ERROR : IER = '',I3)') IER STOP END IF ELSE MH=IH(NF+1)-1 END IF ALLOCATE (IA(4*NF+3+(IFIL+2)*MH),RA(6*NF+1+(IFIL+2)*MH)) NB=0 * * POINTERS FOR AUXILIARY ARRAYS * LGF=1 LS=LGF+NF LXO=LS+NF LGO=LXO+NF LXS=LGO+NF+1 LGS=LXS+NF LHF=LGS+NF LCOL=1 LWN11=LCOL+NF LWN12=LWN11+NF+1 LIW=LWN12+NF+1 LJH=LIW+NF+1 CALL MXVICP(IH(NF+1)-1,JH,IA(LJH)) CALL PNEC(NF,NB,(IFIL+2)*MH,X,IA,RA,RA,RA(LGF),RA(LHF),IH, & IA(LJH),RA(LS),RA(LXO),RA(LGO),RA(LXS),RA(LGS),IA(LCOL), & IA(LWN11),IA(LWN12),IA(LIW),RPAR(1),RPAR(2),RPAR(3),RPAR(4), & RPAR(5),RPAR(6),RPAR(7),GMAX,F,IPAR(1),IPAR(2),IPAR(3),IPAR(4), & IPAR(5),IPAR(6),IPRNT,ITERM) DEALLOCATE (IA,RA) RETURN END ************************************************************************ * SUBROUTINE PNECS ALL SYSTEMS 97/01/22 * PURPOSE : * EASY TO USE SUBROUTINE FOR BOX CONSTRAINED MINIMIZATION OF * LARGE-SCALE FUNCTIONS WITH SPARSE HESSIAN MATRICES * * PARAMETERS : * II NF NUMBER OF VARIABLES. * IU MH NUMBER OF NONZERO ELEMENTS IN THE UPPER PART OF THE HESSIAN * MATRIX. * RI X(NF) VECTOR OF VARIABLES. * II IX(NF) VECTOR CONTAINING TYPES OF BOUNDS. IX(I)=0-VARIABLE * X(I) IS UNBOUNDED. IX(I)=1-LOWER BOUND XL(I).LE.X(I). * IX(I)=2-UPPER BOUND X(I).LE.XU(I). IX(I)=3-TWO SIDE BOUND * XL(I).LE.X(I).LE.XU(I). IX(I)=5-VARIABLE X(I) IS FIXED. * RI XL(NF) VECTOR CONTAINING LOWER BOUNDS FOR VARIABLES. * RI XU(NF) VECTOR CONTAINING UPPER BOUNDS FOR VARIABLES. * II IH(NF+1) POINTERS OF DIAGONAL ELEMENTS OF THE MATRIX H. * IU JH(MH) INDICES OF NONZERO ELEMENTS OF THE MATRIX H * TOGETHER WITH AN ADDITIONAL SPACE USED FOR THE NUMERICAL * DIFFERENTIATION. * II IPAR(7) INTEGER PAREMETERS: * IPAR(1) MAXIMUM NUMBER OF ITERATIONS. * IPAR(2) MAXIMUM NUMBER OF FUNCTION EVALUATIONS. * IPAR(3) MAXIMUM NUMBER OF GRADIENT EVALUATIONS. * IPAR(4) ESTIMATION INDICATOR. IPAR(4)=0-MINIMUM IS NOT * ESTIMATED. IPAR(4)=1-MINIMUM IS ESTIMATED BY THE VALUE * RPAR(6). * IPAR(5) METHOD FOR COMPUTING A TRUST REGION STEP. IPAR(5)=1-THE * STEIHAUG-TOINT METHOD. IPAR(5)=2-THE SHIFTED STEIHAUG-TOINT * METHOD WITH FIVE LANCZOS STEPS. IPAR(5)>2-THE SHIFTED * STEIHAUG-TOINT METHOD WITH IPAR(5) LANCZOS STEPS. * IPAR(6) TYPE OF PRECONDITIONING. IPAR(6)=1-PRECONDITIONING IS * NOT USED. IPAR(6)=2-PRECONDITIONING BY THE INCOMPLETE * GILL-MURRAY DECOMPOSITION. IPAR(6)=3-PRECONDITIONING BY THE * INCOMPLETE GILL-MURRAY DECOMPOSITION WITH A PRELIMINARY * SOLUTION OF THE PRECONDITIONED SYSTEM WHICH IS USED IF IT * SATISFIES THE TERMINATION CRITERION. * IPAR(7) NUMBER DEFINING THE SPACE FOR FILL-IN (THE SIZE OF THIS * SPACE IS IFIL TIMES THE STANDARD SIZE). THE DEFAULT VALUE IS * IFIL=1. THE DEFAULT VALUE HAS TO BE INCREASED IF ITERM IS * LESS OR EQUAL TO -40. * RI RPAR(9) REAL PARAMETERS: * RPAR(1) MAXIMUM STEPSIZE. * RPAR(2) TOLERANCE FOR THE CHANGE OF VARIABLES. * RPAR(3) TOLERANCE FOR THE CHANGE OF FUNCTION VALUES. * RPAR(4) TOLERANCE FOR THE FUNCTION FALUE. * RPAR(5) TOLERANCE FOR THE GRADIENT NORM. * RPAR(6) ESTIMATION OF THE MINIMUM FUNCTION VALUE. * RPAR(7) INITIAL TRUST-REGION RADIUS. * RPAR(8) THIS PARAMETER IS NOT USED IN THE SUBROUTINE PNEC. * RPAR(9) THIS PARAMETER IS NOT USED IN THE SUBROUTINE PNEC. * RO F VALUE OF THE OBJECTIVE FUNCTION. * RO GMAX MAXIMUM PARTIAL DERIVATIVE. * II ISPAS INPUT SPARSE STRUCTURE. ISPAS=1-STANDARD COORDINATE * FORM. ISPAS=2-SPARSE STRUCTURE COMPRESSED BY ROWS. * II IPRNT PRINT SPECIFICATION. IPRNT=0-NO PRINT. * ABS(IPRNT)=1-PRINT OF FINAL RESULTS. * ABS(IPRNT)=2-PRINT OF FINAL RESULTS AND ITERATIONS. * IPRNT>0-BASIC FINAL RESULTS. IPRNT<0-EXTENDED FINAL * RESULTS. * IO ITERM VARIABLE THAT INDICATES THE CAUSE OF TERMINATION. * ITERM=1-IF ABS(X-XO) WAS LESS THAN OR EQUAL TO TOLX IN * MTESX (USUALLY TWO) SUBSEQUENT ITERATIONS. * ITERM=2-IF ABS(F-FO) WAS LESS THAN OR EQUAL TO TOLF IN * MTESF (USUALLY TWO) SUBSEQUENT ITERATIONS. * ITERM=3-IF F IS LESS THAN OR EQUAL TO TOLB. * ITERM=4-IF GMAX IS LESS THAN OR EQUAL TO TOLG. * ITERM=6-IF THE TERMINATION CRITERION WAS NOT SATISFIED, * BUT THE SOLUTION OBTAINED IS PROBABLY ACCEPTABLE. * ITERM=11-IF NIT EXCEEDED MIT. ITERM=12-IF NFV EXCEEDED MFV. * ITERM=13-IF NFG EXCEEDED MFG. ITERM<0-IF THE METHOD FAILED. * VALUES ITERM<=-40 DETECT A LACK OF SPACE. IN THIS CASE, * PARAMETER IPAR(7) HAS TO BE INCREASED. * * VARIABLES IN COMMON /STAT/ (STATISTICS) : * IO NRES NUMBER OF RESTARTS. * IO NDEC NUMBER OF MATRIX DECOMPOSITIONS. * IO NIN NUMBER OF INNER ITERATIONS. * IO NIT NUMBER OF ITERATIONS. * IO NFV NUMBER OF FUNCTION EVALUATIONS. * IO NFG NUMBER OF GRADIENT EVALUATIONS. * IO NFH NUMBER OF HESSIAN EVALUATIONS. * * SUBPROGRAMS USED : * S PNEC DISCRETE NEWTON METHOD WITH ITERATIVE SOLUTION OF THE * TRUST-REGION SUBPROBLEM. * S PFSED3 COMPRESSED SPARSE STRUCTURE OF THE HESSIAN MATRIX IS * COMPUTED FROM THE COORDINATE FORM. * S MXVICP COPYING OF AN INTEGER VECTOR. * * EXTERNAL SUBROUTINES : * SE OBJ COMPUTATION OF THE VALUE OF THE OBJECTIVE FUNCTION. * CALLING SEQUENCE: CALL OBJ(NF,X,FF) WHERE NF IS THE NUMBER * OF VARIABLES, X(NF) IS THE VECTOR OF VARIABLES AND FF IS * THE VALUE OF THE OBJECTIVE FUNCTION. * SE DOBJ COMPUTATION OF THE GRADIENT OF THE OBJECTIVE FUNCTION. * CALLING SEQUENCE: CALL DOBJ(NF,X,GF) WHERE NF IS THE NUMBER * OF VARIABLES, X(NF) IS THE VECTOR OF VARIABLES AND GF(NF) * IS THE GRADIENT OF THE OBJECTIVE FUNCTION. * SUBROUTINE PNECS(NF,MH,X,IX,XL,XU,IH,JH,IPAR,RPAR,F,GMAX,ISPAS, & IPRNT,ITERM) INTEGER NF,MH,IX(*),IH(*),JH(*),IPAR(7),ISPAS,IPRNT,ITERM DOUBLE PRECISION X(*),XL(*),XU(*),RPAR(8),F,GMAX INTEGER NB,LGF,LHF,LS,LXO,LGO,LXS,LGS,LCOL,LWN11,LWN12,LIW,LJH, & IFIL,IER INTEGER NRES,NDEC,NIN,NIT,NFV,NFG,NFH COMMON /STAT/ NRES,NDEC,NIN,NIT,NFV,NFG,NFH INTEGER IA(:) DOUBLE PRECISION RA(:) ALLOCATABLE IA,RA IFIL=IPAR(7) IF (IFIL.LE.0) IFIL=1 IF (ISPAS.LE.1) THEN CALL PFSED3(NF,MH,IH,JH,IER) IF (IER.NE.0) THEN WRITE (6,'(''INPUT ERROR : IER = '',I3)') IER STOP END IF ELSE MH=IH(NF+1)-1 END IF ALLOCATE (IA(4*NF+3+(IFIL+2)*MH),RA(6*NF+1+(IFIL+2)*MH)) NB=1 * * POINTERS FOR AUXILIARY ARRAYS * LGF=1 LS=LGF+NF LXO=LS+NF LGO=LXO+NF LXS=LGO+NF+1 LGS=LXS+NF LHF=LGS+NF LCOL=1 LWN11=LCOL+NF LWN12=LWN11+NF+1 LIW=LWN12+NF+1 LJH=LIW+NF+1 CALL MXVICP(IH(NF+1)-1,JH,IA(LJH)) CALL PNEC(NF,NB,(IFIL+2)*MH,X,IX,XL,XU,RA(LGF),RA(LHF),IH, & IA(LJH),RA(LS),RA(LXO),RA(LGO),RA(LXS),RA(LGS),IA(LCOL), & IA(LWN11),IA(LWN12),IA(LIW),RPAR(1),RPAR(2),RPAR(3),RPAR(4), & RPAR(5),RPAR(6),RPAR(7),GMAX,F,IPAR(1),IPAR(2),IPAR(3),IPAR(4), & IPAR(5),IPAR(6),IPRNT,ITERM) DEALLOCATE (IA,RA) RETURN END ************************************************************************ * SUBROUTINE PNEC ALL SYSTEMS 01/09/22 * PURPOSE : * GENERAL SUBROUTINE FOR LARGE-SCALE BOX CONSTRAINED MINIMIZATION OF * FUNCTIONS WITH SPARSE HESSIAN MATRICES BASED ON THE DISCRETE NEWTON * METHOD WITH ITERATIVE SOLUTION OF THE TRUST-REGION SUBPROBLEM. * * PARAMETERS : * II NF NUMBER OF VARIABLES. * II NB CHOICE OF SIMPLE BOUNDS. NB=0-SIMPLE BOUNDS SUPPRESSED. * NB>0-SIMPLE BOUNDS ACCEPTED. * II MMAX MAXIMUM DIMENSION OF THE SPARSE TABLEAU. * RU X(NF) VECTOR OF VARIABLES. * II IX(NF) VECTOR CONTAINING TYPES OF BOUNDS. IX(I)=0-VARIABLE * X(I) IS UNBOUNDED. IX(I)=1-LOWER BOUND XL(I).LE.X(I). * IX(I)=2-UPPER BOUND X(I).LE.XU(I). IX(I)=3-TWO SIDE BOUND * XL(I).LE.X(I).LE.XU(I). IX(I)=5-VARIABLE X(I) IS FIXED. * RI XL(NF) VECTOR CONTAINING LOWER BOUNDS FOR VARIABLES. * RI XU(NF) VECTOR CONTAINING UPPER BOUNDS FOR VARIABLES. * RU GF(NF) GRADIENT OF THE OBJECTIVE FUNCTION. * RA HF(MMAX) NONZERO ELEMENTS OF THE SPARSE HESSIAN MATRIX * TOGETHER WITH AN ADDITIONAL SPACE USED FOR THE NUMERICAL * DIFFERENTIATION. * II IH(NF+1) POINTERS OF DIAGONAL ELEMENTS OF THE MATRIX H. * IU JH(MMAX) INDICES OF NONZERO ELEMENTS OF THE MATRIX H * TOGETHER WITH AN ADDITIONAL SPACE USED FOR THE NUMERICAL * DIFFERENTIATION. * RO S(NF) DIRECTION VECTOR. * RU XO(NF) VECTORS OF VARIABLES DIFFERENCE. * RI GO(NF) GRADIENTS DIFFERENCE. * RA XS(NF) AUXILIARY VECTOR. * RA GS(NF) AUXILIARY VECTOR. * IA COL(NF) AUXILIARY ARRAY. * IA WN11(NF+1) AUXILIARY VECTOR. * IA WN12(NF+1) AUXILIARY VECTOR. * IA IW(NF+1) AUXILIARY VECTOR. * RI XMAX MAXIMUM STEPSIZE. * RI TOLX TOLERANCE FOR CHANGE OF VARIABLES. * RI TOLF TOLERANCE FOR CHANGE OF FUNCTION VALUES. * RI TOLB TOLERANCE FOR THE FUNCTION VALUE. * RI TOLG TOLERANCE FOR THE GRADIENT NORM. * RI FMIN ESTIMATION OF THE MINIMUM FUNCTION VALUE. * RI XDEL TRUST REGION STEPSIZE. * RO GMAX MAXIMUM PARTIAL DERIVATIVE. * RO F VALUE OF THE OBJECTIVE FUNCTION. * II MIT MAXIMUM NUMBER OF ITERATIONS. * II MFV MAXIMUM NUMBER OF FUNCTION EVALUATIONS. * II MFG MAXIMUM NUMBER OF GRADIENT EVALUATIONS. * II IEST ESTIMATION INDICATOR. IEST=0-MINIMUM IS NOT ESTIMATED. * IEST=1-MINIMUM IS ESTIMATED BY THE VALUE RPAR(8) * II MOS1 METHOD FOR COMPUTING A TRUST REGION STEP. MOS1=1-THE * STEIHAUG-TOINT METHOD. MOS1=2-THE SHIFTED STEIHAUG-TOINT * METHOD WITH FIVE LANCZOS STEPS. MOS1>2-THE SHIFTED * STEIHAUG-TOINT METHOD WITH MOS1 LANCZOS STEPS. * II MOS2 TYPE OF PRECONDITIONING. MOS2=1-PRECONDITIONING IS NOT * USED. MOS2=2-PRECONDITIONING BY THE INCOMPLETE GILL-MURRAY * DECOMPOSITION. MOS2=3-PRECONDITIONING BY THE INCOMPLETE * GILL-MURRAY DECOMPOSITION WITH A PRELIMINARY SOLUTION OF * THE PRECONDITIONED SYSTEM WHICH IS USED IF IT SATISFIES * THE TERMINATION CRITERION. * II IPRNT PRINT SPECIFICATION. IPRNT=0-NO PRINT. * ABS(IPRNT)=1-PRINT OF FINAL RESULTS. * ABS(IPRNT)=2-PRINT OF FINAL RESULTS AND ITERATIONS. * IPRNT>0-BASIC FINAL RESULTS. IPRNT<0-EXTENDED FINAL * RESULTS. * IO ITERM VARIABLE THAT INDICATES THE CAUSE OF TERMINATION. * ITERM=1-IF ABS(X-XO) WAS LESS THAN OR EQUAL TO TOLX IN * MTESX (USUALLY TWO) SUBSEQUENT ITERATIONS. * ITERM=2-IF ABS(F-FO) WAS LESS THAN OR EQUAL TO TOLF IN * MTESF (USUALLY TWO) SUBSEQUENT ITERATIONS. * ITERM=3-IF F IS LESS THAN OR EQUAL TO TOLB. * ITERM=4-IF GMAX IS LESS THAN OR EQUAL TO TOLG. * ITERM=6-IF THE TERMINATION CRITERION WAS NOT SATISFIED, * BUT THE SOLUTION OBTAINED IS PROBABLY ACCEPTABLE. * ITERM=11-IF NIT EXCEEDED MIT. ITERM=12-IF NFV EXCEEDED MFV. * ITERM=13-IF NFG EXCEEDED MFG. ITERM<0-IF THE METHOD FAILED. * VALUES ITERM<=-40 DETECT A LACK OF SPACE. * * VARIABLES IN COMMON /STAT/ (STATISTICS) : * IO NRES NUMBER OF RESTARTS. * IO NDEC NUMBER OF MATRIX DECOMPOSITIONS. * IO NIN NUMBER OF INNER ITERATIONS. * IO NIT NUMBER OF ITERATIONS. * IO NFV NUMBER OF FUNCTION EVALUATIONS. * IO NFG NUMBER OF GRADIENT EVALUATIONS. * IO NFH NUMBER OF HESSIAN EVALUATIONS. * * SUBPROGRAMS USED : * S PCBS04 ELIMINATION OF BOX CONSTRAINT VIOLATIONS. * S PDSGM4 DIRECTION DETERMINATION USING THE STEIHAUG-TOINT AND * SHIFTED STEIHAUG-TOINT TRUST-REGION METHOD. * S PF1HS2 NUMERICAL COMPUTATION OF THE HESSIAN MATRIX USING * DIFFERENCES OF GRADIENTS. * S PS0G01 STEPSIZE SELECTION USING TRUST REGION. * S PYADC0 ADDITION OF A BOX CONSTRAINT. * S PYFUT1 TEST ON TERMINATION. * S PYPTSH DETERMINATION OF GROUPS FOR NUMERICAL DIFFERENTIATION. * S PYRMC0 DELETION OF A BOX CONSTRAINT. * S PYTRCD COMPUTATION OF PROJECTED DIFFERENCES FOR THE VARIABLE * METRIC UPDATE. * S PYTRCG COMPUTATION OF THE PROJECTED GRADIENT. * S PYTRCS COMPUTATION OF THE PROJECTED DIRECTION VECTOR. * S PYTSCH CORRECTION OF THE HESSIAN MATRIX. * S MXSSMI SPARSE SYMMETRIC MATRIX IS REPLACED BY THE UNIT MATRIX. * S MXUDIR VECTOR AUGMENTED BY THE SCALED VECTOR. * S MXVCOP COPYING OF A VECTOR. * S MXVINB PROJECTION OF A SPARSE SYMMETRIC MATRIX TO SATISFY BOX * CONSTRAINTS. * S MXVINE RESTORATION OF A SPARSE SYMMETRIC MATRIX OBTAINED BY * MXVINB * * EXTERNAL SUBROUTINES : * SE OBJ COMPUTATION OF THE VALUE OF THE OBJECTIVE FUNCTION. * CALLING SEQUENCE: CALL OBJ(NF,X,FF) WHERE NF IS THE NUMBER * OF VARIABLES, X(NF) IS THE VECTOR OF VARIABLES AND FF IS * THE VALUE OF THE OBJECTIVE FUNCTION. * SE DOBJ COMPUTATION OF THE GRADIENT OF THE OBJECTIVE FUNCTION. * CALLING SEQUENCE: CALL DOBJ(NF,X,GF) WHERE NF IS THE NUMBER * OF VARIABLES, X(NF) IS THE VECTOR OF VARIABLES AND GF(NF) * IS THE GRADIENT OF THE OBJECTIVE FUNCTION. * * METHOD : * DISCRETE NEWTON METHOD WITH TRUST-REGION STRATEGIES BASED ON * CONJUGATE GRADIENT ITERATIONS. * SUBROUTINE PNEC(NF,NB,MMAX,X,IX,XL,XU,GF,HF,IH,JH,S,XO,GO,XS,GS, & COL,WN11,WN12,IW,XMAX,TOLX,TOLF,TOLB,TOLG,FMIN,XDEL,GMAX,F, & MIT,MFV,MFG,IEST,MOS1,MOS2,IPRNT,ITERM) INTEGER NF,NB,MMAX,IX(*),IH(*),JH(*),COL(*),WN11(*),WN12(*), & IW(*),MIT,MFV,MFG,IEST,MOS2,MOS1,IPRNT,ITERM DOUBLE PRECISION X(*),XL(*),XU(*),GF(*),HF(*),S(*),XO(*),GO(*), & XS(*),GS(*),XMAX,TOLX,TOLF,TOLB,TOLG,FMIN,XDEL,GMAX,F INTEGER IDECF,ITERD,ITERS,ITERH,KD,LD,NTESX,NTESF,MTESX,MTESF, & MRED,KIT,IREST,KBF,MES1,MES2,MES3,MOS3,MAXST,IDIR,ISYS,ITES, & KTERS,IRES1,IRES2,NRED,INEW,IOLD,I,M,MH,N DOUBLE PRECISION R,RO,FO,FP,P,PO,PP,GNORM,GNORMO,SNORM,RMAX,FMAX, & DMAX,UMAX,ETA0,ETA2,ETA9,EPS4,EPS5,EPS8,EPS9,BET1,BET2,GAM1, & GAM2,DEL1,XDELO INTEGER NRES,NDEC,NIN,NIT,NFV,NFG,NFH COMMON /STAT/ NRES,NDEC,NIN,NIT,NFV,NFG,NFH IF (ABS(IPRNT).GT.1) WRITE(6,'(1X,''ENTRY TO PNEC :'')') * * INITIATION * KBF=0 IF (NB.GT.0) KBF=2 NRES=0 NDEC=0 NIN=0 NIT=0 NFV=0 NFG=0 NFH=0 IDIR=0 ISYS=0 ITES=1 MTESX=2 MTESF=2 ITERM=0 ITERD=0 ITERS=2 ITERH=0 KTERS=0 IREST=0 IRES1=999 IRES2=0 IDECF=1 MRED=10 MES1=3 MES2=2 MES3=1 MOS3=1 ETA0=1.0D-15 ETA2=1.0D-8 ETA9=1.0D 120 EPS4=0.10D 0 EPS5=0.90D 0 EPS8=1.00D 0 EPS9=1.00D-8 BET1=0.05D 0 BET2=0.75D 0 GAM1=2.0D 0 GAM2=1.0D 6 DEL1=0.95D 0 RMAX=ETA9 DMAX=ETA9 FMAX=1.0D 60 IF (IEST.LE.0) FMIN=-1.0D 60 IF (IEST.GT.0) IEST= 1 IF (XMAX.LE.0.0D 0) XMAX=1.0D 16 XDEL=MIN(XDEL,XMAX) IF (TOLX.LE.0.0D 0) TOLX=1.0D-16 IF (TOLF.LE.0.0D 0) TOLF=1.0D-14 IF (TOLG.LE.0.0D 0) TOLG=1.0D-6 IF (TOLB.LE.0.0D 0) TOLB=FMIN+1.0D-16 IF (MIT.LE.0) MIT=5000 IF (MFV.LE.0) MFV=5000 IF (MFG.LE.0) MFG=10000 IF (MOS1.LE.0) MOS1=2 IF (MOS1.EQ.2) MOS1=5 IF (MOS2.LE.0) MOS2=2 KD= 2 LD=-1 KIT=0 FO=FMIN * * INITIAL OPERATIONS WITH SIMPLE BOUNDS * IF (KBF.GT.0) THEN DO 2 I = 1,NF IF ((IX(I).EQ.3.OR.IX(I).EQ.4) .AND. XU(I).LE.XL(I)) THEN XU(I) = XL(I) IX(I) = 5 ELSE IF (IX(I).EQ.5 .OR. IX(I).EQ.6) THEN XL(I) = X(I) XU(I) = X(I) IX(I) = 5 END IF 2 CONTINUE END IF M=IH(NF+1)-1 IF (KBF.GT.0) THEN CALL PCBS04(NF,X,IX,XL,XU,EPS9,KBF) CALL PYADC0(NF,N,X,IX,XL,XU,INEW) CALL MXVINE(IH(NF+1)-1,JH) END IF CALL PYPTSH(NF,MMAX,IH,JH,COL,S,XO,GO,WN11,WN12,GF,ITERM) MH=0 CALL MXVINS(M,0,JH(M+1)) IF (ITERM.LT.0) GO TO 11180 CALL OBJ(NF,X,F) NFV=NFV+1 CALL DOBJ(NF,X,GF) NFG=NFG+1 12140 CONTINUE CALL PF1HS2(NF,MH,MMAX,X,IX,S,HF,IH,JH,GO,GF,COL,WN11,WN12,XS, & F,ETA0,KBF,ITERM,ISYS) IF (ISYS.EQ.0) GO TO 12160 CALL DOBJ(NF,X,GO) NFG=NFG+1 GO TO 12140 12160 CONTINUE * * START OF THE ITERATION WITH TESTS FOR TERMINATION. * 11120 CONTINUE CALL PYTRCG(NF,NF,IX,GF,UMAX,GMAX,KBF,IOLD) IF (ABS(IPRNT).GT.1) & WRITE (6,'(1X,''NIT='',I5,2X,''NFV='',I5,2X,''NFG='',I5,2X, & ''F='', G16.9,2X,''G='',E10.3)') NIT,NFV,NFG,F,GMAX CALL PYFUT1(NF,F,FO,UMAX,GMAX,DMAX,TOLX,TOLF,TOLB,TOLG,KD, & NIT,KIT,MIT,NFV,MFV,NFG,MFG,NTESX,MTESX,NTESF,MTESF,ITES, & IRES1,IRES2,IREST,ITERS,ITERM) IF (ITERM.NE.0) GO TO 11180 IF (KBF.GT.0) CALL PYRMC0(NF,N,IX,GF,EPS8,UMAX,GMAX,RMAX,IOLD, & IREST) 11130 CONTINUE IF (IREST.GT.0) THEN CALL MXSSMI(NF,HF,IH) IDECF=0 LD=MIN(LD,1) IF (KIT.LT.NIT) THEN NRES=NRES+1 KIT = NIT ELSE ITERM=-10 IF (ITERS.LT.0) ITERM=ITERS-5 IF (GMAX.LE.1.0D 2*TOLG) ITERM=-ITERM END IF END IF IF (ITERM.NE.0) GO TO 11180 IF (KBF.GT.0) THEN CALL MXVINB(M,IX,JH) CALL PYTSCH(NF,IX,HF,IH,JH,KBF) END IF * * DIRECTION DETERMINATION * CALL PDSGM4(NF,MMAX,IX,GF,HF,IH,JH,S,XO,GO,XS,GS,IW,XMAX, & XDEL,GNORM,GNORMO,SNORM,FMIN,F,P,PP,ETA0,ETA2,DEL1,KD,KBF, & MOS1,MOS2,MOS3,IEST,IDECF,NDEC,NIT,NIN,ITERD,ITERM) * * TEST ON LOCALLY CONSTRAINED STEP AND PREPARATION OF STEPSIZE * SELECTION * IF (ITERD.LT.0) THEN ITERM=ITERD ELSE IF (SNORM.LE.0.0D 0) THEN IREST=MAX(IREST,1) ELSE IREST=0 END IF IF (IREST.EQ.0) THEN RMAX=XMAX/SNORM END IF END IF IF (ITERM.NE.0) GO TO 11180 IF (IREST.NE.0) GO TO 11130 IF (NIT.EQ.1) KIT=NIT CALL PYTRCS(NF,X,IX,XO,XL,XU,GF,GO,S,RO,FP,FO,F,PO,P,RMAX,ETA9, & KBF) IF (RMAX.EQ.0.0D 0) GO TO 11175 11150 CONTINUE CALL PS0G01(R,F,FO,PO,PP,XDEL,XDELO,XMAX,RMAX,SNORM, & BET1,BET2,GAM1,GAM2,EPS4,EPS5,KD,LD,IDIR,ITERS,ITERD, & MAXST,NRED,MRED,KTERS,MES1,MES2,MES3,ISYS) IF (ISYS.EQ.0) GO TO 11154 CALL MXUDIR(NF,R,S,XO,X,IX,KBF) CALL PCBS04(NF,X,IX,XL,XU,EPS9,KBF) CALL OBJ(NF,X,F) NFV=NFV+1 GO TO 11150 11154 CONTINUE KD=2 IF (ITERS.LE.0) THEN R=0.0D 0 F=FO P=PO CALL MXVCOP(NF,XO,X) IF (ITERS.LT.0) THEN ITERM=-6 IF (GMAX.LE.1.0D 2*TOLG) ITERM=-ITERM GO TO 11180 END IF IF (IDIR.EQ.0) IREST=MAX(IREST,1) LD=KD GO TO 11130 END IF IF (KD.GT.LD) THEN CALL DOBJ(NF,X,GF) NFG=NFG+1 14140 CONTINUE CALL PF1HS2(NF,MH,MMAX,X,IX,S,HF,IH,JH,GO,GF,COL,WN11,WN12,XS, & F,ETA0,KBF,ITERM,ISYS) IF (ISYS.EQ.0) GO TO 14160 CALL DOBJ(NF,X,GO) NFG=NFG+1 GO TO 14140 14160 CONTINUE END IF CALL PYTRCD(NF,X,IX,XO,GF,GO,R,F,FO,P,PO,DMAX,KBF,KD,LD,ITERS) 11175 CONTINUE IF (KBF.GT.0) CALL PYADC0(NF,N,X,IX,XL,XU,INEW) GO TO 11120 11180 CONTINUE IF (IPRNT.GT.1.OR.IPRNT.LT.0) & WRITE(6,'(1X,''EXIT FROM PNEC :'')') IF (IPRNT.NE.0) & WRITE (6,'(1X,''NIT='',I5,2X,''NFV='',I5,2X,''NFG='',I5,2X, & ''F='', G16.9,2X,''G='',E10.3,2X,''ITERM='',I3)') NIT,NFV,NFG, & F,GMAX,ITERM IF (IPRNT.LT.0) & WRITE (6,'(1X,''X='',5(G14.7,1X):/(3X,5(G14.7,1X)))') & (X(I),I=1,NF) RETURN END