C REMARK ON ALGORITHM 630, COLLECTED ALGORITHMS FROM ACM. C THIS WORK PUBLISHED IN TRANSACTIONS ON MATHEMATICAL SOFTWARE, C VOL. 11, NO. 2, PP. 103-119. PROGRAM TESTBB * C============================= S T A T U S ============================= C C SYSTEM DEPENDENCE: ISOLATED INTO 5 SMALL SUBROUTINES C C CONVERSION FOR SINGLE OR DOUBLE PRECISION: REQUIRED (SEE CONVRT) C C IGNORE LINES BEGINNING WITH "C!!!!" . C C THIS VERSION IS IN D O U B L E PRECISION. C!!!! THIS VERSION IS IN S I N G L E PRECISION. C C REVISION STATUS: DATE AUTHOR VERSION C C OCT. 30, 1988 A. BUCKLEY 1.1 C C======================== D E S C R I P T I O N ======================== C C THIS IS A ROUTINE PROVIDED TO TEST BBVSCG AND BBLNIR AFTER C INSTALLATION ON A PARTICULAR SYSTEM. IT CALLS BBVSCG AND BBLNIR C TO MINIMIZE A COLLECTION OF 10 TEST FUNCTIONS WHICH ARE PROVIDED C IN ZZFNS. C C IT ALSO SERVES AS A MODEL TO ILLUSTRATE THE USE OF SOME OF THE C FEATURES OF THE MINIMIZATION ALGORITHM IMPLEMENTATION. C C FOR AN EXAMPLE OF THE CODING OF A TEST FUNCTION, SEE THE ROUTINE C ZZFNS. TO SEE HOW TO CHANGE THE INTEGER COMMUNICATION PARAMETERS C (ESSENTIALLY LIKE ENUMERATED TYPES OF PASCAL), SEE THE EXAMPLES C BELOW AT FORTRAN LABEL 19. C C EACH FUNCTION IS MIMIMIZED SEVERAL TIMES. BOTH BBVSCG AND C BBLNIR ARE CALLED. TESTS INVOLVE ANALYTIC AND DIFFERENCED DERI- C VATIVES AND USE BOTH FORWARD AND REVERSE COMMUNICATION. BOTH THE C CONJUGATE GRADIENT AND QUASI-NEWTON CODES ARE TRIED AS WELL AS C THE NOCEDAL UPDATES. C C----- SUMMARY OF TESTS. C C 1. SIMPLE CALL TO BBVSCG, ANALYTIC DERIVATIVES, FORWARD CALLS. C 2. SIMPLE CALL TO BBVSCG, ANALYTIC DERIVATIVES, REVERSE CALLS. C 3. SIMPLE CALL TO BBVSCG, FINITE DIFFERENCES, FORWARD CALLS. C 4. SIMPLE CALL TO BBVSCG, DERIVATIVE TESTING, FORWARD CALLS. C C 5. DIRECT CALL TO BBLNIR, ANALYTIC, FORWARD, METH = 2. C 6. DIRECT CALL TO BBLNIR, ANALYTIC, REVERSE, METH = -2. C C 7. SIMPLE CALL TO BBVSCG, ANALYTIC DERIVATIVES, NOCEDAL UPDATES. C C======================= E N T R Y P O I N T S ======================= C C ONLY THE NATURAL ENTRY POINT TESTBB C C======================== S U B R O U T I N E S ======================== C C BBVSCG C C BBLNIR ( AND ENTRY POINTS THEREIN ) C C BBDFLT C C ZZEVAL ( AND ENTRY POINT ZZECHK THEREIN ) C C ZZFNS ( ALSO THROUGH ENTRY POINTS ZZFSET AND ZZFPAR) C ZZDATE C ZZDTTM C ZZLENG C ZZLCUC C ZZSHFT C ZZTIME C C ZZSECS C C BARD70, BIGGS6, BOX663, CRGLVY, ENGVL2! A SUBSET OF A FULL COLL- C PENAL1, PENAL2, PWSING, ROSENB, SCHMVT! ECTION OF TEST FUNCTIONS C C========================= P A R A M E T E R S ========================= * * INTEGER NINTS, NLOGS, NREALS, NTRACF PARAMETER ( NINTS = 14, NLOGS = 32, NREALS = 2, NTRACF = 15 ) * INTEGER XDRVMD, XNORM, XSCALE, XLTRCU PARAMETER ( XDRVMD = 1, XNORM = 2, XSCALE = 3, XLTRCU = 4 ) * INTEGER XETRCU, XPTRCU, XTTRCU PARAMETER ( XETRCU = 5, XPTRCU = 6, XTTRCU = 7 ) * INTEGER XMETH, XQUADN, XALPS1, XSCGMM PARAMETER ( XMETH = 8, XQUADN = 9, XALPS1 = 10, XSCGMM = 11 ) * INTEGER XHTEST, XUPDTT, XSTSTP PARAMETER ( XHTEST = 12,XUPDTT = 13, XSTSTP = 14 ) * INTEGER XTRACE PARAMETER ( XTRACE = 1 ) * INTEGER XTRF, XTRG, XTTRCE, XTRTST PARAMETER ( XTRF = 16, XTRG = 17, XTTRCE = 18, XTRTST = 19 ) * INTEGER XGRAD, XPOINT, XTGRAD PARAMETER ( XGRAD = 20, XPOINT = 21, XTGRAD = 22 ) * INTEGER XTSTEP, XTSHXG, XTFUNC, XRELF PARAMETER ( XTSTEP = 23,XTSHXG = 24, XTFUNC = 25, XRELF = 26 ) * INTEGER XRELG, XFQUAD, XDIAGL PARAMETER ( XRELG = 27, XFQUAD = 28, XDIAGL = 29 ) * INTEGER XSHNNO, XFRMRS, XFRCEF PARAMETER ( XSHNNO = 30,XFRMRS = 31, XFRCEF = 32 ) * INTEGER XRO, XBETA PARAMETER ( XRO = 1, XBETA = 2 ) * C-----DEFINE AMOUNT OF WORKING STORAGE. * INTEGER EXTRA, MXN, LWORK PARAMETER ( EXTRA = 20, MXN = 40, LWORK = (MXN*(MXN+7))/2 ) * C-----SET UNITS FOR TERMINAL INPUT/OUTPUT FOR CONTROL. C OUTPT IS FOR OUTPUT FROM TEST; SAVED ON A FILE. * 11 INTEGER TRMINP, TRMOUT, OUTPT PARAMETER ( TRMINP = 5, TRMOUT = 6, OUTPT = 8 ) * C-----DEFINE TOTAL NUMBER OF ALLOWABLE FUNCTION ARGUMENTS. * * INTEGER FNO PARAMETER ( FNO = 10 ) * C-----DEFINE TOTAL NUMBER OF TESTS AND PROBLEMS. * INTEGER TESTS, NPROBS, TTESTS PARAMETER ( TESTS = 7, NPROBS = 10, TTESTS = TESTS*NPROBS ) * C-----FOR ZZEVAL: NUMBER OF FUNCTION EVALUATIONS. * INTEGER MAX PARAMETER ( MAX = 300 ) * C-----FOR ZZPRNT: FIRST AND LAST FUNCTION VALUES PRINTED. * INTEGER PFREQ PARAMETER ( PFREQ = -1000 ) * C-----FOR ZZTERM: ACCURACY REQUIREMENT. * DOUBLE PRECISION ACC C!!!! REAL ACC PARAMETER ( ACC = 5.D-4 ) * C-----FOR BBLNIR: ESTIMATE OF FUNCTION REDUCTION. * DOUBLE PRECISION DECRF C!!!! REAL DECRF PARAMETER ( DECRF = -1.0D0 ) * C---- MISCELLANEOUS CONSTANTS (NOT ALL USED). * DOUBLE PRECISION ZERO, ONE, TWO, THREE C!!!! REAL ZERO, ONE, TWO, THREE PARAMETER ( ZERO = 0D0, ONE = 1D0, TWO = 2D0, THREE = 3D0) * DOUBLE PRECISION FOUR, FIVE, SIX, SEVEN C!!!! REAL FOUR, FIVE, SIX, SEVEN PARAMETER ( FOUR = 4D0, FIVE = 5D0,SIX = 6D0, SEVEN = 7D0) * DOUBLE PRECISION EIGHT, NINE, TEN C!!!! REAL EIGHT, NINE, TEN PARAMETER ( EIGHT = 8D0, NINE = 8D0, TEN = 10D0 ) * DOUBLE PRECISION TENTH, FIFTH, HALF C!!!! REAL TENTH, FIFTH, HALF PARAMETER ( TENTH = .1D0, FIFTH = .2D0, HALF = .5D0 ) * DOUBLE PRECISION RPT9, RPT8, RD29 C!!!! REAL RPT9, RPT8, RD29 PARAMETER ( RPT9 = .9D0, RPT8 = .8D0, RD29 = 1D0/29D0 ) DOUBLE PRECISION R11, R12, R13, R14 C!!!! REAL R11, R12, R13, R14 PARAMETER ( R11 = 11D0, R12 = 12D0, R13 = 13D0,R14 = 14D0) * DOUBLE PRECISION R15, R16, R17, R18 C!!!! REAL R15, R16, R17, R18 PARAMETER ( R15 = 15D0, R16 = 16D0, R17 = 17D0,R18 = 18D0) * DOUBLE PRECISION R19, R20, R25, R29 C!!!! REAL R19, R20, R25, R29 PARAMETER ( R19 = 19D0, R20 = 20D0, R25 = 25D0,R29 = 29D0) * DOUBLE PRECISION R32, R36, R40, R42 C!!!! REAL R32, R36, R40, R42 PARAMETER ( R32 = 32D0, R36 = 36D0, R40 = 40D0,R42 = 42D0) * DOUBLE PRECISION R45, R49 C!!!! REAL R45, R49 PARAMETER ( R45 = 45D0, R49 = 49D0 ) * DOUBLE PRECISION R50, R56, R84, R90 C!!!! REAL R50, R56, R84, R90 PARAMETER ( R50 = 50D0, R56 = 56D0, R84 = 84D0,R90 = 90D0) * DOUBLE PRECISION R100, R180, R200 C!!!! REAL R100, R180, R200 PARAMETER ( R100 = 100D0, R180 = 180D0, R200 = 200D0 ) * DOUBLE PRECISION R256, R360, R400 C!!!! REAL R256, R360, R400 PARAMETER ( R256 = 256D0, R360 = 360D0, R400 = 400D0 ) * DOUBLE PRECISION R600, R681, R991 C!!!! REAL R600, R681, R991 PARAMETER ( R600 = 600D0, R681 = 681D0, R991 = 991D0 ) * DOUBLE PRECISION R1162, R2324 C!!!! REAL R1162, R2324 PARAMETER ( R1162 = 1162D0, R2324 = 2324D0 ) * DOUBLE PRECISION R10000, R40000 C!!!! REAL R10000, R40000 PARAMETER ( R10000 = 10000D0, R40000 = 40000D0 ) DOUBLE PRECISION R1PD6, R2PDM6 C!!!! REAL R1PD6, R2PDM6 PARAMETER ( R1PD6 = 1D6, R2PDM6 = 2D-6 ) * DOUBLE PRECISION RP04, RP01, R1PZ1 C!!!! REAL RP04, RP01, R1PZ1 PARAMETER ( RP04 = 4D-2, RP01 = .01D0, R1PZ1 = 1.0001D0 ) * DOUBLE PRECISION R1P2 C!!!! REAL R1P2 PARAMETER ( R1P2 = 1.2D0 ) * DOUBLE PRECISION R1P5, R2P5, R2P625 C!!!! REAL R1P5, R2P5, R2P625 PARAMETER ( R1P5 = 1.5D0, R2P5 = 2.5D0, R2P625 = 2.625D0 ) * DOUBLE PRECISION R10P1, R19P8, R20P2 C!!!! REAL R10P1, R19P8, R20P2 PARAMETER ( R10P1 = 10.1D0,R19P8 = 19.8D0,R20P2 = 20.2D0 ) * DOUBLE PRECISION R2D3, R4D3, R7D3 C!!!! REAL R2D3, R4D3, R7D3 PARAMETER ( R2D3 = 2D0/3D0,R4D3 = 4D0/3D0,R7D3 = 7D0/3D0 ) * DOUBLE PRECISION R2P25 C!!!! REAL R2P25 PARAMETER ( R2P25 = 2.25D0 ) * C======================= D E C L A R A T I O N S ======================= * C---- DEFINE THE ARRAYS FOR REVISING CONTROL PARAMETERS. * INTEGER INTS(NINTS) LOGICAL LOGS(NLOGS) * DOUBLE PRECISION REALS(NREALS) C!!!! REAL REALS(NREALS) * C---- DEFINE THE IW, RW, DW ARRAYS AND THE INNER PRODUCT ROUTINE. * INTEGER IW(EXTRA) REAL RW(EXTRA) DOUBLE PRECISION DW(EXTRA), ZZINNR EXTERNAL ZZINNR * C---- PLACES TO HOLD SOME STATISTICS AND STUFF FROM THE TESTS. * INTEGER ICNTS(TTESTS), FCNTS(TTESTS), FNCT, GRCT, ITCT * DOUBLE PRECISION FVALS(TTESTS), ACCTIM C!!!! REAL FVALS(TTESTS), ACCTIM * INTEGER DIM(NPROBS), IFNC(NPROBS), INDX(NPROBS), COMPNT(NPROBS) * C---- VALUES USED TO REDEFINE COMMUNICATION PARAMETERS. * INTEGER ANAL, DIFF, TEST, DOFG INTEGER NORMFG, NORMAL, RCSTRT, RCRPT, RCNOFG, PSTHRU INTEGER DONE, RCF, RCFG, RCG, NOSTOR, IPMIN, IPUNDF INTEGER BDMETH, LSFAIL, NODESC, XSFUNC, RABORT, USERV, PSBACK * C---- VARIOUS DECLARATIONS NEEDED TO RUN THE TESTS. * INTEGER N, STATUS, METH, ERROR, I, M, ID, IX ,UNIT INTEGER IG, IH, DERVMD, DINDX, DCOMP, CASE, CONTRL, FROM, TO INTEGER ITERS, FUNCS, FREQ, PDONE, TFNCS, TITERS * DOUBLE PRECISION X(MXN), G(MXN), WORK(LWORK+EXTRA), DERERR(NPROBS) C!!!! REAL X(MXN), G(MXN), WORK(LWORK+EXTRA), DERERR(NPROBS) * DOUBLE PRECISION FX, DERR, TIME, AVERRS(NPROBS), AVERR C!!!! REAL FX, DERR, TIME, AVERRS(NPROBS), AVERR * DOUBLE PRECISION FARG(FNO), RPAR1(NPROBS), RPAR2(NPROBS), ACCT C!!!! REAL FARG(FNO), RPAR1(NPROBS), RPAR2(NPROBS), ACCT * CHARACTER*52 TITLE(TESTS) CHARACTER*41 DATE CHARACTER*4 QUITS * EXTERNAL ZZFNS * C=============================== S A V E =============================== C C THERE ARE NO SAVE VARIABLES. C C============================= C O M M O N ============================= C C THERE ARE NO COMMON BLOCKS. C C=============================== D A T A =============================== * C---- DIMENSIONS OF THE TEST PROBLEMS. * DATA DIM / 3, 6, 3, 4, 3, 10, 10, 40, 2, 3 / * C---- NUMBERS OF THE TEST PROBLEMS. * DATA IFNC / 26, 11, 16, 27, 25, 40, 41, 3, 1, 24 / * C---- ARGUMENTS NEEDED FOR THE TEST PROBLEMS. * DATA RPAR1 / 0.D0, 13.D0, 10.D0, 2*0.D0, 2*1.D-5, 3*0D0/ * DATA RPAR2 / 5*0.D0, 2*1.D0, 3*0.D0 / * DATA FARG / FNO * 1.D0 / * C---- OUTPUT IDENTIFICATION. * DATA TITLE - -/ ' CALL BBVSCG, ANALYTIC MODE, FORWARD CALLS.' , - ' CALL BBVSCG, ANALYTIC MODE, REVERSE CALLS.' , - ' CALL BBVSCG, DIFFERENCING, FORWARD CALLS.' , - ' CALL BBVSCG, TESTING MODE, FORWARD CALLS.' , - ' CALL BBLNIR, ANALYTIC MODE, FORWARD CALLS; METH= 2.' , - ' CALL BBLNIR, ANALYTIC MODE, REVERSE CALLS; METH=-2.' , - ' CALL BBVSCG, ANALYTIC MODE, NOCEDAL UPDATES.' / * C========================== E X E C U T I O N ========================== * C WE BEGIN BY REDEFINING THE INTEGER (ENUMERATED) CONTROL PARAMETERS C USED FOR INTER-PROGRAM COMMUNICATION. THE VALUES SET HERE ARE C DIFFERENT THAN THE DEFAULTS. NORMALLY THEY ARE NOT CHANGED UNLESS C THE DEFAULTS CONFLICT WITH THE USER'S PROGRAMS. THEY ARE CHANGED C HERE ONLY FOR ILLUSTRATIVE PURPOSES. OF COURSE, WE MUST BE CAREFUL C TO USE THE VALUES OF THE CODES AS WE DEFINE THEM. SEE, FOR EXAMPLE C LABELS 5200 AND 5333. * C ENTRY BBLDDF ( ANAL, DIFF, TEST, SFIRST ) 19 ANAL = 2 DIFF = 3 TEST = 4 CALL BBLDDF ( ANAL, DIFF, TEST, -1 ) * C ENTRY BBLFDF ( DOF, DOG, DOFG, NONE ) DOFG = 3 CALL BBLFDF ( 1, 2, DOFG, 0 ) C REDEFINE CODES IN TEST ROUTINE ZZFNS WHICH IS PROVIDED AS WELL. CALL ZZFFDF ( 1, 2, DOFG, 0 ) * C ENTRY BBLIDF ( SNRMFG, SNORML, SRCSTR, SRCRPT, SRCNFG, SPSTHR ) CALL BBLIDF ( -1, -2, 1, 2, 0, 99 ) NORMFG = -1 NORMAL = -2 RCSTRT = 1 RCRPT = 2 RCNOFG = 0 PSTHRU = 99 * C ENTRY BBLSDF ( SDONE, SRCF, SRCFG, SRCG, SNSTOR, SIPMIN, C - SIPUNF, SBDMTH, SLSFAL, SNODSC, SXSFNC, SRABRT, SUSERV, SPSBCK) CALL BBLSDF (0, 1, 2, 3, -1, -2, -3, -4, -5, -6, -7, -8, -9, 99) DONE = 0 RCF = 1 RCFG = 2 RCG = 3 NOSTOR = -1 IPMIN = -2 IPUNDF = -3 BDMETH = -4 LSFAIL = -5 NODESC = -6 XSFUNC = -7 RABORT = -8 USERV = -9 PSBACK = 99 * C ENTRY BBLRDF ( SOK, SABORT, SLIMIT, SNOF, SNOG, SNOFG ) CALL BBLRDF ( 0, 1, 2, 3, 4, 5 ) * C---- INITIALIZE TIMING. * CALL ZZSECS (TIME) ACCTIM = ZERO * OPEN ( OUTPT, FILE = 'RESULTS' ) CALL ZZDTTM (DATE) WRITE ( OUTPT, 99993 ) ' STARTING TEST AT ', DATE( 1:10), ' ON ', - DATE(12:41) * C INITIALIZE COUNTS. * ERROR = 0 PDONE = 0 TFNCS = 0 TITERS = 0 FREQ = PFREQ * DO 20 I = 1, TTESTS FVALS(I) = ZERO FCNTS(I) = 0 ICNTS(I) = 0 20 CONTINUE * C---- ASK FOR CHANGES TO CONTROL PARAMETERS. EOF SUFFICIENT IF C NONE REQUIRED. * CALL BBRVAL ( TRMOUT, TRMINP ) * C REDEFINE THE OUTPUT UNIT BY GETTING THE CURRENT VALUES C FROM BBVALS AND RESETTING IT BY CALLING BBSVAL. THIS MAKES C SURE *ALL* OUTPUT GOES TO UNIT OUTPT. THIS IS GENERALLY C NOT NECESSARY. OUTPUT IS SENT TO A DEFAULT VALUE. * UNIT = OUTPT 22 CALL BBVALS ( INTS, LOGS, REALS ) INTS(XPTRCU) = UNIT INTS(XLTRCU) = UNIT INTS(XETRCU) = UNIT CALL BBSVAL ( INTS, LOGS, REALS ) * C---- DO MINIMIZATIONS. RUN EACH OF THE TEST TYPES. * FROM = 0 TO = NPROBS CONTRL = 0 * DO 9000 M = 1,TESTS * C WRITE TITLE, UNLESS NO OUTPUT REQUESTED. * IF ( FREQ .NE. 0 ) WRITE ( UNIT, 99999 ) M, TITLE(M) * C HERE WE HAVE A CHANCE TO CHOOSE A SUBSET OF THE PROBLEMS C TO TEST. * 50 IF ( CONTRL .NE. -3 ) THEN CONTRL = -3 WRITE(TRMOUT,*) ' CONTROL: 0 QUIT' WRITE(TRMOUT,*) ' -1 SKIP TO NEXT SET,' WRITE(TRMOUT,*) ' -2 FINISH THIS SET' WRITE(TRMOUT,*) ' -3 (OR EOF) FINISH FULL RUN' WRITE(TRMOUT,*) ' N > 0 DO PROBLEM #N' READ(TRMINP,'(BN,I2)', END=59 ) CONTRL ENDIF * 59 IF ( CONTRL .GT. 0 ) THEN FROM = CONTRL TO = CONTRL ELSE IF ( CONTRL .EQ. 0 ) THEN GOTO 10000 ELSE IF ( CONTRL .EQ. -1 ) THEN FROM = 0 TO = NPROBS GOTO 8500 ELSE IF ( CONTRL .EQ. -2 ) THEN FROM = MOD(TO,NPROBS)+1 TO = NPROBS ELSE IF ( CONTRL .EQ. -3 ) THEN FROM = MOD(TO,NPROBS)+1 TO = NPROBS ENDIF * C START TIMING. * CALL ZZSECS(TIME) ACCTIM = ACCTIM - TIME * C REPEAT FOR EACH TEST FUNCTION SELECTED. * DO 6000 I = FROM,TO * PDONE = PDONE + 1 IF (PDONE .GT. TTESTS) THEN WRITE ( UNIT, * ) ' TOO MANY TESTS: STOPPING.' GOTO 90000 ENDIF * IF ( FREQ .NE. 0 ) WRITE ( UNIT, 99994 ) I * C ---SET FUNCTION NUMBER IN ZZFNC AND SET FUNCTION ARGUMENTS. C THIS ILLUSTRATES HOW TO DEFINE THE AMOUNT OF EXTRA C STORAGE AVAILABLE. IT IS ONLY NEEDED FOR TEST 7. * IF ( I .EQ. 7 ) THEN CALL ZZFSET ( IFNC(I), EXTRA ) ELSE CALL ZZFSET ( IFNC(I), 1 ) ENDIF * FARG(1) = RPAR1(I) FARG(2) = RPAR2(I) * CALL ZZFPAR ( FARG ) * C ---SET DIMENSION. * N = DIM(I) * C ---SET STARTING POINT. * * GOTO ( 100,200,300,400,500,600,700,800,900,1000) I * 100 X(1) = ONE X(2) = ONE X(3) = ONE * GOTO 5000 * 200 X(1) = ONE X(2) = TWO X(3) = ONE X(4) = ONE X(5) = ONE X(6) = ONE * GOTO 5000 * 300 X(1) = ZERO X(2) = TEN X(3) = R20 * GOTO 5000 * 400 X(1) = ONE X(2) = TWO X(3) = TWO X(4) = TWO * GOTO 5000 * 500 X(1) = ONE X(2) = TWO X(3) = ZERO * GOTO 5000 * 600 DO 650 ID = 1,N X(ID) = ID 650 CONTINUE * GOTO 5000 * 700 DO 750 ID = 1,N X(ID) = HALF 750 CONTINUE * GOTO 5000 * 800 DO 850 ID = 1,N/4 X(4*ID-3) = THREE X(4*ID-2) = - ONE X(4*ID-1) = ZERO X(4*ID ) = ONE 850 CONTINUE * GOTO 5000 * 900 X(1) = -R1P2 X(2) = ONE * GOTO 5000 * 1000 X(1) = HALF X(2) = HALF X(3) = HALF * GOTO 5000 * C ---SET UP CALLS TO MINIMIZE. TESTS 5 AND 6 CALL C BBLNIR DIRECTLY. THE OTHERS CALL BBVSCG. * 5000 GOTO ( 5100, 5200, 5300, 5400, 5500, 5600, 5425 ) M * C TEST 1 THE SIMPLEST CALL. 5100 DERVMD = ANAL STATUS = NORMAL GOTO 5450 * C TEST 2 USING REVERSE COMMUNICATION. 5200 DERVMD = DIFF STATUS = RCSTRT * C MUST INITIALIZE BEFORE CALL TO ZZEVAL. CALL BBDFLT ( FREQ, MAX ) CASE = DOFG CALL ZZEVAL ( ZZFNS, N, X, FX, G, CASE, IW, RW, DW ) * GOTO 5450 * C TEST 3 HERE WE SEE HOW TO REDEFINE THE METHOD OF COMPUTING C DERIVATIVES. FOR CONVENIENCE, A SET OF NAMED INTEGER C INDICES IS PROVIDED FOR ACCESSING THE APPROPRIATE C ENTRIES OF EACH ARRAY. THESE ARE DOCUMENTED IN BBVALS. * 5300 CONTINUE * 5333 DERVMD = DIFF CALL BBVALS ( INTS, LOGS, REALS ) INTS(XDRVMD) = DERVMD CALL BBSVAL ( INTS, LOGS, REALS ) STATUS = NORMAL * GOTO 5450 * C TEST 4 TEST THE DERIVATIVE CALCULATIONS. 5400 DERVMD = TEST CALL BBVALS ( INTS, LOGS, REALS ) INTS(XDRVMD) = DERVMD CALL BBSVAL ( INTS, LOGS, REALS ) STATUS = NORMAL * GOTO 5450 * C TEST 7 HERE WE USE THE SAME TECHNIQUE TO RESET C THE MINIMIZATION METHOD, THE UPDATE C STRATEGY AND THE DERIVATIVE MODE. THIS C TESTS NOCEDAL'S UPDATE STRATEGY. 5425 DERVMD = ANAL * CALL BBVALS ( INTS, LOGS, REALS ) INTS(XMETH ) = 2 INTS(XUPDTT) = 2 INTS(XDRVMD) = DERVMD CALL BBSVAL ( INTS, LOGS, REALS ) STATUS = NORMAL * GOTO 5450 * 5450 CONTINUE * ITERS = FREQ FUNCS = MAX ACCT = ACC * 5460 CALL BBVSCG ( ZZFNS, N, X, FX, G, ACCT, STATUS, - ITERS, FUNCS, WORK, LWORK ) * C CHECK FOR REVERSE COMMUNICATION TEST IF ( M .EQ. 2 .AND. STATUS .EQ. RCFG ) THEN * C RE-EVALUATE F AND G. CASE = DOFG CALL ZZEVAL ( ZZFNS, N, X, FX, G, CASE, IW, RW, DW ) * C AND RE-ENTER BBVSCG. STATUS = RCRPT * GOTO 5460 * ELSE IF ( M .EQ. 4 ) THEN C GET TESTING INFO. CALL ZZECHK ( DERR, AVERR, DCOMP, DINDX ) * DERERR(I) = DERR AVERRS(I) = AVERR COMPNT(I) = DCOMP INDX (I) = DINDX * ENDIF * GOTO 5900 * C TEST 5 AGAIN WE SEE HOW TO REDEFINE THE METHOD OF COMPUTING C DERIVATIVES. FOR CONVENIENCE, A SET OF NAMED INTEGER C INDICES IS PROVIDED FOR ACCESSING THE APPROPRIATE C ENTRIES OF EACH ARRAY. THESE ARE DOCUMENTED IN BBVALS. C HERE WE ARE CALLING BBLNIR DIRECTLY. * 5500 DERVMD = ANAL CALL BBVALS ( INTS, LOGS, REALS ) INTS(XDRVMD) = DERVMD CALL BBSVAL ( INTS, LOGS, REALS ) STATUS = NORMAL C REDEFINE THE ALGORITHM STRATEGY. METH = 2 * GOTO 5650 * C TEST 6 LIKE TEST 5, BUT WITH REVERSE COMMUNICATION AND C A DIFFERENT METHOD. 5600 DERVMD = ANAL STATUS = RCSTRT METH = -2 * GOTO 5650 * C DIRECT CALL TO BBLNIR. FIRST GET CURRENT SETTINGS. 5650 CALL BBVALS ( INTS, LOGS, REALS ) * C THEN INITIALIZE ZZEVAL. CALL ZZESET ( LOGS(XTRF), LOGS(XTRG), LOGS(XTRTST), - INTS(XETRCU) ) CALL ZZESRT ( INTS(XSCALE), DERVMD, MAX ) C THEN INITIALIZE ZZPRNT. CALL ZZP1ST ( INTS(XPTRCU), LOGS(XGRAD), LOGS(XPOINT), - FREQ ) CALL ZZP2ST ( INTS(XPTRCU), LOGS(XGRAD), LOGS(XPOINT), 0 ) * C THEN INITIALIZE ZZTERM. QUITS = 'FFFF' IF ( LOGS(XTGRAD) ) QUITS(1:1) = 'T' IF ( LOGS(XTSTEP) ) QUITS(2:2) = 'T' IF ( LOGS(XTSHXG) ) QUITS(3:3) = 'T' IF ( LOGS(XTFUNC) ) QUITS(4:4) = 'T' CALL ZZTSET (INTS(XNORM), QUITS, LOGS(XTTRCE), INTS(XTTRCU)) * C INITIALIZE BBLNIR THROUGH ENTRY POINT. HERE ALL THE C VALUES ARE TAKEN FROM THOSE STORED IN BBVALS. IN C MOST APPLICATIONS, VALUES WOULD BE DIRECTLY INPUT BY C THE USER FOR THOSE VALUES HE WISHED TO CHANGE. CALL BBLSET ( METH, INTS(XQUADN), INTS(XALPS1),INTS(XSTSTP), - INTS(XSCGMM), INTS(XHTEST),INTS(XUPDTT), - REALS(XRO), REALS(XBETA), - LOGS(XFQUAD), LOGS(XDIAGL),LOGS(XSHNNO), - LOGS(XFRMRS), LOGS(XFRCEF), - LOGS(XRELF), LOGS(XRELG), - INTS(XLTRCU), LOGS(XTRACE) ) * C DEFINE THE WORKING ARRAY SIZES. ID = 1 IX = ID + N IG = IX + N IH = IG + N * C INITIAL FUNC/GRAD VALUES FOR REVERSE COMMUNICATION. IF ( STATUS .EQ. RCSTRT ) THEN CASE = DOFG CALL ZZEVAL ( ZZFNS,N,X,FX, G, CASE, IW, RW, DW ) ENDIF * 5690 ACCT = ACC CALL BBLNIR(ZZFNS, N, X, FX, DECRF, G, ACCT, STATUS, ZZINNR, - WORK(ID), WORK(IX), WORK(IG), WORK(IH), LWORK-3*N, - IW, RW, DW ) * C CHECK FOR REVERSE COMMUNICATION TEST IF ( M .EQ. 6 .AND. STATUS .EQ. RCRPT ) THEN * C RE-EVALUATE F AND G. CASE = DOFG CALL ZZEVAL( ZZFNS, N, X, FX, G, CASE, IW, RW, DW ) * C AND RE-ENTER BBVSCG. STATUS = RCRPT GOTO 5690 * ENDIF * C ALL TESTS: ADD NUMBER OF ERRORS. * 5900 IF ( STATUS .NE. DONE ) THEN ERROR = ERROR + 1 ENDIF * C GET STATISTICAL COUNTS. CALL ZZEGET( FNCT, GRCT, ACCTIM ) CALL ZZPGET( ACCTIM, ITCT ) * FVALS(PDONE) = FX FCNTS(PDONE) = FNCT ICNTS(PDONE) = ITCT * TFNCS = TFNCS + FNCT TITERS = TITERS + ITCT * IF ( FREQ .NE. 0 ) WRITE ( UNIT, 99997 ) STATUS * 6000 CONTINUE CALL ZZSECS(TIME) ACCTIM = ACCTIM + TIME * 8000 IF ( M .EQ. 4 ) THEN IF ( FREQ .NE. 0 ) WRITE ( UNIT, 99998 ) - (DERERR(I),COMPNT(I),INDX(I),AVERRS(I),I = 1,NPROBS) ENDIF * 8500 IF ( TO .NE. NPROBS ) GOTO 50 * 9000 CONTINUE * 10000 WRITE ( UNIT, 99991 ) WRITE ( UNIT, 99992 ) (I,ICNTS(I),FVALS(I),FCNTS(I),I=1,PDONE) * CALL ZZSECS (TIME) WRITE ( UNIT, 99995 ) ACCTIM, TIME WRITE ( UNIT, 99996 ) PDONE, ERROR, TFNCS, TITERS WRITE ( TRMOUT, * ) ' ' WRITE ( TRMOUT, * ) ' ' WRITE ( TRMOUT, * ) ' ' WRITE ( TRMOUT, 99995 ) ACCTIM, TIME * CALL ZZDTTM (DATE) WRITE( TRMOUT, 99993 ) ' TEST ENDED AT ', DATE(1:10), ' ON ', - DATE(12:41) WRITE( UNIT, 99993 ) ' TEST ENDED AT ', DATE(1:10), ' ON ', - DATE(12:41) * C=============================== E X I T =============================== * 90000 STOP * C============================ F O R M A T S ============================ * 99991 FORMAT (// ' RN ITS FUNCT. VALUE FNS |', - ' RN ITS FUNCT. VALUE FNS |', - ' RN ITS FUNCT. VALUE FNS'/ - ' ------------------------|', - '-------------------------|', - '------------------------') * 99992 FORMAT ( (2(I3,I4,1X,E12.6,I4,' |'), (I3,I4,1X,E12.6,I4) )) * 99993 FORMAT ( 4A // ) * 99994 FORMAT ( /' FUNCTION #', I2 / ) * 99995 FORMAT (/' TIME USED WAS ', F12.3, ' SECONDS WITHOUT PROMPTS;'/ - ' ', F12.3, ' TOTAL.' ) * 99996 FORMAT ( // ' *****TEST FINISHED**** PROBLEMS DONE ',I3, - '; NUMBER OF ERRORS IS ',I2,'.' / - ' TOTAL FUNCTION CALLS = ',I4, - ' TOTAL ITERATIONS = ', I4 // ) * 99997 FORMAT ( /' ************* RUN COMPLETE, STATUS = ', I3, '.'/ ) * 99998 FORMAT ( //' TESTING MODE DERIVATIVE ESTIMATION ERRORS' // - ' MAX ERROR COMPONENT ITERATE AV. DECIMALS ' // - 10 ( E10.2, I7, 7X, I5, F9.2 /) ) * 99999 FORMAT ( '1 BEGINNING RUN #', I1, ':', A ) * C================================ E N D ================================ * END SUBROUTINE ZZDATE (CHDATE) * C============== A R G U M E N T D E C L A R A T I O N S ============== * CHARACTER *(*) CHDATE * C============================= S T A T U S ============================= C C SYSTEM DEPENDENCE: SYSTEM ROUTINE FOR DATE. C C THE CURRENT VERSION MAY BE PROCESSED TO PRODUCE C CODE FOR ONE OF THE FOLLOWING: C C THIS VERSION IS FOR SUN4 C CONVERSION FOR SINGLE OR DOUBLE PRECISION: NOT REQUIRED C C REVISION STATUS: DATE AUTHOR VERSION C C JUN. 16, 1988 A. BUCKLEY 1.2 C C======================== D E S C R I P T I O N ======================== C C THIS ROUTINE MUST CALL A SYSTEM ROUTINE TO GET THE CURRENT DATE. C ZZDATE MUST RETURN THIS DATE IN THE CHARACTER VARIABLE CHDATE C IN THE FORM C C (YY+MM+DD) C C AS REQUIRED BY THE ROUTINE ZZDTTM. CHDATE MUST OF LENGTH 10. C ONLY THE 6 CHARACTERS YY MM DD ARE ACTUALLY USED. THE OTHERS C CAN BE ANYTHING, I.E. ONLY THE POSITION OF THE YY MM DD MATTERS. C C THIS MUST CALL A SYSTEM ROUTINE TO GET THE DATE. C TO IMPLEMENT THIS ON ANOTHER SYSTEM, ONE MAY EITHER C C (A) INCORPORATE AN ALTERNATE VERSION OF ZZDATE; C C (B) USE THE "DUMMY" VERSION OF THIS ROUTINE WITH THE SINGLE C EXECUTABLE STATEMENT CHDATE='( + + )', IN WHICH CASE NO C DATE INFORMATION WILL APPEAR IN THE OUTPUT. C C======================= E N T R Y P O I N T S ======================= C C ONLY THE NATURAL ENTRY POINT ZZDATE C C======================== S U B R O U T I N E S ======================== C C SYSTEM DATE ROUTINE. C C========================= P A R A M E T E R S ========================= C C THERE ARE NO PARAMETERS DEFINED. C C================= L O C A L D E C L A R A T I O N S ================= * CHARACTER * 24 UNXDAT * CHARACTER * 3 NAME (12), TEMP * INTEGER I * C=============================== S A V E =============================== * SAVE NAME * C======================= E Q U I V A L E N C E S ======================= C C THERE ARE NO EQUIVALENCES. C C============================= C O M M O N ============================= C C THERE ARE NO COMMON BLOCKS. C C=============================== D A T A =============================== * DATA NAME /'JAN','FEB','MAR','APR','MAY','JUN', - 'JUL','AUG','SEP','OCT','NOV','DEC' / * C========================== E X E C U T I O N ========================== * C---- INITIALIZE CHDATE * CHDATE = '( + + )' * CALL FDATE(UNXDAT) * CHDATE(2:3) = UNXDAT(23:24) * CHDATE(8:9) = UNXDAT(9:10) * TEMP = UNXDAT(5:7) CALL ZZLCUC(TEMP) * DO 100 I = 1,12 * IF ( TEMP .EQ. NAME(I) ) THEN * WRITE ( CHDATE(5:6), '(I2.2)' ) I * GOTO 90000 * ENDIF * 100 CONTINUE * C=============================== E X I T =============================== * 90000 RETURN * C============================ F O R M A T S ============================ C C THERE ARE NONE. C C================================ E N D ================================ * END SUBROUTINE ZZDTTM ( CHDATE ) * C============== A R G U M E N T D E C L A R A T I O N S ============== * CHARACTER *(*) CHDATE * C============================= S T A T U S ============================= C C SYSTEM DEPENDENCE: NONE C C CONVERSION FOR SINGLE OR DOUBLE PRECISION: NOT REQUIRED C C REVISION STATUS: DATE AUTHOR VERSION C C JUN. 1, 1985 A. BUCKLEY 1.0 C C======================== D E S C R I P T I O N ======================== C C THIS ROUTINE RETURNS IN CHDATE A 41-CHARACTER DATE OF THE FORM C GIVEN IN MODEL(BELOW). IT USES THE TIME AND DATE AS OBTAINED C FROM THE OPERATING SYSTEM (VIA THE ROUTINES ZZTIME AND ZZDATE) C AND CONVERTS THEM TO THE FORM OF THE MODEL GIVEN BELOW. C IT ASSUMES THAT THE ROUTINES ZZTIME AND ZZDATE RETURN 10 C CHARACTER STRINGS, RESPECTIVELY, OF THE FORM: C C TIME: (HH+MM+SS) C DATE: (YY+MM+DD) C C NOTE THAT EXCESS BLANKS IN THE DATE ARE ELIMINATED. C IF CHDATE IS MORE THAN 41 CHARACTERS IN LENGTH, ONLY THE C LEFTMOST 41 WILL BE ALTERED. IF IT IS LESS THAN 41 IN C LENGTH, ONLY THE LEFTMOST CHARACTERS OF THE DATE WILL BE C RETURNED. C C======================= E N T R Y P O I N T S ======================= C C ONLY THE NATURAL ENTRY POINT ZZDTTM C C======================== S U B R O U T I N E S ======================== C C ZZDATE USER ROUTINE TO GET A DATE. C ZZTIME USER ROUTINE TO GET THE TIME OF DAY. C ZZLENG USER ROUTINE TO GET STRING LENGTH. C ZZSHFT USER ROUTINE TO SHIFT A STRING. C C MIN, INT, LEN, MOD, REAL ...INTRINSIC C C========================= P A R A M E T E R S ========================= * INTEGER PTHOUR, PTMIN, PTAMPM PARAMETER ( PTHOUR = 1, PTMIN = 4, PTAMPM = 7 ) * INTEGER PTMON, PTDAY, PTYEAR, PTDAYN PARAMETER ( PTMON = 24, PTDAY = 34, PTYEAR = 40, PTDAYN = 13 ) * CHARACTER*(*) MODEL PARAMETER ( MODEL ='00:00 A.M., WEDNESDAY, SEPTEMBER 00, 1999') * C================= L O C A L D E C L A R A T I O N S ================= * INTEGER KMON, TO, K, DAYNO, MODLEN * INTEGER ZZLENG * CHARACTER *10 TEMP CHARACTER *41 TDATE CHARACTER * 9 MONTHS(12), DAYS(0:6) * C=============================== S A V E =============================== * SAVE MONTHS, DAYS * C======================= E Q U I V A L E N C E S ======================= C C THERE ARE NO EQUIVALENCES. C C============================= C O M M O N ============================= C C THERE ARE NO COMMON BLOCKS. C C=============================== D A T A =============================== * DATA MONTHS( 1), MONTHS( 2)/'JANUARY ','FEBRUARY '/ DATA MONTHS( 3), MONTHS( 4)/'MARCH ','APRIL '/ DATA MONTHS( 5), MONTHS( 6)/'MAY ','JUNE '/ DATA MONTHS( 7), MONTHS( 8)/'JULY ','AUGUST '/ DATA MONTHS( 9), MONTHS(10)/'SEPTEMBER','OCTOBER '/ DATA MONTHS(11), MONTHS(12)/'NOVEMBER ','DECEMBER '/ * DATA DAYS(0) / 'SUNDAY ' / DATA DAYS(1) / 'MONDAY ' / DATA DAYS(2) / 'TUESDAY ' / DATA DAYS(3) / 'WEDNESDAY' / DATA DAYS(4) / 'THURSDAY ' / DATA DAYS(5) / 'FRIDAY ' / DATA DAYS(6) / 'SATURDAY ' / * C========================== E X E C U T I O N ========================== * TDATE = MODEL MODLEN = LEN(TDATE) * CALL ZZDATE(TEMP) * IF ( TEMP(8:8) .EQ. '0' ) THEN TEMP(8:8) = ' ' ENDIF * TDATE ( PTDAY : PTDAY+1 ) = TEMP(8:9) TDATE ( PTYEAR : PTYEAR+1 ) = TEMP(2:3) * READ ( TEMP(8:9), '(I2)' ) DAYNO * READ ( TEMP(2:3), '(I2)' ) K * K = K + 1900 * READ ( TEMP(5:6), '(I2)' ) KMON * TDATE(PTMON:PTMON+8) = MONTHS(KMON) * TO = ZZLENG ( MONTHS(KMON) ) * IF ( TO .NE. 9 ) THEN * CALL ZZSHFT ( TDATE, PTMON+9, PTMON+TO, MODLEN ) * ENDIF * IF ( KMON .EQ. 1 .OR. KMON .EQ. 2 ) THEN * KMON = KMON + 13 K = K - 1 * ELSE * KMON = KMON + 1 * ENDIF * DAYNO = DAYNO + INT ( REAL(KMON) * 30.6001 ) DAYNO = DAYNO + INT ( REAL( K ) * 365.25 ) * DAYNO = MOD ( DAYNO+5, 7 ) * CALL ZZTIME(TEMP) * TDATE(PTMIN:PTMIN+1) = TEMP(5:6) * READ ( TEMP(2:3), '(I2)' ) K * IF ( K .GE. 13 ) THEN * K = K-12 * TDATE(PTAMPM:PTAMPM) = 'P' * ELSE IF ( K .EQ. 12 ) THEN * TDATE(PTAMPM:PTAMPM) = 'P' * ELSE IF ( K .EQ. 0 ) THEN * K = K + 12 * TDATE(PTAMPM:PTAMPM) = 'A' * ELSE * TDATE(PTAMPM:PTAMPM) = 'A' * ENDIF * WRITE ( TDATE(PTHOUR:PTHOUR+1), '(I2)' ) K * TDATE(PTDAYN:PTDAYN+8) = DAYS(DAYNO) * K = ZZLENG (DAYS(DAYNO)) * IF ( K .NE. 9 ) THEN C ==> SHIFT OVER BLANKS. * CALL ZZSHFT ( TDATE, PTDAYN+9, PTDAYN+K, MODLEN ) * ENDIF * GOTO 90000 * C=============================== E X I T =============================== * 90000 MODLEN = MIN ( MODLEN, LEN(CHDATE) ) * CHDATE(1:MODLEN) = TDATE * RETURN * C============================ F O R M A T S ============================ C C THERE ARE NO FORMATS USED. C C================================ E N D ================================ * END SUBROUTINE ZZFNS ( IFG, N, X, F, G, IW, DUMMY, WORK ) C!!!! SUBROUTINE ZZFNS ( IFG, N, X, F, G, IW, WORK, DUMMY ) * C============== A R G U M E N T D E C L A R A T I O N S ============== * INTEGER N, IFG, IW(*) * DOUBLE PRECISION F, X(N), G(N), WORK(*) C!!!! REAL F, X(N), G(N), WORK(*) * C *** NOTE THAT THESE ARE **DELIBERATELY** OPPOSITE TO OTHER PAIRS. REAL DUMMY(*) C!!!! DOUBLE PRECISION DUMMY(*) * C============================= S T A T U S ============================= C C SYSTEM DEPENDENCE: NONE C C CONVERSION FOR SINGLE OR DOUBLE PRECISION: REQUIRED (SEE CONVRT) C C IGNORE LINES BEGINNING WITH "C!!!!" . C C THIS VERSION IS IN D O U B L E PRECISION. C!!!! THIS VERSION IS IN S I N G L E PRECISION. C C REVISION STATUS: DATE AUTHOR VERSION C C JUN. 1, 1985 A. BUCKLEY 1.0 C C======================== D E S C R I P T I O N ======================== C C THIS TEST FUNCTION EVALUATES ONE OF THE STANDARD TEST C FUNCTIONS PROVIDED WITH TESTPACK. THE ARGUMENTS IN THE CALLING C SEQUENCE HAVE PRECISELY THE SAME MEANING AS IN THE ROUTINE ZZEVAL. C C THE TEST FUNCTION TO USE IS SELECTED BY CALLING THE ENTRY C POINT ZZFSET ( FUNCNO ). THE VALUE OF THE INTEGER, FUNCNO, C SPECIFIES WHICH OF THE TEST FUNCTIONS IS TO BE USED; THE FUNCTION C IS CHOSEN USING A COMPUTED GOTO. C C SOME OF THE FUNCTIONS NEED SPECIAL ARGUMENTS (OTHER THAN THE C VALUE OF X); THESE ARE PROVIDED THROUGH THE ENTRY POINT ZZFPAR. A C MAXIMUM OF FIVE ARGUMENTS ARE PROVIDED. IF THE MAXIMUM NUMBER OF C ARGUMENTS IS TO BE INCREASED, THE PARAMETER FNO SHOULD BE C INCREASED. IT MUST AGREE WITH THE VALUE USED IN ZZTP. C C ALL FUNCTION ARGUMENTS ARE REAL. INTEGER VALUES MAY BE PASSED C BY ASSIGNING THE INTEGER VALUE TO A REAL ARGUMENT AND THEN USING C NINT TO RECOVER THE INTEGER VALUE. C C THE AMOUNT OF SPACE AVAILABLE IN THE ARRAY WORK IS DEFINED C BY CALLING THE ENTRY POINT ZZFSET. THIS MEANS THAT IT DOES NOT C HAVE TO BE PROVIDED IN THE CALL TO ZZFNS OR IN THE CALL TO ZZEVAL. C IT IS ALSO EASIER SINCE IT SELDOM CHANGES. C C======================= E N T R Y P O I N T S ======================= C C ZZFNS THE NATURAL ENTRY POINT. C ZZFSET THE ENTRY POINT TO SELECT A PARTICULAR FUNCTION. C IT ALSO SETS THE SIZE OF WORKING STORAGE AVAILABLE. C ZZFPAR AN ENTRY TO DEFINE ARGUMENTS NEEDED BY TEST FUNCTIONS. C C======================== S U B R O U T I N E S ======================== C C PREDEFINED FUNCTIONS : SIN, COS, TAN, ACOS, ATAN, ABS, MAX, NINT C EXP, LOG, MIN, MOD, SIGN, SQRT, REAL(DBLE) C C STATEMENT FUNCTION: RD C C========================= P A R A M E T E R S ========================= * DOUBLE PRECISION ZERO, ONE, TWO, THREE C!!!! REAL ZERO, ONE, TWO, THREE PARAMETER ( ZERO = 0D0, ONE = 1D0, TWO = 2D0, THREE = 3D0) * DOUBLE PRECISION FOUR, FIVE, SIX, SEVEN C!!!! REAL FOUR, FIVE, SIX, SEVEN PARAMETER ( FOUR = 4D0, FIVE = 5D0,SIX = 6D0, SEVEN = 7D0) * DOUBLE PRECISION EIGHT, NINE, TEN C!!!! REAL EIGHT, NINE, TEN PARAMETER ( EIGHT = 8D0, NINE = 8D0, TEN = 10D0 ) * DOUBLE PRECISION TENTH, FIFTH, HALF C!!!! REAL TENTH, FIFTH, HALF PARAMETER ( TENTH = .1D0, FIFTH = .2D0, HALF = .5D0 ) * DOUBLE PRECISION RPT9, RPT8, RD29 C!!!! REAL RPT9, RPT8, RD29 PARAMETER ( RPT9 = .9D0, RPT8 = .8D0, RD29 = 1D0/29D0 ) DOUBLE PRECISION R11, R12, R13, R14 C!!!! REAL R11, R12, R13, R14 PARAMETER ( R11 = 11D0, R12 = 12D0, R13 = 13D0,R14 = 14D0) * DOUBLE PRECISION R15, R16, R17, R18 C!!!! REAL R15, R16, R17, R18 PARAMETER ( R15 = 15D0, R16 = 16D0, R17 = 17D0,R18 = 18D0) * DOUBLE PRECISION R19, R20, R25, R29 C!!!! REAL R19, R20, R25, R29 PARAMETER ( R19 = 19D0, R20 = 20D0, R25 = 25D0,R29 = 29D0) * DOUBLE PRECISION R32, R36, R40, R42 C!!!! REAL R32, R36, R40, R42 PARAMETER ( R32 = 32D0, R36 = 36D0, R40 = 40D0,R42 = 42D0) * DOUBLE PRECISION R45, R49 C!!!! REAL R45, R49 PARAMETER ( R45 = 45D0, R49 = 49D0 ) * DOUBLE PRECISION R50, R56, R84, R90 C!!!! REAL R50, R56, R84, R90 PARAMETER ( R50 = 50D0, R56 = 56D0, R84 = 84D0,R90 = 90D0) * DOUBLE PRECISION R100, R180, R200 C!!!! REAL R100, R180, R200 PARAMETER ( R100 = 100D0, R180 = 180D0, R200 = 200D0 ) * DOUBLE PRECISION R256, R360, R400 C!!!! REAL R256, R360, R400 PARAMETER ( R256 = 256D0, R360 = 360D0, R400 = 400D0 ) * DOUBLE PRECISION R600, R681, R991 C!!!! REAL R600, R681, R991 PARAMETER ( R600 = 600D0, R681 = 681D0, R991 = 991D0 ) * DOUBLE PRECISION R1162, R2324 C!!!! REAL R1162, R2324 PARAMETER ( R1162 = 1162D0, R2324 = 2324D0 ) * DOUBLE PRECISION R10000, R40000 C!!!! REAL R10000, R40000 PARAMETER ( R10000 = 10000D0, R40000 = 40000D0 ) DOUBLE PRECISION R1PD6, R2PDM6 C!!!! REAL R1PD6, R2PDM6 PARAMETER ( R1PD6 = 1D6, R2PDM6 = 2D-6 ) * DOUBLE PRECISION RP04, RP01, R1PZ1 C!!!! REAL RP04, RP01, R1PZ1 PARAMETER ( RP04 = 4D-2, RP01 = .01D0, R1PZ1 = 1.0001D0 ) * DOUBLE PRECISION R1P2 C!!!! REAL R1P2 PARAMETER ( R1P2 = 1.2D0 ) * DOUBLE PRECISION R1P5, R2P5, R2P625 C!!!! REAL R1P5, R2P5, R2P625 PARAMETER ( R1P5 = 1.5D0, R2P5 = 2.5D0, R2P625 = 2.625D0 ) * DOUBLE PRECISION R10P1, R19P8, R20P2 C!!!! REAL R10P1, R19P8, R20P2 PARAMETER ( R10P1 = 10.1D0,R19P8 = 19.8D0,R20P2 = 20.2D0 ) * DOUBLE PRECISION R2D3, R4D3, R7D3 C!!!! REAL R2D3, R4D3, R7D3 PARAMETER ( R2D3 = 2D0/3D0,R4D3 = 4D0/3D0,R7D3 = 7D0/3D0 ) * DOUBLE PRECISION R2P25 C!!!! REAL R2P25 PARAMETER ( R2P25 = 2.25D0 ) * INTEGER ALPHA, BETA, GAMMA PARAMETER ( ALPHA = 5, BETA = 14, GAMMA = 3 ) * * INTEGER JUSTF, BOTH, JUSTG, NOOP PARAMETER ( JUSTF = 1, BOTH = 0, JUSTG = -1, NOOP = 2 ) * INTEGER FNO PARAMETER ( FNO = 10 ) C THE RETURN CODES TO BE USED BY THE FUNCTION EVALUATION ROUTINE C TO INDICATE TO THE MINIMIZATION ROUTINE WHETHER OR NOT THE CALL C WAS SUCCESSFUL. * INTEGER COK, CABORT, CLIMIT PARAMETER ( COK = 0, CABORT = -1, CLIMIT = -2 ) * INTEGER CNOF, CNOG, CNOFG PARAMETER ( CNOF = -3, CNOG = -4, CNOFG = -5 ) * C================= L O C A L D E C L A R A T I O N S ================= * INTEGER OK, ABORT, LIMIT, NOF, NOG, NOFG INTEGER SOK, SABORT, SLIMIT, SNOF, SNOG, SNOFG * INTEGER FUNCNO, SIZE, SDOF, SDOFG, SDOG, SNOOP INTEGER SETFNC, SETSIZ, F1 * INTEGER I, J, K, JLO, JHI, IEG, IY INTEGER I1, I2, IC, IS, IV, IA, IB, M, L INTEGER NPTS, NEQ, IP, IQ, NSYS, IMESH, IH, TT INTEGER IB0, IB1, IDB0, IDB1, ILAST, BACK, RET, P * INTEGER ALPHA2, NOVER2, DOF, DOG, DOFG, NEITHR * LOGICAL EVEN, GFIRST, FIRST, ERROR, FONLY, GONLY LOGICAL DONE, PROB13 * DOUBLE PRECISION FARG ( FNO ), SARG ( FNO ), ZZMPAR, HUGE C!!!! REAL FARG ( FNO ), SARG ( FNO ), ZZMPAR, HUGE * * C--------- VARIABLES FOR THE TEST FUNCTIONS. * * DOUBLE PRECISION X1, X2, X3, X4, X5, X6 C!!!! REAL X1, X2, X3, X4, X5, X6 * DOUBLE PRECISION X7, X8, X9, X10, X11 C!!!! REAL X7, X8, X9, X10, X11 * DOUBLE PRECISION G1, G2, G3, G4, G5, G6 C!!!! REAL G1, G2, G3, G4, G5, G6 * DOUBLE PRECISION G7, G8, G9, G10, G11 C!!!! REAL G7, G8, G9, G10, G11 * DOUBLE PRECISION W1, W2, W3, W4, W5, W6 C!!!! REAL W1, W2, W3, W4, W5, W6 * DOUBLE PRECISION W7, W8, W9, W10, W11, W12 C!!!! REAL W7, W8, W9, W10, W11, W12 * DOUBLE PRECISION R, S, T, R1, BIGGST, SMLLST C!!!! REAL R, S, T, R1, BIGGST, SMLLST * DOUBLE PRECISION R2, R3, RI, RK , SK, TI C!!!! REAL R2, R3, RI, RK , SK, TI * DOUBLE PRECISION XI, XK, YI, PI, U, SUM C!!!! REAL XI, XK, YI, PI, U, SUM * DOUBLE PRECISION XP1, XM1, R2P, RD, TPI, TPIS C!!!! REAL XP1, XM1, R2P, RD, TPI, TPIS * DOUBLE PRECISION HJ, HJJ, DELTX, SINX C!!!! REAL HJ, HJJ, DELTX, SINX * DOUBLE PRECISION U1, U2, RF1, RF2, RF3, RF4, DH, DHH C!!!! REAL U1, U2, RF1, RF2, RF3, RF4, DH, DHH * * C--------- DATA ARRAYS FOR TESTPACK FUNCTIONS * * DOUBLE PRECISION AL (50), ARGASY (15), BARD7Y (15) C!!!! REAL AL (50), ARGASY (15), BARD7Y (15) * DOUBLE PRECISION HIM32A (7), HIM32B (7) C!!!! REAL HIM32A (7), HIM32B (7) * DOUBLE PRECISION KOWOSU (11), KOWOSY (11) C!!!! REAL KOWOSU (11), KOWOSY (11) * DOUBLE PRECISION ORBETA (33), OD (33), MEY (16) C!!!! REAL ORBETA (33), OD (33), MEY (16) * DOUBLE PRECISION OSB1Y (33), OSB2Y (65) C!!!! REAL OSB1Y (33), OSB2Y (65) * INTEGER A (50), B (56) * C=============================== S A V E =============================== * SAVE FUNCNO, FARG, SIZE, GFIRST, FIRST, PI, R2P, BIGGST, SMLLST SAVE NPTS, NEQ, IP, IQ, NSYS, IA, IB, DH, DHH, HUGE SAVE IMESH, IH, IB0, IB1, IDB0, IDB1, PROB13, DONE, M SAVE DOF, DOG, DOFG, NEITHR SAVE OK, ABORT, LIMIT, NOF, NOG, NOFG * C--------- SAVE DATA ARRAYS FOR THE TEST FUNCTIONS. * * SAVE ARGASY, BARD7Y, AL, HIM32A, HIM32B, KOWOSU, KOWOSY SAVE ORBETA, MEY , OD, OSB1Y , OSB2Y , A , B * C======================= E Q U I V A L E N C E S ======================= C C THERE ARE NO EQUIVALENCES. C C============================= C O M M O N ============================= C C THERE ARE NO COMMON BLOCKS. C C=============================== D A T A =============================== * DATA FUNCNO /1/, SIZE /1/ * DATA FARG / FNO * 1.D0 / * DATA FIRST/.TRUE./, GFIRST/.TRUE./ * DATA DOF/JUSTF/, DOFG/ BOTH/, DOG/JUSTG/, NEITHR/NOOP/ DATA OK/ COK/, ABORT/CABORT/, LIMIT/CLIMIT/ DATA NOF/ CNOF/, NOG/ CNOG/, NOFG/ CNOFG/ * C--------- DATA FOR TESTPACK FUNCTION ARGAUS * DATA ARGASY( 1), ARGASY( 2), ARGASY( 3), ARGASY( 4), ARGASY( 5) - / 9.000 D-4, 4.400 D-3, 1.750 D-2, 5.400 D-2, 1.295 D-1 / * DATA ARGASY( 6), ARGASY( 7), ARGASY( 8), ARGASY( 9), ARGASY(10) - / 2.420 D-1, 3.521 D-1, 3.989 D-1, 3.521 D-1, 2.420 D-1 / * DATA ARGASY(11), ARGASY(12), ARGASY(13), ARGASY(14), ARGASY(15) - / 1.295 D-1, 5.400 D-2, 1.750 D-2, 4.400 D-3, 9.000 D-4 / * * C--------- DATA FOR TESTPACK FUNCTION BARD70 * DATA BARD7Y( 1), BARD7Y( 2), BARD7Y( 3), BARD7Y( 4), BARD7Y( 5) - / .14 D0 , .18 D0 , .22 D0 , .25 D0 , .29 D0 / * DATA BARD7Y( 6), BARD7Y( 7), BARD7Y( 8), BARD7Y( 9), BARD7Y(10) - / .32 D0 , .35 D0 , .39 D0 , .37 D0 , .58 D0 / * DATA BARD7Y(11), BARD7Y(12), BARD7Y(13), BARD7Y(14), BARD7Y(15) - / .73 D0 , .96 D0 , 1.34 D0 , 2.10 D0 , 4.39 D0 / * C--------- DATA FOR TESTPACK FUNCTION CHNRSN * DATA AL(1), AL(2), AL(3), AL(4), AL(5), AL(6), AL(7), AL(8) - / 1.25D0,1.40D0,2.40D0,1.40D0,1.75D0,1.20D0,2.25D0,1.20D0/ * DATA AL(9) ,AL(10),AL(11),AL(12),AL(13),AL(14),AL(15),AL(16) - / 1.00D0,1.10D0,1.50D0,1.60D0,1.25D0,1.25D0,1.20D0,1.20D0/ * DATA AL(17),AL(18),AL(19),AL(20),AL(21),AL(22),AL(23),AL(24) - / 1.40D0,0.50D0,0.50D0,1.25D0,1.80D0,0.75D0,1.25D0,1.40D0/ * DATA AL(25),AL(26),AL(27),AL(28),AL(29),AL(30) - / 1.60D0,2.00D0,1.00D0,1.60D0,1.25D0,2.75D0/ * DATA AL(31),AL(32),AL(33),AL(34),AL(35),AL(36),AL(37),AL(38) - / 1.25D0,1.25D0,1.25D0,3.00D0,1.50D0,2.00D0,1.25D0,1.40D0/ * DATA AL(39),AL(40),AL(41),AL(42),AL(43),AL(44),AL(45),AL(46) - / 1.80D0,1.50D0,2.20D0,1.40D0,1.50D0,1.25D0,2.00D0,1.50D0/ * DATA AL(47),AL(48),AL(49),AL(50) - / 1.25D0,1.40D0,0.60D0,1.50D0/ * C--------- DATA FOR TESTPACK FUNCTION HIMM32 * DATA HIM32A(1), HIM32A(2), HIM32A(3), HIM32A(4) - / 0.0D0, 4.28D-4, 1.0D-3, 1.61D-3 / * DATA HIM32A(5), HIM32A(6), HIM32A(7) - / 2.09D-3, 3.48D-3, 5.25D-3 / * DATA HIM32B(1), HIM32B(2), HIM32B(3), HIM32B(4) - / 7.391D0, 1.118D1, 1.644D1, 1.62D1 / * DATA HIM32B(5), HIM32B(6), HIM32B(7) - / 2.22D1, 2.402D1, 3.132D1 / * C--------- DATA FOR TESTPACK FUNCTION KOWOSB * DATA KOWOSU(1), KOWOSU(2), KOWOSU(3), KOWOSU(4) - / 4.0D0, 2.0D0, 1.0D0, 0.5D0 / * DATA KOWOSU(5), KOWOSU(6), KOWOSU(7), KOWOSU(8) - / 0.25D0, 0.167D0, 0.125D0, 0.1D0 / * DATA KOWOSU(9), KOWOSU(10), KOWOSU(11) - / 0.0833D0, 0.0714D0, 0.0625D0 / * DATA KOWOSY(1), KOWOSY(2), KOWOSY(3), KOWOSY(4) - / 0.1957D0, 0.1947D0, 0.1735D0, 0.1600D0 / * DATA KOWOSY(5), KOWOSY(6), KOWOSY(7), KOWOSY(8) - / 0.0844D0, 0.0627D0, 0.0456D0, 0.0342D0 / * DATA KOWOSY(9), KOWOSY(10), KOWOSY(11) - / 0.0323D0, 0.0235D0, 0.0246D0 / * C--------- DATA FOR TESTPACK FUNCTION MEYER * DATA MEY(1), MEY(2), MEY(3), MEY(4), MEY(5), MEY(6) - /3.478D4, 2.861D4, 2.365D4, 1.963D4, 1.637D4, 1.372D4/ * DATA MEY(7), MEY(8), MEY(9), MEY(10), MEY(11), MEY(12) - /1.154D4, 9.744D3, 8.261D3, 7.030D3, 6.005D3, 5.147D3/ * DATA MEY(13), MEY(14), MEY(15), MEY(16) - /4.427D3, 3.820D3, 3.307D3, 2.872D3/ * C--------- DATA FOR TESTPACK FUNCTION ORTOIT * DATA ORBETA(1),ORBETA(2),ORBETA(3),ORBETA(4),ORBETA(5) - /1.0D0, 1.5D0, 1.0D0, 0.1D0, 1.5D0/ * DATA ORBETA(6),ORBETA(7),ORBETA(8),ORBETA(9),ORBETA(10) - /2.0D0, 1.0D0, 1.5D0, 3.0D0, 2.0D0/ * DATA ORBETA(11),ORBETA(12),ORBETA(13),ORBETA(14),ORBETA(15) - /1.0D0, 3.0D0, 0.1D0, 1.5D0, 0.15D0/ * DATA ORBETA(16),ORBETA(17),ORBETA(18),ORBETA(19),ORBETA(20) - /2.0D0, 1.0D0, 0.1D0, 3.0D0, 0.1D0/ * DATA ORBETA(21),ORBETA(22),ORBETA(23),ORBETA(24),ORBETA(25) - /1.2D0, 1.0D0, 0.1D0, 2.0D0, 1.2D0/ * DATA ORBETA(26),ORBETA(27),ORBETA(28),ORBETA(29),ORBETA(30) - /3.0D0, 1.5D0, 3.0D0, 2.0D0, 1.0D0/ * DATA ORBETA(31),ORBETA(32),ORBETA(33) - /1.2D0, 2.0D0, 1.0D0/ * DATA OD(1), OD(2), OD(3), OD(4), OD(5), OD(6) - / 5.0D0,5.0D0,5.0D0,2.5D0,6.0D0,6.0D0 / * DATA OD(7), OD(8), OD(9), OD(10), OD(11), OD(12) - / 5.0D0,6.0D0,10.0D0,6.0D0,5.0D0,9.0D0 / * DATA OD(13), OD(14), OD(15), OD(16), OD(17), OD(18) - / 2.0D0,7.0D0,2.5D0,6.0D0,5.0D0,2.0D0 / * DATA OD(19), OD(20), OD(21), OD(22), OD(23), OD(24) - / 9.0D0,2.0D0,5.0D0,5.0D0,2.5D0,5.0D0 / * DATA OD(25), OD(26), OD(27), OD(28), OD(29), OD(30) - / 6.0D0,10.0D0,7.0D0,10.0D0,6.0D0,5.0D0 / * DATA OD(31), OD(32), OD(33) - / 4.0D0,4.0D0,4.0D0 / * DATA A(1),A(2),A(3),A(4),A(5),A(6),A(7),A(8),A(9) - / -31,-1,-2,-4,-6,-8,-10,-12,+11 / DATA A(10),A(11),A(12),A(13),A(14),A(15),A(16),A(17),A(18) - / +13,-14,-16,+9,-18,+5,+20,-21,-19 / DATA A(19),A(20),A(21),A(22),A(23),A(24),A(25),A(26),A(27) - / -23,+7,-25,-28,-29,-32,+3,-33,-35 / DATA A(28),A(29),A(30),A(31),A(32),A(33),A(34),A(35),A(36) - / -36,+30,-37,+38,-39,-40,-41,-44,-46 / DATA A(37),A(38),A(39),A(40),A(41),A(42),A(43),A(44),A(45) - / +42,+45,+48,-50,+26,+34,-43,+15,+17 / DATA A(46),A(47),A(48),A(49),A(50) - / +24,-47,-49,-22,-27 / * DATA B(1),B(2),B(3),B(4),B(5),B(6),B(7),B(8),B(9) - / -1,+2,-3,+4,-5,+6,-7,+8,-9 / DATA B(10),B(11),B(12),B(13),B(14),B(15),B(16),B(17),B(18) - / +10,-11,+12,-13,+14,-15,+16,-17,+18 / DATA B(19),B(20),B(21),B(22),B(23),B(24),B(25),B(26),B(27) - / -19,-20,0,+22,+23,-24,+25,-26,+27 / DATA B(28),B(29),B(30),B(31),B(32),B(33),B(34),B(35),B(36) - / -28,+29,-30,+31,-32,+33,-34,-35,+21 / DATA B(37),B(38),B(39),B(40),B(41),B(42),B(43),B(44),B(45),B(46) - / -36,+37,-38,-39,-40,+41,-42,+43,+44,-50 / DATA B(47),B(48),B(49),B(50),B(51),B(52),B(53),B(54),B(55),B(56) - / +45,+46,-47,-48,-49,0,0,0,0,0 / * C--------- DATA FOR TESTPACK FUNCTION OSBRN1 * DATA OSB1Y(1), OSB1Y(2), OSB1Y(3), OSB1Y(4), OSB1Y(5) -/.844D0, .908D0, .932D0, .936D0, .925D0/ * DATA OSB1Y(6), OSB1Y(7), OSB1Y(8), OSB1Y(9), OSB1Y(10) -/.908D0, .881D0, .850D0, .818D0, .784D0/ * DATA OSB1Y(11), OSB1Y(12), OSB1Y(13), OSB1Y(14), OSB1Y(15) -/.751D0, .718D0, .685D0, .658D0, .628D0/ * DATA OSB1Y(16), OSB1Y(17), OSB1Y(18), OSB1Y(19), OSB1Y(20) -/.603D0, .580D0, .558D0, .538D0, .522D0/ * DATA OSB1Y(21), OSB1Y(22), OSB1Y(23), OSB1Y(24), OSB1Y(25) -/.506D0, .490D0, .478D0, .467D0, .457D0/ * DATA OSB1Y(26), OSB1Y(27), OSB1Y(28), OSB1Y(29), OSB1Y(30) -/.448D0, .438D0, .431D0, .424D0, .420D0/ * DATA OSB1Y(31), OSB1Y(32), OSB1Y(33) -/.414D0, .411D0, .406D0/ * C--------- DATA FOR TESTPACK FUNCTION OSBRN2 * DATA OSB2Y(1), OSB2Y(2), OSB2Y(3), OSB2Y(4), OSB2Y(5) -/1.366D0, 1.191D0, 1.112D0, 1.013D0, .991D0/ * DATA OSB2Y(6), OSB2Y(7), OSB2Y(8), OSB2Y(9), OSB2Y(10) -/.885D0, .831D0, .847D0, .786D0, .725D0/ * DATA OSB2Y(11), OSB2Y(12), OSB2Y(13), OSB2Y(14), OSB2Y(15) -/.746D0, .679D0, .608D0, .655D0, .616D0/ * DATA OSB2Y(16), OSB2Y(17), OSB2Y(18), OSB2Y(19), OSB2Y(20) -/.606D0, .602D0, .626D0, .651D0, .724D0/ * DATA OSB2Y(21), OSB2Y(22), OSB2Y(23), OSB2Y(24), OSB2Y(25) -/.649D0, .649D0, .694D0, .644D0, .624D0/ * DATA OSB2Y(26), OSB2Y(27), OSB2Y(28), OSB2Y(29), OSB2Y(30) -/.661D0, .612D0, .558D0, .533D0, .495D0/ * DATA OSB2Y(31), OSB2Y(32), OSB2Y(33), OSB2Y(34), OSB2Y(35) -/.50D0, .423D0, .395D0, .375D0, .372D0/ * DATA OSB2Y(36), OSB2Y(37), OSB2Y(38), OSB2Y(39), OSB2Y(40) -/.391D0, .396D0, .405D0, .428D0, .429D0/ * DATA OSB2Y(41), OSB2Y(42), OSB2Y(43), OSB2Y(44), OSB2Y(45) -/.523D0, .562D0, .607D0, .653D0, .672D0/ * DATA OSB2Y(46), OSB2Y(47), OSB2Y(48), OSB2Y(49), OSB2Y(50) -/.708D0, .633D0, .668D0, .645D0, .632D0/ * DATA OSB2Y(51), OSB2Y(52), OSB2Y(53), OSB2Y(54), OSB2Y(55) -/.591D0, .559D0, .597D0, .625D0, .739D0/ * DATA OSB2Y(56), OSB2Y(57), OSB2Y(58), OSB2Y(59), OSB2Y(60) -/.710D0, .729D0, .720D0, .636D0, .581D0/ * DATA OSB2Y(61), OSB2Y(62), OSB2Y(63), OSB2Y(64), OSB2Y(65) -/.428D0, .292D0, .162D0, .098D0, .054D0/ * * C========================== E X E C U T I O N ========================== * C--------- FUNCTION DEFINITION * RD (I) = DBLE (I) C!!!! RD (I) = REAL (I) * C--------- SOME ONE TIME ONLY CONSTANTS. * IF ( GFIRST ) THEN PI = ACOS(-ONE) TPI = TWO * PI TPIS = TPI * PI R2P = ONE / TPI HUGE = ZZMPAR(3)/TEN SMLLST = LOG(ZZMPAR(2)*TEN) BIGGST = LOG(HUGE) ENDIF * C--------- SET LOGICAL FLAGS AND SELECT FUNCTION. * FONLY = IFG .EQ. DOF GONLY = IFG .EQ. DOG RET = OK * GOTO( 51, 99, 50, 99, 99, 99, 99, 99, 99, - 99, 8, 99, 99, 99, 99, 10, 99, 99, 99, - 99, 99, 99, 99, 52, 22, 5, 20, 99, 99, - 99, 99, 99, 99, 99, 99, 99, 99, 99, 99, - 44, 45, 99, 99, 99, 99, 99, 99, 99, 99, - 99, 99, 99, 99, 99, 99, 99, 99, 99, 99, - 99, 99, 99, 99, 99, 99, 99, 99 ) - FUNCNO * C THESE DUMMY GO TO'S JUST MAKE POSSIBLE FUTURE CHANGES C MORE CONVENIENT. THE TEST FUNCTIONS APPEAR HERE IN C ALPHABETICAL ORDER. * C BARD70 5 GOTO 3600 * C BIGGS6 8 GOTO 2100 * C BOX663 10 GOTO 2600 * C CRGLVY 20 GOTO 3700 * C ENGVL2 22 GOTO 3500 * C PENAL1 44 GOTO 5000 * C PENAL2 45 GOTO 5100 * C PWSING 50 GOTO 1300 * C ROSENB 51 GOTO 1100 * C SCHMVT 52 GOTO 3400 * 99 GOTO 10000 * * C>>>>> NOTE : IF WE SUPPOSE THAT EACH OF THESE TEST FUNCTIONS HAD C>>>>> BEEN CODED AS A SEPARATE ROUTINE, THEN, UNLESS C>>>>> OTHERWISE SPECIFIED, ALL TEST FUNCTIONS WOULD HAVE C>>>>> HAD AN ARGUMENT LIST AS FOLLOWS: C>>>>> C>>>>> ( CASE, N, X, F, G ) C>>>>> C>>>>> THOSE WHICH WOULD REQUIRE ADDITIONAL ARGUMENTS ARE C>>>>> NOTED BY GIVING A SUITABLE CALLING SEQUENCE. THIS C>>>>> SERVES TO DEFINE THE SPECIAL ARGUMENTS FOR THOSE TEST C>>>>> FUNCTIONS. SEE FOR EXAMPLE PENAL2 AT 5100. * * C--------- TESTPACK FUNCTION ROSENB * * 1100 X1 = X(1) W1 = ONE - X1 W2 = X(2) - X1*X1 * IF ( .NOT. GONLY ) THEN F = R100*W2*W2 + W1*W1 ENDIF * IF ( .NOT. FONLY ) THEN G(1) = -R400*W2*X1 - TWO*W1 G(2) = R200*W2 ENDIF * GOTO 10000 * * * C--------- TESTPACK FUNCTION PWSING * * 1300 IF ( .NOT. GONLY ) THEN F = ZERO ENDIF * IF ( 4 * (N/4) .NE. N ) THEN * IF ( .NOT. FONLY ) THEN DO 1310 I = 1,N G(I) = ZERO 1310 CONTINUE ENDIF * ELSE * DO 1320 I=1,N/4 * J = 4*I * W1 = X(J-3) W2 = X(J-2) W3 = X(J-1) W4 = X(J ) * W5 = W1 + TEN * W2 W6 = W3 - W4 W2 = W2 - TWO * W3 W3 = W2 * W2 *W2 W1 = W1 - W4 W4 = W1 * W1 * W1 * IF ( .NOT. GONLY ) THEN F = F + W5*W5 + FIVE*W6*W6 + W2*W3 + TEN*W1*W4 ENDIF * IF ( .NOT. FONLY ) THEN G(J-3) = TWO * W5 + R40 * W4 G(J-2) = R20 * W5 + FOUR * W3 G(J-1) = TEN * W6 - EIGHT * W3 G(J ) = -TEN * W6 - R40 * W4 ENDIF * 1320 CONTINUE * ENDIF * GOTO 10000 * * * C--------- TESTPACK FUNCTION BIGGS6 ( N, X, F, G, IFG, NINT(FARG(1))) C--------- NINT(FARG(1)) IS M * * 2100 X1 = X(1) X2 = X(2) X3 = X(3) X4 = X(4) X5 = X(5) X6 = X(6) * IF ( .NOT. GONLY ) THEN F = ZERO ENDIF * IF ( .NOT. FONLY ) THEN G1 = ZERO G2 = ZERO G3 = ZERO G4 = ZERO G5 = ZERO G6 = ZERO ENDIF * DO 2110 I = 1, NINT(FARG(1)) T = RD(I) TI = T/TEN IF ( MAX(-T,-TI*FOUR,-TI*X1,-TI*X2,-TI*X5) .LE. BIGGST ) THEN YI = EXP(-TI) - FIVE * EXP(-T) + THREE*EXP(-FOUR*TI) W3 = EXP(-TI*X1) W4 = EXP(-TI*X2) W5 = EXP(-TI*X5) ELSE RET = NOFG GOTO 90000 ENDIF RI = X3*W3 - X4*W4 + X6*W5 - YI * IF ( .NOT. GONLY ) THEN F = F + RI*RI ENDIF * IF ( .NOT. FONLY ) THEN W1 = TI*RI G1 = G1 - W3*W1 G2 = G2 + W4*W1 G3 = G3 + W3*RI G4 = G4 - W4*RI G5 = G5 - W5*W1 G6 = G6 + W5*RI ENDIF * 2110 CONTINUE * IF ( .NOT. FONLY ) THEN G(1) = TWO*X3*G1 G(2) = TWO*X4*G2 G(3) = TWO * G3 G(4) = TWO * G4 G(5) = TWO*X6*G5 G(6) = TWO * G6 ENDIF * GOTO 10000 * * * C--------- TESTPACK FUNCTION BOX663 ( N, X, F, G, IFG, NINT(FARG(1))) C--------- NINT(FARG(1)) IS M * * 2600 IF ( .NOT. GONLY ) THEN F = ZERO ENDIF * IF ( .NOT. FONLY ) THEN G1 = ZERO G2 = ZERO G3 = ZERO ENDIF * DO 2610 I = 1,NINT(FARG(1)) W2 = RD(I) TI = W2/TEN IF ( MAX(-W2,-TI,-TI*X(1),-TI*X(2)) .LE. BIGGST ) THEN W3 = EXP(-TI * X(1)) W4 = EXP(-TI * X(2)) W5 = EXP(-TI) - EXP(-W2) ELSE RET = NOFG GOTO 90000 ENDIF RI = W3 - W4 - W5*X(3) * IF ( .NOT. GONLY ) THEN IF ( ABS(RI) .LE. SQRT(HUGE-MAX(F,ZERO)) ) THEN F = F + RI*RI ELSE RET = NOFG GOTO 90000 ENDIF ENDIF * IF ( .NOT. FONLY ) THEN W2 = TI*RI G1 = G1 - W3*W2 G2 = G2 + W4*W2 G3 = G3 - W5*RI ENDIF * 2610 CONTINUE * IF ( .NOT. FONLY ) THEN G(1) = TWO * G1 G(2) = TWO * G2 G(3) = TWO * G3 ENDIF * GOTO 10000 * * * * C--------- TESTPACK FUNCTION SCHMVT * * 3400 IF ( FIRST ) THEN FIRST = .FALSE. PI = ACOS(-ONE) ENDIF * X1 = X(1) X2 = X(2) X3 = X(3) * W1 = X1 - X2 W2 = X1 + X3 * W3 = ONE + W1*W1 W4 = (PI*X2 + X3) / TWO W5 = (W2/X2) - TWO IF ( -W5**2 .LE. BIGGST ) THEN W6 = EXP(-W5*W5) ELSE RET = NOFG GOTO 90000 ENDIF * IF ( .NOT. GONLY ) THEN F = - ((ONE/W3) + SIN(W4) + W6 ) ENDIF * IF ( .NOT. FONLY ) THEN * W3 = TWO*W1/(W3*W3) W4 = COS(W4)/TWO W6 = TWO*W5*W6/X2 * G(1) = W3 + W6 G(2) = -W3 - PI*W4 - W6*W2/X2 G(3) = -W4 + W6 * ENDIF * GOTO 10000 * C--------- TESTPACK FUNCTION ENGVL2 * * 3500 X1 = X(1) X2 = X(2) X3 = X(3) * W1 = X1*X1 W2 = X1*W1 W3 = X2*X2 W4 = X3*X3 * W5 = X3 - TWO W6 = FIVE*X3 - X1 + ONE W7 = W1 + W3 - ONE * W8 = W7 + W4 W9 = W7 + W5*W5 W10 = X1 + X2 + X3 - ONE W11 = X1 + X2 - X3 + ONE W12 = W2 + THREE*W3 + W6*W6 - R36 * IF ( .NOT. GONLY ) THEN F = W8*W8 + W9*W9 + W10*W10 + W11*W11 + W12*W12 ENDIF * IF ( .NOT. FONLY ) THEN W10 = W8 + W9 G(1) = TWO*(TWO*X1*W10 + TWO*(X1+X2) + W12*(THREE*W1-TWO*W6)) G(2) = TWO*(TWO*X2*W10 + TWO*(X1+X2) + SIX*W12*X2) G(3) = TWO*(TWO*(W8*X3+W5*W9) + TWO*X3 - TWO + TEN*W12*W6) ENDIF * GOTO 10000 * * C--------- TESTPACK FUNCTION BARD70 * * 3600 X1 = X(1) X2 = X(2) X3 = X(3) * IF ( .NOT. GONLY ) THEN F = ZERO ENDIF * IF ( .NOT. FONLY ) THEN G1 = ZERO G2 = ZERO G3 = ZERO ENDIF * DO 3610 I=1,15 * W1 = RD(I) W2 = RD(16-I) W3 = MIN(W1,W2) * W4 = X2*W2 + X3*W3 RI = BARD7Y(I) - (X1 + W1/W4) W4 = W4*W4 * IF ( .NOT. GONLY ) THEN F = F + RI*RI ENDIF * IF ( .NOT. FONLY ) THEN W4 = RI*W1/W4 G1 = G1 - RI G2 = G2 + W2*W4 G3 = G3 + W3*W4 ENDIF 3610 CONTINUE * IF ( .NOT. FONLY ) THEN G(1) = G1*TWO G(2) = G2*TWO G(3) = G3*TWO ENDIF * GOTO 10000 * C--------- TESTPACK FUNCTION CRGLVY * * 3700 X1 = X(1) X2 = X(2) X3 = X(3) X4 = X(4) * W1 = X2 - X3 W2 = X3 - X4 W3 = X4 - ONE * IF ( X1 .LE. BIGGST ) THEN W4 = EXP(X1) ELSE RET = NOFG GOTO 90000 ENDIF W5 = W4 - X2 W6 = TAN(W2) * IF ( .NOT. GONLY ) THEN F = W5**4 + R100*W1**6 + W6**4 + X1**8 + W3*W3 ENDIF * IF ( .NOT. FONLY ) THEN * W2 = COS(W2) W5 = FOUR * W5**3 W1 = R600 * W1**5 W6 = FOUR * W6**3 / (W2*W2) * G(1) = W4*W5 + EIGHT*X1**7 G(2) = -W5 + W1 G(3) = -W1 + W6 G(4) = -W6 + TWO*W3 ENDIF * GOTO 10000 * * * * C--------- TESTPACK FUNCTION PENAL1 ( N, X, F, G, IFG, C FARG(1), FARG(2) ) C--------- FARG(1) IS A C--------- FARG(2) IS B * 5000 RF1 = FARG ( 1 ) RF2 = FARG ( 2 ) * W1 = - ONE / FOUR W2 = ZERO * DO 5010 J = 1, N W3 = X(J) W1 = W1 + W3*W3 W3 = W3 - ONE W2 = W2 + W3*W3 5010 CONTINUE * IF ( .NOT. GONLY ) THEN F = RF1*W2 + RF2 *W1*W1 ENDIF * IF ( .NOT. FONLY ) THEN W1 = FOUR*RF2*W1 W2 = TWO*RF1 DO 5020 J = 1, N W3 = X(J) G(J) = W2 * (W3 - ONE) + W3*W1 5020 CONTINUE ENDIF * GOTO 10000 * * C--------- TESTPACK FUNCTION PENAL2 ( N, X, F, G, IFG, C FARG(1), FARG(2), WORK, SIZE) C--------- FARG(1) IS A C--------- FARG(2) IS B * 5100 RF1 = FARG ( 1 ) RF2 = FARG ( 2 ) * IF ( SIZE .LT. 2 * N ) THEN F = ZERO DO 5110 K = 1, N G(K) = ZERO 5110 CONTINUE GO TO 10000 ENDIF * W1 = EXP(TENTH) W2 = EXP(-TENTH) W3 = ZERO * I1 = 0 I2 = N * DO 5120 K = 1, N W4 = X(K) W3 = W3 + RD( N - K + 1 ) * W4 * W4 IF ( TENTH*W4 .LE. BIGGST ) THEN W5 = EXP (TENTH * W4) ELSE RET = NOFG GOTO 90000 ENDIF * IF ( K .EQ. 1 ) THEN W6 = ZERO W7 = ONE * ELSE W7 = W9 * W1 W10 = W5 + W8 - (W7 + W9) W11 = W5 - W2 * IF ( .NOT. FONLY ) THEN WORK(I1+K) = W10 WORK(I2+K) = W11 ENDIF * IF ( .NOT. GONLY ) THEN W6 = W6 + W10*W10 + W11*W11 ENDIF * ENDIF * W8 = W5 W9 = W7 * 5120 CONTINUE * W1 = X(1) - FIFTH W2 = W3 - ONE * IF ( .NOT. GONLY ) THEN F = RF1 * W6 + RF2* ( W1*W1 + W2*W2 ) ENDIF * IF ( .NOT. FONLY ) THEN W3 = FIFTH * RF1 W2 = FOUR * RF2 * W2 * DO 5130 K = 1, N * C ---NOTE THAT W8 DOES NOT NEED TO BE PRE-DEFINED WHEN K = 1. * W4 = X(K) IF ( TENTH*W4 .LE. BIGGST ) THEN W5 = EXP(TENTH * W4) ELSE RET = NOFG GOTO 90000 ENDIF W6 = W8 W7 = WORK(I2+K) * IF ( K .LT. N ) THEN W8 = WORK(I1+K+1) IF ( K .EQ. 1 ) THEN G(1) = W3 * W5 * ( W8 ) - + W2 * W4 * RD(N) + W1 * TWO * RF2 * ELSE G(K) = W3 * W5 * ( W6 + W7 + W8 ) - + W2 * W4 * RD( N - K + 1 ) * ENDIF * ELSE G(N) = W3 * W5 * ( W6 + W7 ) - + W2 * W4 * ENDIF * 5130 CONTINUE * ENDIF * GOTO 10000 * * 10000 GOTO 90000 * C>>>>>>>>>>>>>>>>>>>>>>>>>> E N T R Y ZZFFDF <<<<<<<<<<<<<<<<<<<<<<<<<< * ENTRY ZZFFDF ( SDOF, SDOG, SDOFG, SNOOP ) * DOF = SDOF DOG = SDOG DOFG = SDOFG NEITHR = SNOOP * RETURN * C>>>>>>>>>>>>>>>>>>>>>>>>>> E N T R Y ZZFRDF <<<<<<<<<<<<<<<<<<<<<<<<<< * ENTRY ZZFRDF ( SOK, SABORT, SLIMIT, SNOF, SNOG, SNOFG ) * OK = SOK ABORT = SABORT LIMIT = SLIMIT NOF = SNOF NOG = SNOG NOFG = SNOFG * RETURN * C>>>>>>>>>>>>>>>>>>>>>>>>>> E N T R Y ZZFSET <<<<<<<<<<<<<<<<<<<<<<<<<< * ENTRY ZZFSET ( SETFNC, SETSIZ ) * FUNCNO = SETFNC SIZE = SETSIZ FIRST = .TRUE. * RETURN * C>>>>>>>>>>>>>>>>>>>>>>>>>> E N T R Y ZZFPAR <<<<<<<<<<<<<<<<<<<<<<<<<< * ENTRY ZZFPAR ( SARG ) * DO 80000 I = 1, FNO * FARG ( I ) = SARG ( I ) * 80000 CONTINUE * RETURN * C=============================== E X I T =============================== * 90000 IFG = RET GFIRST = .FALSE. * RETURN * C============================ F O R M A T S ============================ C C THERE ARE NO FORMATS USED. C C================================ E N D ================================ * END SUBROUTINE ZZLCUC ( LINE ) * C============== A R G U M E N T D E C L A R A T I O N S ============== * CHARACTER * (*) LINE * C============================= S T A T U S ============================= C C SYSTEM DEPENDENCE: THE ACTIONS DEPEND UPON THE UNDERLYING C CHARACTER SET. MOST VERSIONS ASSUME AN ASCII C OR EBCDIC CHARACTER SET WITH EACH CHARACTER C OCCUPYING ONE BYTE. THE CYBER VERSION ASSUMES 6- C 12 BIT CODES FOR UPPER/LOWER CASE. ALL THAT HAP- C PENS IS THAT THE ESCAPE CODES @ AND ^ ARE REMOVED. C C THE CURRENT VERSION MAY BE PROCESSED TO PRODUCE C CODE FOR ONE OF THE FOLLOWING: C C THIS VERSION IS FOR SUN4 C CONVERSION FOR SINGLE OR DOUBLE PRECISION: NOT REQUIRED C C REVISION STATUS: DATE AUTHOR VERSION C C JUN. 16, 1988 A. BUCKLEY 1.2 C C======================== D E S C R I P T I O N ======================== C C LINE CONTAINS THE CHARACTERS TO BE CONVERTED TO UPPER CASE. C NOTE THAT ON SOME SYSTEMS (NOTABLY THE OLD ARCHITECTURE C CDC CYBERS), UPPER AND LOWER CASE CHARACTERS HAVE DIFFERENT C LENGTHS, SO THE LINE LENGTH MAY BE CHANGED BY THIS ROUTINE. C C THE ROUTINE ASSUMES THAT THE VALUE C C ICHAR ( UPPER-CASE OF LETTER ) - ICHAR ( LOWER-CASE LETTER ) C C IS FIXED, I.E. IT IS THE SAME FOR EACH LETTER =A,B,...,Z. C C THIS SHOULD BE TRUE ON MOST EBCDIC OR ASCII SYSTEMS. C C======================= E N T R Y P O I N T S ======================= C C ONLY THE NATURAL ENTRY POINT ZZLCUC C C======================== S U B R O U T I N E S ======================== C C INDEX ...GENERIC C ZZLENG ...NON-BLANK LENGTH OF A LINE C C CHAR, ICHAR ...INTRINSIC C C========================= P A R A M E T E R S ========================= * CHARACTER * (*) LCALPH PARAMETER ( LCALPH = 'abcdefghijklmnopqrstuvwxyz' ) * C================= L O C A L D E C L A R A T I O N S ================= * LOGICAL FIRST * CHARACTER CH * INTEGER I, ZZLENG, WLEN, SHIFT * C=============================== S A V E =============================== SAVE FIRST, SHIFT * C======================= E Q U I V A L E N C E S ======================= C C THERE ARE NO EQUIVALENCES. C C============================= C O M M O N ============================= C C THERE ARE NO COMMON BLOCKS. C C=============================== D A T A =============================== DATA FIRST/.TRUE./ * C========================== E X E C U T I O N ========================== * WLEN = ZZLENG(LINE) * IF ( FIRST ) THEN SHIFT = ICHAR( 'A' ) - ICHAR( 'a' ) FIRST = .FALSE. ENDIF * DO 100 I = 1, WLEN * CH = LINE( I:I ) * IF ( INDEX( LCALPH, CH ) .NE. 0 ) THEN * LINE( I:I ) = CHAR( ICHAR( CH ) + SHIFT ) * ENDIF * 100 CONTINUE * GO TO 90000 * C=============================== E X I T =============================== * 90000 RETURN * C============================ F O R M A T S ============================ C C THERE ARE NONE. C C================================ E N D ================================ * END INTEGER FUNCTION ZZLENG (LINE) * C============== A R G U M E N T D E C L A R A T I O N S ============== * CHARACTER*(*) LINE * C============================= S T A T U S ============================= C C SYSTEM DEPENDENCE: NONE C C CONVERSION FOR SINGLE OR DOUBLE PRECISION: NOT REQUIRED C C REVISION STATUS: DATE AUTHOR VERSION C C JUN. 1, 1985 A. BUCKLEY 1.0 C C======================== D E S C R I P T I O N ======================== C C THIS ROUTINE DETERMINES THE POSITION OF THE LAST NONBLANK C CHARACTER IN THE STRING LINE. IF THE LINE IS ENTIRELY C BLANK, THEN ZZLENG IS SET TO 0. C C======================= E N T R Y P O I N T S ======================= C C ONLY THE NATURAL ENTRY POINT ZZLENG C C======================== S U B R O U T I N E S ======================== C C LEN ...INTRINSIC C C========================= P A R A M E T E R S ========================= * CHARACTER*(*) BLANK PARAMETER ( BLANK = ' ' ) * C================= L O C A L D E C L A R A T I O N S ================= * INTEGER I * C=============================== S A V E =============================== C C THERE ARE NO SAVE VARIABLES. C C======================= E Q U I V A L E N C E S ======================= C C THERE ARE NO EQUIVALENCES. C C============================= C O M M O N ============================= C C THERE ARE NO COMMON BLOCKS. C C=============================== D A T A =============================== C C THERE ARE NO DATA VALUES. C C========================== E X E C U T I O N ========================== * ZZLENG = 0 * DO 1000 I = LEN(LINE), 1, -1 * IF ( LINE(I:I) .NE. BLANK ) THEN ZZLENG = I GOTO 90000 ENDIF * 1000 CONTINUE * C=============================== E X I T =============================== * 90000 RETURN * C============================ F O R M A T S ============================ C C THERE ARE NO FORMATS USED. C C================================ E N D ================================ * END SUBROUTINE ZZSHFT (STRING, FROM, TO, NUMBER ) * C============== A R G U M E N T D E C L A R A T I O N S ============== * INTEGER FROM, TO, NUMBER * CHARACTER *(*) STRING * C============================= S T A T U S ============================= C C SYSTEM DEPENDENCE: NONE C C CONVERSION FOR SINGLE OR DOUBLE PRECISION: NOT REQUIRED C C REVISION STATUS: DATE AUTHOR VERSION C C JUN. 1, 1985 A. BUCKLEY 1.0 C C======================== D E S C R I P T I O N ======================== C C THIS ROUTINE PERFORMS A SHIFT OF CHARACTERS WITHIN STRING. THE C NUMBER OF CHARACTERS SHIFTED IS NUMBER AND THEY ARE SHIFTED SO C THAT THE CHARACTER IN POSITION FROM IS MOVED TO POSITION TO. C CHARACTERS IN THE TO POSITION ARE OVERWRITTEN. BLANKS REPLACE C CHARACTERS IN THE FROM POSITION. SHIFTING MAY BE LEFT OR RIGHT, C AND THE FROM AND TO POSITIONS MAY OVERLAP. CARE IS TAKEN NOT C TO ALTER OR USE ANY CHARACTERS BEYOND THE DEFINED LIMITS C OF THE STRING. C C======================= E N T R Y P O I N T S ======================= C C ONLY THE NATURAL ENTRY POINT ZZSHFT C C======================== S U B R O U T I N E S ======================== C C LEN MIN MAX ...INTRINSIC C C========================= P A R A M E T E R S ========================= * CHARACTER *(*) BLANK PARAMETER ( BLANK = ' ' ) * C================= L O C A L D E C L A R A T I O N S ================= * INTEGER N, SHIFT, INCR, I, IS, IE, IBS, ETO, EFROM, K, SLEN * C=============================== S A V E =============================== C C THERE ARE NO SAVE VARIABLES. C C======================= E Q U I V A L E N C E S ======================= C C THERE ARE NO EQUIVALENCES. C C============================= C O M M O N ============================= C C THERE ARE NO COMMON BLOCKS. C C=============================== D A T A =============================== C C THERE ARE NO DATA VALUES. C C========================== E X E C U T I O N ========================== * SLEN = LEN (STRING) N = NUMBER - 1 SHIFT = FROM - TO * IF ( FROM .NE. TO ) THEN * IF ( TO .LE. FROM ) THEN * INCR = 1 IS = MIN( FROM+MAX(0,1-TO), SLEN+1 ) IE = MIN( FROM+N, SLEN ) IBS = MAX( IE-SHIFT+1, MAX(FROM,0) ) * ELSE * INCR = -1 ETO = TO + N EFROM = FROM + N IS = MAX( EFROM - MAX(0,ETO-SLEN) , 0 ) IE = MAX(FROM , 0) IBS = MIN( TO-1 , MIN(EFROM,SLEN) ) * ENDIF * DO 1000 I=IS,IE,INCR K = I - SHIFT STRING(K:K) = STRING(I:I) 1000 CONTINUE * DO 2000 I=IBS,IE,INCR STRING(I:I) = BLANK 2000 CONTINUE * ENDIF * GOTO 90000 * C=============================== E X I T =============================== * 90000 RETURN * C============================ F O R M A T S ============================ C C THERE ARE NO FORMATS USED. C C================================ E N D ================================ * END SUBROUTINE ZZTIME ( CHTIME ) * C============== A R G U M E N T D E C L A R A T I O N S ============== * CHARACTER *(*) CHTIME * C============================= S T A T U S ============================= C C SYSTEM DEPENDENCE: SYSTEM DEPENDENT ROUTINE FOR TIME OF DAY. C C THE CURRENT VERSION MAY BE PROCESSED TO PRODUCE C CODE FOR ONE OF THE FOLLOWING: C C THIS VERSION IS FOR SUN4 C C CONVERSION FOR SINGLE OR DOUBLE PRECISION: NOT REQUIRED C C REVISION STATUS: DATE AUTHOR VERSION C C JUN. 16, 1988 A. BUCKLEY 1.2 C C======================== D E S C R I P T I O N ======================== C C THIS ROUTINE MUST CALL A SYSTEM ROUTINE TO GET THE CURRENT TIME. C IT MUST RETURN THIS TIME IN THE CHARACTER VARIABLE CHTIME C IN THE FORM C C (HH+MM+SS) C C AS REQUIRED BY THE ROUTINE ZZDTTM. CHTIME MUST BE OF LENGTH 10. C ONLY THE SIX CHARACTERS HH MM SS ARE ACTUALLY USED. THE OTHERS C CAN BE ANYTHING, BUT THE HH, MM AND SS MUST BE IN THE CORRECT C POSITION. HH SHOULD BE BETWEEN 0 AND 23, I.E. BASED ON THE 24 HOUR C CLOCK. OF COURSE MM AND SS MUST BE IN THE RANGE 0 TO 59. C C THIS VERSION CALLS A SYSTEM DEPENDENT ROUTINE. C TO IMPLEMENT THIS ON ANOTHER SYSTEM, ONE MAY EITHER C C (A) INCORPORATE AN ALTERNATE VERSION OF ZZTIME; C C (B) USE THE "DUMMY" VERSION OF THIS ROUTINE WITH THE SINGLE C EXECUTABLE STATEMENT CHTIME=' ', IN WHICH CASE NO C TIME INFORMATION WILL APPEAR IN THE OUTPUT. C C======================= E N T R Y P O I N T S ======================= C C ONLY THE NATURAL ENTRY POINT ZZTIME C C======================== S U B R O U T I N E S ======================== C C SYSTEM ROUTINE TO GET TIME OF DAY. C C========================= P A R A M E T E R S ========================= C C THERE ARE NO PARAMETERS DEFINED. C C================= L O C A L D E C L A R A T I O N S ================= * CHARACTER * 24 UNXTIM * C=============================== S A V E =============================== C C THERE ARE NO SAVE VARIABLES. C C======================= E Q U I V A L E N C E S ======================= C C THERE ARE NO EQUIVALENCES. C C============================= C O M M O N ============================= C C THERE ARE NO COMMON BLOCKS. C C=============================== D A T A =============================== C C THERE ARE NO DATA STATEMENTS. C C========================== E X E C U T I O N ========================== * CALL FDATE (UNXTIM) * CHTIME(2:9) = UNXTIM(12:19) * C=============================== E X I T =============================== * 90000 RETURN * C============================ F O R M A T S ============================ C C THERE ARE NONE. C C================================ E N D ================================ * END SUBROUTINE ZZEVAL ( ZZUFNC, N, X, F, G, INDIC, IW, RW, DW ) * C============== A R G U M E N T D E C L A R A T I O N S ============== * EXTERNAL ZZUFNC * INTEGER INDIC, N, IW(*) * DOUBLE PRECISION F, X(N), G(N), ZZUFNC C!!!! REAL F, X(N), G(N) * DOUBLE PRECISION DW(*) REAL RW(*) * C============================= S T A T U S ============================= C C SYSTEM DEPENDENCE: NONE C C CONVERSION FOR SINGLE OR DOUBLE PRECISION: REQUIRED (SEE CONVRT) C C IGNORE LINES BEGINNING WITH "C!!!!" . C C THIS VERSION IS IN D O U B L E PRECISION. C!!!! THIS VERSION IS IN S I N G L E PRECISION. C C REVISION STATUS: DATE AUTHOR VERSION C C MAY. 21, 1987 A. BUCKLEY 1.2 C C======================== D E S C R I P T I O N ======================== C C THIS SUBROUTINE EVALUATES A TEST FUNCTION AT THE GIVEN C POINT "X". IT RETURNS THE VALUE OF THE FUNCTION AND / OR C THE VALUE OF THE GRADIENT AT X. IT ALLOWS THE APPLICATION OF A C NONLINEAR SCALING TO THE FUNCTION IF DESIRED (SEE FSCALE BELOW). C IT ALSO ALLOWS THE USE OF FINITE DIFFERENCES (SEE SDRVMD BELOW). C IT CAN ALSO ACT AS A NOOP, I.E. AS A DO NOTHING ROUTINE; (SEE C INDIC BELOW). C C-----ON ENTRY: C C ZZUFNC THE NAME OF THE FUNCTION TO EVALUATE. THERE C MUST BE A SUBROUTINE PROVIDED OF THE FORM C C SUBROUTINE ZZUFNC(INDIC,N,X,F,G,IW,RW,DW) C C (WHERE N, X, F, G, INDIC, IW, RW AND DW HAVE THE C SAME MEANING AS IN THIS SUBROUTINE ZZEVAL.) C C N THE DIMENSION OF THE PROBLEM, I.E. THE C NUMBER OF VARIABLES IN THE FUNCTION ZZUFNC. C C X CONTAINS THE VALUE OF THE N-COORDINATES X[1],...,X[N] C AT WHICH TO EVALUATE THE FUNCTION. C C INDIC = DOF ONLY EVALUATE THE FUNCTION. C = DOFG EVALUATE BOTH. C = DOG ONLY EVALUATE THE GRADIENT. C = NONE ACTUALLY, IF INDIC HAS ANY VALUE OTHER THAN C ONE OF THE FIRST THREE, THEN JUST CALL ZZUFNC C WITH THIS SAME CODE FOR INDIC; I.E. ZZEVAL SHOULD DO C NOTHING. THIS IS INTENDED FOR THE CONVENIENCE OF THE C WRITER OF ZZUFNC. C C NOTE THAT THE VALUES OF THESE CODES CAN BE REDEFINED C THROUGH THE ENTRY POINT ZZFDEF BELOW. DEFAULT VALUES C ARE GIVEN IN THE PARAMETER SECTION BELOW. C C IW THESE ARE 3 WORK ARRAYS WHICH ARE NOT USED AT ALL BY C RW ZZEVAL, BUT WHICH ARE JUST PASSED TO THE USER'S C DW ROUTINE ZZUFNC TO BE USED AS DESIRED. WITH THESE ARRAYS C AVAILABLE, IT IS OFTEN NOT NECESSARY TO USE REVERSE C COMMUNICATION. NOTE THAT THERE IS ONE AVAILABLE OF C EACH BASIC NUMERIC TYPE. C C-----ON EXIT: C C F CONTAINS THE FUNCTION VALUE (WITH THE SCALING C APPLIED IF REQUIRED). C C G CONTAINS THE GRADIENT VALUE (WITH THE SCALING C APPLIED IF REQUIRED). C C NEITHER F NOR G IS REFERENCED UNLESS ITS VALUE IS REQUESTED. C C INDIC = OK THE REQUEST MADE ON THE CALL WAS COMPLETED SATIS- C FACTORILY. F AND/OR G ARE AVAILABLE AS REQUESTED. C ABORT THE MINIMIZATION ROUTINE WHICH CALLED ZZEVAL IS C HEREBY REQUESTED TO EXIT IMMEDIATELY TO THE ROUTINE C WHICH CALLED IT. THIS CAN BE USED BY THE ROUTINE C ZZUFNC TO TRIGGER PREMATURE TERMINATION DUE TO C CIRCUMSTANCES OF WHICH THE MINIMIZATION ROUTINE MAY C NOT BE AWARE. C LIMIT TERMINATE THE MINIMIZATION; THE PRESET LIMIT ON THE C NUMBER OF FUNCTION EVALUATIONS ALLOWED HAS BEEN C EXCEEDED. SEE MAXFN BELOW. C NOF THE FUNCTION VALUE COULD NOT BE DETERMINED. C NOG THE GRADIENT VALUE COULD NOT BE DETERMINED. C NOFG NEITHER F NOR G COULD BE EVALUATED. C C THESE CODES CAN BE REDEFINED THROUGH AN ENTRY POINT BELOW, C AND HAVE DEFAULT VALUES SPECIFIED IN THE PARAMETER SECTION. C C-----SET THROUGH ENTRY POINT CALLS. C C ZZESRT ( FSCALE, SDRVMD, MAXFN ) M A N D A T O R Y C C THIS IS CALLED BEFORE MINIMIZING EACH FUNCTION. THIS CALL C IS M A N D A T O R Y. C C FSCALE CONTROLS THE NONLINEAR SCALING OF ZZUFNC. C C = 0 NO EFFECT. C C = K>0 THIS ROUTINE COMPUTES AND RETURNS FF( ZZUFNC(X) ), C WHERE FF IS THE K-TH OF THE NONLINEAR FUNCTIONS C OF ONE VARIABLE DEFINED IN THE ROUTINE ZZSCAL. C C NOTE THAT FOR CERTAIN SCALINGS, IF YOU CALL ZZEVAL C JUST FOR A GRADIENT VALUE, IT MAY BE NECESSARY TO C REQUEST A FUNCTION VALUE AS WELL IN ORDER TO DO THE C SCALING. THAT FUNCTION VALUE WILL NOT BE PASSED BACK. C THOSE WHICH DO NOT REQUIRE F FOR THE SCALING ARE THOSE C WITH K = 1,2,..,REQF - 1. FOR K = REQF,..., THE VALUE C OF F IS NECESSARY. C C SDRVMD THIS SPECIFIES THE METHOD BY WHICH DERIVATIVES ARE C TO BE COMPUTED, WHEN REQUESTED. THE CHOICE IS BETWEEN C C CANAL USE ANALYTIC FORMULAE WHICH MUST BE CODED AND C AVAILABLE IN THE USER ROUTINE ZZUFNC. C C CDIFF USE FINITE DIFFERENCE APPROXIMATIONS. IN THIS CASE, C THE USER ROUTINE MAY IGNORE CALLS WITH INDIC <> JUSTF, C AND NEED ONLY BE ABLE TO COMPUTE FUNCTION VALUES. C FURTHER COMMENTS APPEAR IN THE DISCUSSION OF FINITE C DIFFERENCE COMPUTATIONS (BELOW). C C CTEST IN THIS CASE BOTH ANALYTIC AND FINITE DIFFERENCES C ARE COMPUTED. THEY ARE THEN COMPARED AND A RECORD C IS KEPT TO SEE TO WHAT EXTENT THEY DISAGREE. A C RECORD OF THE LEVEL OF AGREEMENT IS AVAILABLE C THROUGH THE ENTRY POINT ZZECHK GIVEN BELOW. A MORE C COMPLETE DESCRIPTION IS ALSO GIVEN WHERE ZZECHK IS C DISCUSSED BELOW. C C CFIRST THIS CASE IS PRECISELY THE SAME AS FOR CTEST, WITH C THE SOLE EXCEPTION THAT THE TESTING ONLY TAKES PLACE C ON THE FIRST CALL TO ZZEVAL. C C THE INTEGER VALUES OF THE CODES FOR CANAL, ETC ARE C SET IN THE PARAMETER SECTION BELOW. THEY MAY BE C RESET VIA THE ENTRY POINT ZZEDEF DESCRIBED BELOW. C C MAXFN THE MAXIMUM VALUE ALLOWED FOR THE COUNT IFNCT. C C <= 0 ON ENTRY SPECIFIES NO MAXIMUM, I.E. MAXFN IS C IGNORED. C C = K>0 SPECIFIES THE MAXIMUM NUMBER OF TIMES THAT ZZUFNC C MAY BE CALLED. IF THE FUNCTION EVALUATION COUNT C IN IFNCT IS GREATER THAN OR EQUAL TO MAXFN ON C ENTRY TO ZZEVAL, THEN THE FUNCTION IS NOT C EVALUATED AND THE RETURN CODE INDIC IS SET AS C ABOVE. NOTE THAT THE COUNT IN IFNCT DOES N O T C INCLUDE FUNCTION EVALUATIONS USED FOR COMPUTING C FINITE DIFFERENCE GRADIENTS. C C C THE NEXT FOUR PARAMETERS ARE NOT IN THE CALLING SEQUENCE OF C ZZESRT, BUT THEY ARE INITIALIZED WHEN ZZESRT IS CALLED. C C IFNCT COUNTS THE NUMBER OF TIMES THE ROUTINE IS CALLED C TO EVALUATE THE FUNCTION. IT IS INITIALIZED TO 0 C DURING THE CALL TO ZZESRT. C C IGRCT COUNTS THE NUMBER OF TIMES THE ROUTINE IS CALLED C TO EVALUATE THE GRADIENT. IT IS INITIALIZED TO 0 C DURING THE CALL TO ZZESRT. C C FTIME RECORDS THE TIME ACCUMULATED IN EVALUATING THE C FUNCTION AND/OR THE GRADIENT. IT IS PRESET TO ZERO C WHEN ZZESRT IS CALLED. THE TIME USED IN THE FINAL C SCALING IS INCLUDED IN THE TIMING WHEN THE VALUE OF C FSCALE IS NON-ZERO. TIMING COMMENCES ON ENTRY TO C ZZEVAL, AND ENDS JUST BEFORE RETURN FROM ZZEVAL. C C ERR THE ESTIMATE OF THE ERROR BETWEEN THE ANALYTIC C AND DIFFERENCE VALUES FOR THE GRADIENT IS RECORDED C IN A SET OF VARIABLES ERR, SERR, DCNT, INDEX AND C GCNT, SO THESE ARE INITIALIZED TO 0. C C ZZESET ( TRF, TRG, ITRUN ) C C THIS IS CALLED BEFORE USING ZZEVAL (THESE VALUES ALSO HAVE C INTERNALLY SET DEFAULT VALUES GIVEN IN [..], SO THE CALL TO C ZZESET IS NOT MANDATORY.) C C TRF = TRUE IF THE FUNCTION VALUE IS TO BE PRINTED [FALSE] C TRG = TRUE IF THE GRADIENT VALUE IS TO BE PRINTED [FALSE] C C ITRUN THE UNIT NUMBER FOR OUTPUT OF COMPUTED VALUES [6] C C NOTE THAT AN ERROR MESSAGE IS PRINTED WHEN THE MAXIMUM NUMBER C OF FUNCTION EVALUATIONS IS EXCEEDED, PROVIDED EITHER TRF OR C TRG IS TRUE. C C ZZEDDF ( SANAL, SDIFF, STEST, SFIRST ) C C THIS MAY BE CALLED BEFORE USING ZZEVAL, AS FOR ZZESET. THIS C ALLOWS THE CODES FOR ANAL, ETC., TO BE REDEFINED. ALL HAVE C DEFAULTS, SO THIS CALL IS NOT MANDATORY. C C SANAL THE INTEGER VALUE FOR THE CODE FOR USING ANALYTIC C DERIVATIVES [CANAL]. C SDIFF THE INTEGER VALUE FOR THE CODE FOR USING FINITE C DIFFERENCES TO APPROXIMATE DERIVATIVES [CDIFF]. C STEST THE INTEGER VALUE FOR THE CODE FOR USING BOTH CANAL C AND CDIFF AND DOING A TEST FOR AGREEMENT [CTEST]. C SFIRST THE INTEGER VALUE FOR THE CODE FOR USING BOTH CANAL C AND CDIFF ON THE FIRST ITERATION ONLY. C C ZZEFDF ( SDOF, SDOG, SDOFG, SNONE ) C C THIS MAY BE CALLED BEFORE USING ZZEVAL, JUST AS FOR ZZEDEF. C C DOF THE CODE INDICATING THAT JUST THE FUNCTION VALUE IS C DESIRED. [JUSTF] C DOG THE CODE INDICATING THAT JUST THE GRADIENT VALUE IS C DESIRED. [JUSTG] C DOFG THE CODE INDICATING THAT BOTH THE FUNCTION AND GRADIENT C VALUES ARE DESIRED. [BOTH] C SNONE THE CODE INDICATING THAT NO ACTION IS TO BE TAKEN AND C THAT ZZUFNC SHOULD BE CALLED WITH NO OTHER PROCESSING. C C ZZERDF ( OK, LIMIT, ABORT, NOF, NOG, NOFG ) C C THIS MAY BE CALLED BEFORE USING ZZEVAL, JUST AS FOR ZZEDEF. C C OK THIS CODE INDICATES THAT THE REQUEST WAS SUCCESSFULLY DONE. C ABORT THIS MEANS THAT THE CALLING ROUTINE SHOULD IMMEDIATELY C TERMINATE THE MINIMIZATION AND RETURN TO THE ROUTINE WHICH C CALLED IT. C LIMIT THIS MEANS THAT THE ALLOWED NUMBER OF FUNCTION EVALUATIONS C HAS BEEN EXCEEDED. C NOF THIS MEANS THAT ZZEVAL WAS UNABLE TO SUCCESSFULLY EVALUATE C THE FUNCTION. C NOG THIS MEANS THAT ZZEVAL WAS UNABLE TO SUCCESSFULLY EVALUATE C THE GRADIENT. C NOFG THIS MEANS THAT ZZEVAL WAS UNABLE TO OBTAIN EITHER A C FUNCTION OR GRADIENT VALUE. C C-----AVAILABLE THROUGH ENTRY POINT CALLS AFTER A FUNCTION HAS BEEN C MINIMIZED: C C ZZEGET ( FNCT, GRCT, TIME ) C C THIS MAY BE CALLED AFTER MINIMIZING A FUNCTION TO OBTAIN SOME C SIMPLE STATISTICS WHICH HAVE BEEN ACCUMULATED SINCE THE LAST C CALL TO ZZESRT. THESE ARE C C FNCT THE NUMBER OF CALLS TO EVALUATE THE FUNCTION, I.E. C CALLS WITH INDIC = JUSTF OR BOTH. C C GRCT THE NUMBER OF CALLS TO EVALUATE THE GRADIENT, I.E. C CALLS WITH INDIC = JUSTG OR BOTH. C C TIME THE AMOUNT OF CPU TIME SPENT IN ZZEVAL. C C ZZECHK ( ERR, AVERR, INDX, ITERAT ) C C C THIS MAY ALSO BE CALLED AFTER A SEQUENCE OF CALLS TO ZZEVAL. C IT GIVES AN ESTIMATE OF THE AGREEMENT BETWEEN ANALYTIC AND C FINITE DIFFERENCE DERIVATIVES. OF COURSE THESE VALUES ARE ONLY C DEFINED IF SDRVMD = CTEST. C C ERR WILL BE RETURNED AS AN ESTIMATE OF THE LARGEST ERROR C WHICH OCCURRED AND AVERR IS AN ESTIMATE OF THE AVERAGE NUMBER OF C DECIMAL DIGITS OF AGREEMENT BETWEEN THE COMPONENTS OF THE ANALYTIC C AND DIFFERENCE DERIVATIVES. C C TO BE SPECIFIC, WHEN IN TEST MODE, EACH COMPONENT OF THE C ANALYTIC DERIVATIVE IS COMPUTED, AND THAT IS RETURNED IN G AS THE C GRADIENT. AS WELL, FOR EACH COMPONENT, A FINITE DIFFERENCE C APPROXIMATION IS COMPUTED (AS DESCRIBED BELOW) AND THE RELATIVE C DIFFERENCE BETWEEN THAT AND THE ANALYTIC COMPONENT IS DETERMINED. C THIS QUANTITY IS MONITORED, AND THE LARGEST SUCH VALUE IS C RECORDED. IN ADDITION, INDX INDICATES IN WHICH COMPONENT OF THAT C GRADIENT THE ERROR OCCURRED AND ITERAT TELLS WHICH GRADIENT C EVALUATION WAS IN PROGRESS WHEN THE ERROR OCCURRED; I.E. ITERAT C JUST RECORDS THE CURRENT VALUE OF IGRCT. NOTE THAT INDX C AND ITERAT ONLY REFER TO THE POINT AT WHICH THE LARGEST C ERROR OCCURRED. C C IF THE FUNCTION AND GRADIENT EVALUATIONS ARE CORRECT, ONE C WOULD NORMALLY EXPECT THE RELATIVE ERROR TO BE OF THE ORDER OF C 10**-(T/2), WHERE T IS THE NUMBER OF FIGURES OF RELATIVE C ACCURACY OF THE MACHINE IN USE. HOWEVER, AS THE MINIMUM IS C APPROACHED AND THE GRADIENT COMPONENTS GENERALLY BECOME VERY C SMALL, THIS RELATIVE ACCURACY MAY BE MUCH WORSE THAN EXPECTED. C THEREFORE WE ALSO MAINTAIN AN ESTIMATE OF THE AVERAGE AGREEMENT. C HERE, FOR EACH COMPONENT OF EACH GRADIENT COMPUTATION, WE COMPUTE C THE BASE 10 LOG OF THE RELATIVE ACCURACY; THIS IS ROUGHLY THE C NUMBER OF SIGNIFICANT FIGURES OF AGREEMENT BETWEEN THE TWO VALUES. C THIS QUANTITY IS MONITORED AND AVERR IS RETURNED AS THE AVERAGE C VALUE OF THE NUMBER OF SIGNIFICANT FIGURES OF AGREEMENT. C C WHEN FUNCTION AND GRADIENT COMPUTATIONS ARE CORRECT, ERR WILL C GENERALLY BE AT LEAST AS SMALL AS 10**(-T/2), ALTHOUGH IT CAN BE C MORE LIKE 10**(-T/4). GROSS BLUNDERS WILL USUALLY GIVE ERR A C VALUE VERY NEAR TO 1, BUT NOT ALWAYS. IF ALL IS WELL, AVERR WILL C USUALLY BE ABOUT T/2; BLUNDERS WILL OFTEN RESULT IN AVERR BEING C NEAR 0 OR 1. C C-----FINITE DIFFERENCE COMPUTATIONS C C FOR FIRST DERIVATIVES, SIMPLE FORWARD DIFFERENCES ARE USED. C TO ESTIMATE THE I-TH COMPONENT OF THE GRADIENT OF F, WE COMPUTE C C ( F(X + H*E[I]) - F(X) ) / H, C C WHERE H = EPS * ABS(X[I]). WHEN X[I] = 0, WE JUST CHOOSE H = EPS. C HERE EPS IS THE ROOT OF ETA, WHERE ETA DEFINES THE RELATIVE C MACHINE ACCURACY. THIS IS USED WHEN SDRVMD = CDIFF OR CTEST. C C WHEN SDRVMD = CTEST, MORE INFORMATION IS REQUIRED; THUS WE C ALSO COMPUTE F(X + SQRT(H)*E[I]). THIS MEANS THAT WHEN IN TEST C MODE, TWICE AS MANY FUNCTION EVALUATIONS ARE NEEDED. THIS IS C REQUIRED TO ELIMINATE SCALING EFFECTS IN THE ESTIMATE OF FIGURES C OF AGREEMENT. C C======================= E N T R Y P O I N T S ======================= C C ZZEVAL THE NATURAL ENTRY POINT. C C ZZESRT AN ENTRY TO INITIALIZE FOR TESTING EACH FUNCTION. C ZZESET AN ENTRY TO INITIALIZE CONTROL PARAMETERS. C C ZZEDDF AN ENTRY TO REDEFINE CODES FOR DERIVATIVE CALCULATIONS. C ZZEFDF AN ENTRY TO REDEFINE CODES FOR FUNCTION EVALUATIONS. C ZZERDF AN ENTRY TO REDEFINE CODES FOR RETURN CODES. C C ZZEGET AN ENTRY TO RETURN FUNCTION/GRADIENT COUNTS AND TIME. C ZZECHK RETURNS THE ERROR VALUE IF IN TESTMODE. C C======================== S U B R O U T I N E S ======================== C C ZZSECS ...FOR FUNCTION TIMING. C ZZMPAR ...FOR MACHINE PARAMETERS. C ZZUFNC ...THE USER ROUTINE. C ZZSCAL ...PERFORMS SCALING. C C SQRT, MAX, ABS ...INTRINSIC FUNCTIONS. C LOG10, MIN, SIGN ...INTRINSIC FUNCTIONS. C C========================= P A R A M E T E R S ========================= * DOUBLE PRECISION ZERO, ONE, TWO, THREE C!!!! REAL ZERO, ONE, TWO, THREE PARAMETER ( ZERO = 0D0, ONE = 1D0, TWO = 2D0, THREE = 3D0) * DOUBLE PRECISION FOUR, FIVE, SIX, SEVEN C!!!! REAL FOUR, FIVE, SIX, SEVEN PARAMETER ( FOUR = 4D0, FIVE = 5D0,SIX = 6D0, SEVEN = 7D0) * DOUBLE PRECISION EIGHT, NINE, TEN C!!!! REAL EIGHT, NINE, TEN PARAMETER ( EIGHT = 8D0, NINE = 8D0, TEN = 10D0 ) DOUBLE PRECISION R11, R12, R13, R14 C!!!! REAL R11, R12, R13, R14 PARAMETER ( R11 = 11D0, R12 = 12D0, R13 = 13D0,R14 = 14D0) * DOUBLE PRECISION R15, R16, R17, R18 C!!!! REAL R15, R16, R17, R18 PARAMETER ( R15 = 15D0, R16 = 16D0, R17 = 17D0,R18 = 18D0) * DOUBLE PRECISION R19, R20, R25, R29 C!!!! REAL R19, R20, R25, R29 PARAMETER ( R19 = 19D0, R20 = 20D0, R25 = 25D0,R29 = 29D0) * DOUBLE PRECISION R32, R36, R40, R42 C!!!! REAL R32, R36, R40, R42 PARAMETER ( R32 = 32D0, R36 = 36D0, R40 = 40D0,R42 = 42D0) * DOUBLE PRECISION R45, R49 C!!!! REAL R45, R49 PARAMETER ( R45 = 45D0, R49 = 49D0 ) * DOUBLE PRECISION R50, R56, R84, R90 C!!!! REAL R50, R56, R84, R90 PARAMETER ( R50 = 50D0, R56 = 56D0, R84 = 84D0,R90 = 90D0) * DOUBLE PRECISION R100, R180, R200 C!!!! REAL R100, R180, R200 PARAMETER ( R100 = 100D0, R180 = 180D0, R200 = 200D0 ) * DOUBLE PRECISION R256, R360, R400 C!!!! REAL R256, R360, R400 PARAMETER ( R256 = 256D0, R360 = 360D0, R400 = 400D0 ) * DOUBLE PRECISION R600, R681, R991 C!!!! REAL R600, R681, R991 PARAMETER ( R600 = 600D0, R681 = 681D0, R991 = 991D0 ) * DOUBLE PRECISION R1162, R2324 C!!!! REAL R1162, R2324 PARAMETER ( R1162 = 1162D0, R2324 = 2324D0 ) * DOUBLE PRECISION R10000, R40000 C!!!! REAL R10000, R40000 PARAMETER ( R10000 = 10000D0, R40000 = 40000D0 ) * INTEGER XEPS, XSMALL, XBIG PARAMETER ( XEPS = 1, XSMALL = 2, XBIG = 3 ) * INTEGER REQF PARAMETER ( REQF = 2 ) * C DEFINE THE DERIVATIVE CODES * * INTEGER CANAL, CDIFF, CTEST, CFIRST PARAMETER ( CANAL = 1, CDIFF = 2, CTEST = 3, CFIRST = 4 ) * C DEFINE THE FUNCTION CODES * * INTEGER JUSTF, BOTH, JUSTG, NOOP PARAMETER ( JUSTF = 1, BOTH = 0, JUSTG = -1, NOOP = 2 ) * C DEFINE THE RETURN CODES * C THE RETURN CODES TO BE USED BY THE FUNCTION EVALUATION ROUTINE C TO INDICATE TO THE MINIMIZATION ROUTINE WHETHER OR NOT THE CALL C WAS SUCCESSFUL. * INTEGER COK, CABORT, CLIMIT PARAMETER ( COK = 0, CABORT = -1, CLIMIT = -2 ) * INTEGER CNOF, CNOG, CNOFG PARAMETER ( CNOF = -3, CNOG = -4, CNOFG = -5 ) * C================= L O C A L D E C L A R A T I O N S ================= * LOGICAL TRF, TRG, FIRST, FONLY, GONLY, VALID, BAD, TRTEST, FCALL * INTEGER I, IFNCT, FSCALE, IGRCT, SDRVMD INTEGER ITRUN, MAXFN, DERVMD INTEGER CASE, CALLS, COUNT, INDEX, GCNT, DCNT * DOUBLE PRECISION FT, FV, FTIME, TT, SCALE, SERR, RH C!!!! REAL FT, FV, FTIME, TT, SCALE, SERR, RH * DOUBLE PRECISION FVAL, FVAL2, ERR, ETA, EPS, H, ZZMPAR, TERR C!!!! REAL FVAL, FVAL2, ERR, ETA, EPS, H, ZZMPAR, TERR * C-----DECLARATIONS FOR ENTRY POINT DUMMY ARGUMENTS. * INTEGER DITRUN, FSCAL, MAXM, FNCT, GRCT INTEGER DDERV, INDX, ITERAT * INTEGER ANAL, DIFF, TEST, DOF, DOG, DOFG, NONE, TFIRST INTEGER SANAL, SDIFF, STEST, SDOF, SDOG, SDOFG, SNONE, SFIRST * INTEGER OK, ABORT, LIMIT, NOF, NOG, NOFG INTEGER SOK, SABORT, SLIMIT, SNOF, SNOG, SNOFG * LOGICAL DTRF, DTRG, DTRTST * DOUBLE PRECISION TIME, ERROR, AVERR C!!!! REAL TIME, ERROR, AVERR * C=============================== S A V E =============================== * SAVE ITRUN, FSCALE, IFNCT, IGRCT, SERR, DCNT SAVE TRF, TRG, FTIME, SDRVMD, MAXFN, TRTEST SAVE FIRST, ERR, INDEX, GCNT, EPS, ETA, FCALL SAVE ANAL, DIFF, TEST, TFIRST, DOF, DOG, DOFG, NONE SAVE OK, ABORT, LIMIT, NOF, NOG, NOFG * C======================= E Q U I V A L E N C E S ======================= C C THERE ARE NO EQUIVALENCES. C C============================= C O M M O N ============================= C C THERE ARE NO COMMON BLOCKS. C C=============================== D A T A =============================== * DATA FIRST /.TRUE./, FCALL/.TRUE./ * DATA TRF, TRG, TRTEST / 3 * .FALSE. /, ITRUN / 6 / * DATA SDRVMD/CANAL/, FSCALE/0/, MAXFN/0/ * DATA ANAL/CANAL/, DIFF/ CDIFF/, TEST/ CTEST/, TFIRST/CFIRST/ DATA DOF/JUSTF/, DOG/ JUSTG/, DOFG/ BOTH/, NONE/ NOOP/ DATA OK/ COK/, ABORT/CABORT/, LIMIT/CLIMIT/ DATA NOF/ CNOF/, NOG/ CNOG/, NOFG/ CNOFG/ * C========================== E X E C U T I O N ========================== * C-----STATEMENT FUNCTION. * BAD() = CASE .EQ. ABORT .OR. CASE .EQ. LIMIT .OR. CASE .EQ. NOF - .OR. CASE .EQ. NOG .OR. CASE .EQ. NOFG * C-----FIRST TEST FOR NOOP CALL. * VALID = INDIC .EQ. DOF .OR. INDIC .EQ. DOG .OR. INDIC .EQ. DOFG * IF ( .NOT. VALID ) THEN CALL ZZUFNC ( INDIC, N, X, F, G, IW, RW, DW ) GOTO 90500 ENDIF * IF ( MAXFN .GT. 0 .AND. IFNCT .GE. MAXFN ) THEN GOTO 91000 ENDIF * DERVMD = SDRVMD * IF ( FIRST ) THEN * FIRST = .FALSE. * ETA = ZZMPAR(XEPS) EPS = SQRT (ETA) * ENDIF * IF ( FCALL ) THEN IF ( DERVMD .EQ. TFIRST ) DERVMD = TEST FCALL = .FALSE. ELSE IF ( DERVMD .EQ. TFIRST ) DERVMD = ANAL ENDIF * FONLY = INDIC .EQ. DOF GONLY = INDIC .EQ. DOG CASE = INDIC * CALL ZZSECS (TT) FTIME = FTIME - TT * IF ( .NOT. GONLY ) IFNCT = IFNCT + 1 IF ( .NOT. FONLY ) IGRCT = IGRCT + 1 * C-----FIRST COMPUTE REQUIRED FUNCTION AND/OR GRADIENT VALUES. * C DETERMINE NO OF EXTRA CALLS TO USER ROUTINE WHICH WILL BE NEEDED. * IF ( DERVMD .EQ. ANAL .OR. FONLY ) THEN CALLS = 0 ELSE CALLS = N ENDIF * C FORCE FUNCTION EVALUATION IF REQUIRED FOR SCALING. * IF ( FSCALE .GE. REQF .AND. GONLY ) THEN CASE = DOFG ENDIF * C FIRST COMPUTE F(X) --- AND G(X) IF NEEDED. * CALL ZZUFNC ( CASE, N, X, FVAL, G, IW, RW, DW ) * IF ( BAD() ) THEN INDIC = CASE GOTO 90000 ENDIF * IF ( INDIC .NE. DOG ) THEN FT = FVAL ENDIF * C AFTER FIRST CALL, FUNCTION VALUES ONLY. * C-----DO EXTRA CALLS, IF REQUIRED. * DO 1500 COUNT = 1, CALLS * TT = X(COUNT) * IF ( TT .EQ. ZERO ) THEN H = EPS ELSE H = EPS * ABS( TT ) ENDIF * X(COUNT) = TT + H * C COMPUTE F( X + H * E[COUNT] ) * CASE = DOF CALL ZZUFNC ( CASE, N, X, FVAL, G, IW, RW, DW ) * IF ( BAD() ) THEN INDIC = CASE GOTO 90000 ENDIF * X(COUNT) = TT * IF ( DERVMD .EQ. TEST ) THEN * C ---IF TRACE REQUESTED, PRINT ESTIMATED AND ANALYTIC VALUES. * IF ( TRTEST ) WRITE(ITRUN,99995) G(COUNT),COUNT,(FVAL-FT)/H * C ---ESTIMATE ERROR, AND LEAVE COMPUTED C ANALYTIC GRADIENTS IN G. USE F AT C X + A * E[COUNT], FOR A = H AND SQRT(H). * RH = SQRT(H) X(COUNT) = TT + RH * CASE = DOF CALL ZZUFNC ( CASE, N, X, FVAL2, G, IW, RW, DW ) IF ( BAD() ) THEN INDIC = CASE GOTO 90000 ENDIF * X(COUNT) = TT * IF ( ABS(FVAL2-FT) .GT. R100*ETA*ABS(FT) ) THEN * TERR = (FVAL-FT - H*G(COUNT))/(FVAL2-FT - RH*G(COUNT)) * IF (TT .GT. ONE) TERR = TERR / TT * C TRUNCATE TO INTERVAL [ETA,1]. * TERR = MAX( MIN(ONE,ABS(TERR)), ETA ) * C ESTIMATE NUMBER OF FIGURES OF AGREEMENT. * SERR = SERR - LOG10 (TERR) DCNT = DCNT + 1 * IF (TRTEST) WRITE(ITRUN,99994) TERR,-LOG10(TERR) * IF ( TERR .GT. ABS(ERR) ) THEN * INDEX = COUNT GCNT = IGRCT ERR = SIGN (TERR, ERR) ENDIF * ELSE * C FLAG CASE WHERE THERE IS EXCESSIVE CANCELLATION. * ERR = - ABS(ERR) IF (TRTEST) WRITE(ITRUN,99993) * ENDIF * ELSE * C ---ESTIMATE GRADIENTS USING FORWARD FINITE DIFFERENCE C FORMULAE AND STORE IN G. * G(COUNT) = ( FVAL - FT ) / H * ENDIF * 1500 CONTINUE * * C-----DO SCALING: DEFINE FV AND SCALE. NOTE THAT IN SOME INSTANCES C THIS MAY REQUIRE AN EXTRA CALL TO GET THE FUNCTION C VALUE WHEN INDIC = DOG; THIS WAS DONE IN THE CALLS C ABOVE. * IF ( FSCALE .NE. 0 ) THEN CALL ZZSCAL( FT, FV, SCALE, FSCALE, FONLY, GONLY ) ELSE FV = FT SCALE = ONE ENDIF * C-----NOW REVISE THE FUNCTION AND GRADIENT AS NECESSARY. * IF ( .NOT. GONLY ) THEN F = FV ENDIF * IF ( .NOT. FONLY .AND. SCALE .NE. ONE ) THEN DO 5100 I=1,N G(I) = G(I) * SCALE 5100 CONTINUE ENDIF * INDIC = OK * GOTO 90000 * * C>>>>>>>>>>>>>>>>>>>>>>>>>> E N T R Y ZZESRT <<<<<<<<<<<<<<<<<<<<<<<<<< * ENTRY ZZESRT ( FSCAL, DDERV, MAXM ) * FCALL = .TRUE. * FSCALE = FSCAL SDRVMD = DDERV MAXFN = MAXM * IFNCT = 0 IGRCT = 0 FTIME = ZERO * IF ( SDRVMD .EQ. TEST .OR. SDRVMD .EQ. TFIRST ) THEN ERR = ZERO SERR = ZERO DCNT = 0 INDEX = 0 GCNT = 0 ENDIF * RETURN * C>>>>>>>>>>>>>>>>>>>>>>>>>> E N T R Y ZZESET <<<<<<<<<<<<<<<<<<<<<<<<<< * ENTRY ZZESET ( DTRF, DTRG, DTRTST, DITRUN ) * TRF = DTRF TRG = DTRG TRTEST = DTRTST ITRUN = DITRUN * RETURN * C>>>>>>>>>>>>>>>>>>>>>>>>>> E N T R Y ZZEDDF <<<<<<<<<<<<<<<<<<<<<<<<<< * ENTRY ZZEDDF ( SANAL, SDIFF, STEST, SFIRST ) * ANAL = SANAL DIFF = SDIFF TEST = STEST TFIRST = SFIRST * RETURN * C>>>>>>>>>>>>>>>>>>>>>>>>>> E N T R Y ZZEFDF <<<<<<<<<<<<<<<<<<<<<<<<<< * ENTRY ZZEFDF ( SDOF, SDOG, SDOFG, SNONE ) * DOF = SDOF DOG = SDOG DOFG = SDOFG NONE = SNONE * RETURN * C>>>>>>>>>>>>>>>>>>>>>>>>>> E N T R Y ZZERDF <<<<<<<<<<<<<<<<<<<<<<<<<< * ENTRY ZZERDF ( SOK, SABORT, SLIMIT, SNOF, SNOG, SNOFG ) * OK = SOK ABORT = SABORT LIMIT = SLIMIT NOF = SNOF NOG = SNOG NOFG = SNOFG * RETURN * C>>>>>>>>>>>>>>>>>>>>>>>>>> E N T R Y ZZEGET <<<<<<<<<<<<<<<<<<<<<<<<<< * ENTRY ZZEGET ( FNCT, GRCT, TIME ) * FNCT = IFNCT GRCT = IGRCT TIME = FTIME * RETURN * C>>>>>>>>>>>>>>>>>>>>>>>>>> E N T R Y ZZECHK <<<<<<<<<<<<<<<<<<<<<<<<<< * ENTRY ZZECHK ( ERROR, AVERR, INDX, ITERAT ) * ERROR = ERR INDX = INDEX ITERAT = GCNT * AVERR = SERR / DCNT * RETURN * C=============================== E X I T =============================== * 90000 CALL ZZSECS (TT) FTIME = FTIME + TT * IF ( TRF .AND. .NOT. BAD() ) WRITE (ITRUN,99998) F * IF ( TRG .AND. .NOT. BAD() ) THEN WRITE (ITRUN,99997) WRITE (ITRUN,99996) G ENDIF * 90500 RETURN * C ALTERNATE RETURN IF MAXIMUM NUMBER OF FUNCTION EVALUATIONS C EXCEEDED. * 91000 IF ( TRF .OR. TRG ) WRITE ( ITRUN,99999 ) * INDIC = LIMIT * RETURN * C============================ F O R M A T S ============================ * 99993 FORMAT( ' EXCESSIVE ERROR IN GRADIENT ESTIMATION.') * 99994 FORMAT( ' ERROR ESTIMATE IN GRADIENT ESTIMATION: ', G15.7/ - ' ESTIMATED FIGURES OF AGREEMENT: ', G9.2 ) * 99995 FORMAT( ' ANALYTIC GRADIENT ', G22.15, ' (COMPONENT ',I3,')'/ - ' ESTIMATED DERIVATIVE', G22.15 ) * 99996 FORMAT( ' ', 5 G15.8 ) * 99997 FORMAT( ' (ZZEVAL) GRADIENT = ' ) * 99998 FORMAT( ' (ZZEVAL) FUNCTION = ', G26.16 ) * 99999 FORMAT(/' THE NUMBER OF FUNCTION EVALUATIONS ALLOWED HAS ', - 'BEEN EXCEEDED.') * C================================ E N D ================================ * END SUBROUTINE ZZPRNT ( N, X, F, G, NRMG, INCR ) * C============== A R G U M E N T D E C L A R A T I O N S ============== * INTEGER N, INCR * DOUBLE PRECISION F, X(N), G(N), NRMG C!!!! REAL F, X(N), G(N), NRMG * C============================= S T A T U S ============================= C C SYSTEM DEPENDENCE: NONE C C CONVERSION FOR SINGLE OR DOUBLE PRECISION: REQUIRED (SEE CONVRT) C C IGNORE LINES BEGINNING WITH "C!!!!" . C C THIS VERSION IS IN D O U B L E PRECISION. C!!!! THIS VERSION IS IN S I N G L E PRECISION. C C REVISION STATUS: DATE AUTHOR VERSION C C MAY 21, 1987 A. BUCKLEY 1.2 C C======================== D E S C R I P T I O N ======================== C C THIS ROUTINE PRINTS (OPTIONALLY) A POINT X , THE VALUE F OF C SOME FUNCTION AT THE POINT X, ALONG WITH THE NORM OF THE GRADIENT C AT THAT POINT, AND (OPTIONALLY) THE VALUE OF THE GRADIENT C G AT THE POINT X. C C TWO OUTPUT UNITS MAY BE SIMULTANEOUSLY DEFINED; EITHER OR BOTH C MAY BE USED. THE PRINT INTERVAL MAY BE DEFINED INDEPENDENTLY FOR C EACH. IF FORCE IS SET, IT APPLIES TO BOTH, IF BOTH ARE DEFINED. C THE DESCRIPTION BELOW ONLY APPLIES TO ONE UNIT. IDENTICAL COMMENTS C APPLY TO THE OTHER. OUTPUT TO A UNIT MAY BE TURNED OFF BY SETTING C THE LEVEL TO 0, OR THE UNIT NUMBER TO 0. ONLY ONE CALL IS NEEDED C AT EACH ITERATE; PRINTING WILL BE DONE ON EITHER OR BOTH UNITS, C AS NEEDED. C C SOME OF THE CONTROL OF PRINT IS THROUGH VARIABLES WHICH ARE C SET THROUGH AN ENTRY POINT CALL TO ZZPSET. THESE ARE DECLARED C AS SAVE VARIABLES. ZZPSET SHOULD BE CALLED TO INITIALIZE ZZPRNT C EACH TIME A FUNCTION IS TO BE MINIMIZED. THE DESCRIPTION OF THE C CONTROL FOLLOWS. C C---DESCRIPTION OF PARAMETERS. C C N THE DIMENSION OF THE PROBLEM. C X THE CURRENT POINT. C F THE FUNCTION VALUE AT X. C G THE GRADIENT VALUE AT X, IF NEEDED. C NRMG THE NORM OF THE GRADIENT. C INCR SEE (5) BELOW. C C---NOTE: 1. CONTROL IS UNDER PLEV1. (PLEV MEANS PRINT LEVEL) C SEE ENTRY POINT ZZP1ST BELOW. C C LET IP = ABS(PLEV1). THEN IF C C PLEV1 = 0 THERE IS NO OUTPUT. C C PLEV1 < 0 PRINT EVERY IP-TH ITERATION: C C THE ITERATION NUMBER IN ITCT, C THE FUNCTION VALUE IN F, C THE NO. OF FUNC EVALUATIONS IN IFNCT. C THE NO. OF GRAD EVALUATIONS IN IGRCT. C C THESE COUNTS ARE OBTAINED THROUGH A CALL C TO THE ENTRY POINT ZZEGET IN ZZEVAL. C C PLEV1 > 0 PRINT EVERY IP-TH ITERATION, AS FOR C PLEV1 < 0, BUT ALSO PRINT: C C THE POINT X, AND C THE GRADIENT G (SEE POINT 2 BELOW). C C 2. SETTING GRAD1 = FALSE WILL ENSURE THAT THE GRADIENTS ARE C NEVER PRINTED, REGARDLESS OF THE VALUE OF PLEV1. THIS C WOULD BE APPROPRIATE WHEN GRADIENTS ARE NOT AVAILABLE OR C TO PRINT X WITHOUT PRINTING G WHEN PLEV1 > 0. THE C SAME COMMENTS APPLY TO SUPPRESSING X WITH POINT1=FALSE. C C 3. PRPT1 RECORDS THE NUMBER OF THE NEXT ITERATION AT WHICH C TO PRINT. WHEN ZZP1ST IS CALLED, THE ITERATION COUNT C ITCT IS INITIALIZED (TO -1) AND PRPT1 IS SET TO 0; THIS C IS WHY ZZP1ST MUST BE CALLED BEFORE EACH FUNCTION IS C MINIMIZED. ON ENTRY TO ZZPRNT, ITCT IS FIRST INCREMENTED C BY 1. THEN, IF ITCT IS LESS THAN PRPT1, NO ACTION TAKES C PLACE AT ALL. IF ON THE OTHER HAND, ITCT = PRPT1, THEN C PRINTING OF THE APPROPRIATE INFORMATION IS DONE, AND C THEN PRPT1 IS ADVANCED BY ABS(PLEV1) TO MARK THE POINT C AT WHICH NEXT TO PRINT. IF ON ENTRY THE VALUE OF ITCT C IS BEYOND THAT OF PRPT1, PRPT1 IS REPEATEDLY INCREMENTED BY C ABS(PLEV1) UNTIL ONE OF THE FIRST TWO CASES OCCURS. OF C COURSE THIS IS NOT SUPPOSED TO HAPPEN IF ZZP1ST WAS C CALLED TO INITIALIZE ZZPRNT. C C 4. ITCT IS INCREMENTED INTERNALLY, BUT THE COUNTING OF C IFNCT AND IGRCT WILL BE DONE BY ZZEVAL; THE VALUES ARE C OBTAINED BY A CALL TO THE ENTRY POINT ZZEGET. C C 5. SETTING INCR CAN BE USED TO FORCE PRINTING, REGARDLESS C OF THE VALUES OF PRPT1 AND ITCT. THIS IS USEFUL FOR C FORCING PRINTING OF THE FINAL POINT REACHED. IN FACT, C INCR DEFINES THE AMOUNT BY WHICH TO INCREMENT THE INTERNAL C ITERATION COUNTER OF ZZPRNT. THUS, NORMALLY ZZPRNT WILL C BE CALLED WITH INCR = 1. TO FORCE PRINTING, ZZPRNT MAY BE C CALLED WITH INCR = 0; THE POINT IS PRINTED BUT THE ITERA- C TION COUNTER IS NOT ADVANCED. FINALLY, IF ONE WISHES TO C INSIST THAT THE ITERATION COUNTER BE UPDATED CORRECTLY AND C THAT THE POINT BE PRINTED REGARDLESS OF THE VALUE OF PRPT1, C ONE MAY CALL ZZPRNT WITH INCR = -1; BECAUSE INCR <=0, THE C PRINT OF THE POINT WILL BE FORCED, AND IT IS IN FACT C ABS(INCR) THAT IS USED TO UPDATE THE ITERATION COUNT. THIS C PARTICULAR CASE IS USEFUL AT THE FINAL POINT. C C NOTE THAT PRPT1 IS STILL ADVANCED, BUT ONLY IF APPROPRIATE, C I.E. IF PRINTING WOULD HAVE BEEN DONE ANYWAY, AS EXPLAINED C IN 2. C C ALSO, WHEN FORCING IS DONE, THE ROUTINE IS CAREFUL NOT C TO REPEAT A PRINTING REQUEST. IF THE OUTPUT UNIT OR C THE STATUS OF GRAD1 OR THE ITERATION COUNT ITCT IS C DIFFERENT, THEN THE PRINTING IS DONE; OTHERWISE IT C IS CONSIDERED A REPEAT OF A PREVIOUS REQUEST AND IT C IS IGNORED. C C 6. PRTIME IS USED FOR ACCUMULATING THE TIME SPENT IN THE C PRINT ROUTINE. IT IS INITIALIZED TO ZERO AT THE CALL C TO ZZP1ST, AND EACH CALL TO ZZPRNT INCREMENTS PRTIME BY C THE AMOUNT OF TIME SPENT IN THE ROUTINE. C C---ALSO AVAILABLE THROUGH ZZPGET (TIME, ITER). C C THE USER MAY CALL ZZPGET AT ANY TIME TO GET THE AMOUNT OF TIME C SPENT IN THE PRINT ROUTINE AND THE CURRENT ITERATION COUNT. C THESE ARE, RESP., THE ARGUMENTS ITER AND TIME. C C======================= E N T R Y P O I N T S ======================= C C ZZPRNT ...THE NATURAL ENTRY. C ZZP1ST ...TO INITIALIZE CONTROL VARIABLES FOR FIRST UNIT. C ZZP2ST ...TO INITIALIZE CONTROL VARIABLES FOR SECOND UNIT. C ZZPGET ...TO RETURN ITERATION COUNT AND TIME. C C======================== S U B R O U T I N E S ======================== C C ABS ...INTRINSIC FUNCTION. C C ZZSECS ...FOR PRINT TIMING. C ZZEGET ...ENTRY TO ZZEVAL FOR FUNCTION/GRADIENT COUNTS. C C========================= P A R A M E T E R S ========================= * DOUBLE PRECISION ZERO, ONE, TWO, THREE C!!!! REAL ZERO, ONE, TWO, THREE PARAMETER ( ZERO = 0D0, ONE = 1D0, TWO = 2D0, THREE = 3D0) * DOUBLE PRECISION FOUR, FIVE, SIX, SEVEN C!!!! REAL FOUR, FIVE, SIX, SEVEN PARAMETER ( FOUR = 4D0, FIVE = 5D0,SIX = 6D0, SEVEN = 7D0) * DOUBLE PRECISION EIGHT, NINE, TEN C!!!! REAL EIGHT, NINE, TEN PARAMETER ( EIGHT = 8D0, NINE = 8D0, TEN = 10D0 ) * C================= L O C A L D E C L A R A T I O N S ================= * INTEGER IFNCT, IGRCT, ITCT INTEGER PRPT1, PLEV1, UNIT1, LSTIT1, LSTUN1 INTEGER PRPT2, PLEV2, UNIT2, LSTIT2, LSTUN2 * LOGICAL GRAD1, LSTGR1, POINT1, LSTPT1 LOGICAL GRAD2, LSTGR2, POINT2, LSTPT2 LOGICAL FORCE, GOT * DOUBLE PRECISION SECS, PRTIME, DTIME C!!!! REAL SECS, PRTIME, DTIME * C-----DECLARATIONS FOR ENTRY POINT DUMMY ARGUMENTS. * LOGICAL DGRAD1, DPINT1 LOGICAL DGRAD2, DPINT2 * INTEGER DPRUN1, DPRNT1 INTEGER DPRUN2, DPRNT2 * INTEGER ITER * DOUBLE PRECISION TIME C!!!! REAL TIME * C=============================== S A V E =============================== * SAVE PRTIME, ITCT * SAVE PRPT1, GRAD1, PLEV1, UNIT1, POINT1 SAVE PRPT2, GRAD2, PLEV2, UNIT2, POINT2 SAVE LSTIT1, LSTUN1, LSTGR1, LSTPT1 SAVE LSTIT2, LSTUN2, LSTGR2, LSTPT2 * C======================= E Q U I V A L E N C E S ======================= C C THERE ARE NO EQUIVALENCES. C C============================= C O M M O N ============================= C C THERE ARE NO COMMON BLOCKS. C C=============================== D A T A =============================== * DATA PLEV1/0/, PLEV2/0/, UNIT1/6/, UNIT2/0/ * C========================== E X E C U T I O N ========================== * FORCE = INCR .LE. 0 * ITCT = ITCT + ABS(INCR) * CALL ZZSECS (SECS) PRTIME = PRTIME - SECS * GOT = .FALSE. * IF ( FORCE .AND. ( ITCT .EQ. LSTIT1 ) - .AND. ( UNIT1 .EQ. LSTUN1 ) - .AND. ( POINT1 .EQV. LSTPT1 ) - .AND. ( GRAD1 .EQV. LSTGR1 ) ) THEN C DON'T REPEAT AN EARLIER REQUEST. GOTO 2000 * ENDIF * 100 IF ( PLEV1 .NE. 0 .AND. ITCT .GT. PRPT1 ) THEN * PRPT1 = PRPT1 + ABS(PLEV1) GOTO 100 * ENDIF * IF ( (UNIT1 .NE. 0 ) .AND. - (PLEV1 .NE. 0 ) .AND. - (FORCE .OR. (ITCT .EQ. PRPT1)) ) THEN * C -----SAVE INFORMATION DEFINING THIS PRINT REQUEST. * LSTIT1 = ITCT LSTUN1 = UNIT1 LSTGR1 = GRAD1 LSTPT1 = POINT1 * C ------PRINT ITERATION NUMBER, FUNCTION VALUE, NORM OF G, AND C NUMBER OF FUNCTION/GRADIENT EVALUATIONS. * CALL ZZEGET ( IFNCT, IGRCT, DTIME ) GOT = .TRUE. * WRITE ( UNIT1, 99999 ) ITCT,F,IFNCT,NRMG,IGRCT,DTIME * C ------IF PLEV1 > 0 , ALSO PRINT X AND G. * IF ( PLEV1 .GT. 0 ) THEN * IF ( POINT1) THEN WRITE (UNIT1,99998) X ENDIF * IF ( GRAD1 ) THEN WRITE (UNIT1,99997) G ENDIF * ENDIF * C ------UPDATE COUNTER. * IF (ITCT .EQ. PRPT1) PRPT1 = PRPT1 + ABS(PLEV1) * ENDIF * 2000 IF ( FORCE .AND. ( ITCT .EQ. LSTIT2 ) - .AND. ( UNIT2 .EQ. LSTUN2 ) - .AND. ( POINT2 .EQV. LSTPT2 ) - .AND. ( GRAD2 .EQV. LSTGR2 ) ) THEN C DON'T REPEAT AN EARLIER REQUEST. GOTO 4000 * ENDIF * 2200 IF ( PLEV2 .NE. 0 .AND. ITCT .GT. PRPT2 ) THEN * PRPT2 = PRPT2 + ABS(PLEV2) GOTO 2200 * ENDIF * IF ( (UNIT2 .NE. 0 ) .AND. - (PLEV2 .NE. 0 ) .AND. - (FORCE .OR. (ITCT .EQ. PRPT2)) ) THEN * C -----SAVE INFORMATION DEFINING THIS PRINT REQUEST. * LSTIT2 = ITCT LSTUN2 = UNIT2 LSTGR2 = GRAD2 LSTPT2 = POINT2 * C ------PRINT ITERATION NUMBER, FUNCTION VALUE, NORM OF G, AND C NUMBER OF FUNCTION/GRADIENT EVALUATIONS. * IF ( .NOT. GOT ) CALL ZZEGET ( IFNCT, IGRCT, DTIME ) * WRITE ( UNIT2, 99989 ) ITCT,F,IFNCT,NRMG * C ------IF PLEV2 > 0 , ALSO PRINT X AND G. * IF ( PLEV2 .GT. 0 ) THEN * IF ( POINT2) THEN WRITE (UNIT2,99988) X ENDIF * IF ( GRAD2 ) THEN WRITE (UNIT2,99987) G ENDIF * ENDIF * C ------UPDATE COUNTER. * IF (ITCT .EQ. PRPT2) PRPT2 = PRPT2 + ABS(PLEV2) * ENDIF * 4000 GOTO 90000 * C>>>>>>>>>>>>>>>>>>>>>>>>>> E N T R Y ZZP1ST <<<<<<<<<<<<<<<<<<<<<<<<<< * ENTRY ZZP1ST ( DPRUN1, DGRAD1, DPINT1, DPRNT1 ) * UNIT1 = DPRUN1 GRAD1 = DGRAD1 POINT1 = DPINT1 PLEV1 = DPRNT1 * PRPT1 = 0 ITCT = -1 PRTIME = ZERO * RETURN * C>>>>>>>>>>>>>>>>>>>>>>>>>> E N T R Y ZZP2ST <<<<<<<<<<<<<<<<<<<<<<<<<< * ENTRY ZZP2ST ( DPRUN2, DGRAD2, DPINT2, DPRNT2 ) * UNIT2 = DPRUN2 GRAD2 = DGRAD2 POINT2 = DPINT2 PLEV2 = DPRNT2 * PRPT2 = 0 ITCT = -1 PRTIME = ZERO * RETURN * C>>>>>>>>>>>>>>>>>>>>>>>>>> E N T R Y ZZPGET <<<<<<<<<<<<<<<<<<<<<<<<<< * ENTRY ZZPGET ( TIME, ITER ) * TIME = PRTIME ITER = ITCT * RETURN * C=============================== E X I T =============================== * 90000 CALL ZZSECS (SECS) PRTIME = PRTIME + SECS * 91000 RETURN * C============================ F O R M A T S ============================ * 99987 FORMAT(' GRAD: ', 7G9.2 / (1X,8G9.2) ) * 99988 FORMAT(' POINT X:', 7G9.2 / (1X,8G9.2) ) * 99989 FORMAT(' PT #',I3,'; F=',G15.8,'(#',I3,') !!G!!=',E7.2) * 99997 FORMAT(' THE GRADIENT AT THIS POINT IS ', 3G15.8 / (1X,5G15.8) ) * 99998 FORMAT(' THE VARIABLES HAVE THE CURRENT VALUES GIVEN BY ',4X, - G26.16 / (2X,3G26.16) ) * 99999 FORMAT(' ',' ...PT ',I3,'; F=',G23.16,'(#',I3,') !!G!!=', - E7.2, '(#',I3,'); ',F8.3,' SECS' ) * C================================ E N D ================================ * END SUBROUTINE ZZSCAL ( FT, FV, SCALE, FSCALE, FONLY, GONLY ) * C============== A R G U M E N T D E C L A R A T I O N S ============== * INTEGER FSCALE * LOGICAL FONLY, GONLY * DOUBLE PRECISION FT, FV, SCALE C!!!! REAL FT, FV, SCALE * C============================= S T A T U S ============================= C C SYSTEM DEPENDENCE: NONE C C CONVERSION FOR SINGLE OR DOUBLE PRECISION: REQUIRED (SEE CONVRT) C C IGNORE LINES BEGINNING WITH "C!!!!" . C C THIS VERSION IS IN D O U B L E PRECISION. C!!!! THIS VERSION IS IN S I N G L E PRECISION. C C REVISION STATUS: DATE AUTHOR VERSION C C DEC. 15, 1986 A. BUCKLEY 1.0 C C======================== D E S C R I P T I O N ======================== C C THIS SUBROUTINE APPLIES ONE OF SEVERAL SCALINGS (LINEAR OR C NONLINEAR) TO A FUNCTION VALUE. C C-----ON ENTRY: C C FT - THE PRESENT FUNCTION VALUE C C FSCALE - THE CODE FOR THE TYPE OF SCALE DESIRED. WHERE C THE SCALE FUNCTION IS ONE OF THE FOLLOWING: C C 1: F(Z) = 1 + Z C 2: F(Z) = Z*Z C 3: F(Z) = -1 / (1 + Z*Z) C 4: F(Z) = SQRT(1 + Z*Z) C 5: F(Z) = Z*Z*Z C C FONLY - IF TRUE ONLY THE FUNCTION IS EVALUATED. C C GONLY - IF TRUE ONLY THE GRADIENT IS EVALUATED. C C-----ON EXIT: C C FV - THE SCALED FUNCTION VALUE. C C SCALE - GRADIENT SCALING FACTOR. C C======================= E N T R Y P O I N T S ======================= C C NONE ARE USED. C C======================== S U B R O U T I N E S ======================== C C SQRT... INTRINSIC C C========================= P A R A M E T E R S ========================= * DOUBLE PRECISION ZERO, ONE, TWO, THREE C!!!! REAL ZERO, ONE, TWO, THREE PARAMETER ( ZERO = 0D0, ONE = 1D0, TWO = 2D0, THREE = 3D0) * DOUBLE PRECISION FOUR, FIVE, SIX, SEVEN C!!!! REAL FOUR, FIVE, SIX, SEVEN PARAMETER ( FOUR = 4D0, FIVE = 5D0,SIX = 6D0, SEVEN = 7D0) * DOUBLE PRECISION EIGHT, NINE, TEN C!!!! REAL EIGHT, NINE, TEN PARAMETER ( EIGHT = 8D0, NINE = 8D0, TEN = 10D0 ) * C================= L O C A L D E C L A R A T I O N S ================= C C THERE ARE NO LOCAL DECLARATIONS. C C=============================== S A V E =============================== C C THERE ARE NO VARIABLES TO BE SAVED. C C======================= E Q U I V A L E N C E S ======================= C C THERE ARE NO EQUIVALENCES. C C============================= C O M M O N ============================= C C THERE ARE NO COMMON BLOCKS. C C=============================== D A T A =============================== C C NO DATA STATEMENTS ARE USED. C C========================== E X E C U T I O N ========================== * GOTO (2100,2200,2300,2400,2500), FSCALE * C -----FF(Z) = 1 + F(Z) -------FSCALE = 1. * 2100 IF ( .NOT. GONLY ) FV = FT + ONE IF ( .NOT. FONLY ) SCALE = ONE GOTO 90000 * C -----FF(Z) = Z*Z ------------FSCALE = 2. * 2200 IF ( .NOT. GONLY ) FV = FT * FT IF ( .NOT. FONLY ) SCALE = TWO * FT GOTO 90000 * C -----FF(Z) = -1/(1+Z**2) --- FSCALE = 3. * 2300 FV = -ONE / ( ONE + FT**2 ) IF ( .NOT. FONLY ) SCALE = TWO * FT * FV**2 GOTO 90000 * C -----FF(Z) = SQRT(1+Z**2) -- FSCALE = 4. * 2400 FV = SQRT(ONE + FT**2) IF ( .NOT. FONLY ) SCALE = FT/FV GOTO 90000 * C -----FF(Z) = Z*Z*Z --------- FSCALE = 5. * 2500 IF ( .NOT. GONLY ) FV = FT*FT*FT IF ( .NOT. FONLY ) SCALE = THREE*FT*FT GOTO 90000 * C=============================== E X I T =============================== * 90000 RETURN * C================================ E N D ================================ * END SUBROUTINE ZZTER