PROGRAM TML ************************************************************************ * TML v1.0, Mar 2003 * * This code implements the method: * * "Topographical Multi Level Single Linkage" * * due to A.H.G. Rinnooy Kan and G.T. Timmer * described in the following two articles: * * 1) A.H.G. Rinnooy Kan and G.T.Timmer, * Stochastic global optimization methods * Part I: Clustering methods * Mathematical Programming, Vol 39, (1987) pp 27-56. * * 2) A.H.G. Rinnooy Kan and G.T.Timmer, * Stochastic global optimization methods * Part II: Multi level methods * Mathematical Programming, Vol 39, (1987) pp 57-78. * * with modifications suggested by Aimo Torn and Sami Viitanen * described in: * * A. Torn and S. Viitanen, * Topographical global optimization using pre-sampled points, * Journal of Global Optimization, Vol 5, (1994) pp 267-276. * * and modifications due to M.M. Ali and C. Storey * described in: * * M.M. Ali and C. Storey, * Topographical multilevel single linkage, * Journal of Global Optimization, Vol 5, (1994) pp 349-358. * * Newer versions of this code are available from: * http://nrt.cs.uoi.gr/merlin ************************************************************************ IMPLICIT DOUBLE PRECISION (A-H,O-Z) C --------------------------------------------------------------------- C Maximum number of minima. PARAMETER ( MXNM = 10000 ) C Maximum number of parameters. PARAMETER ( MXV = 200 ) C Maximum number of sample points. PARAMETER ( MXNS = 500 ) C Maximum number of nearest neighbors. PARAMETER ( MXNNN = 8 ) C LP=.TRUE. enables informative printout. LP=.FALSE. disables it. LOGICAL LP PARAMETER ( LP = .FALSE. ) C --------------------------------------------------------------------- C Input file. CHARACTER*(*) INFIL PARAMETER ( INFIL = 'TML.DAT' ) CHARACTER*(*) CMDFIL PARAMETER ( CMDFIL = 'command' ) C --------------------------------------------------------------------- C Workspaces: C Maximum total points (minima + sample). PARAMETER ( MXNT = MXNM+MXNS ) C Total set of points, set of minima. DIMENSION S(MXV,MXNT), XMINZ(MXV,MXNM) DIMENSION SS(MXV) C Corresponding values of the objective function for the above sets. DIMENSION FVAL(MXNT), FMINZ(MXNM) C Indices (pointers) for the Graph-Minima and for the local starters. DIMENSION IGMIN(MXNT), ISTART(MXNT) C Indices for the nearest neighbors, Euclidean distances. DIMENSION IREGIO(MXNNN), DIS(MXNT) DIMENSION POINT(MXV), XLL(MXV), XRL(MXV), IXAT(MXV), GRAD(MXV) DIMENSION ICODE(4) COMMON /MODE3/ MODEV, IFFORM, INMCL C --------------------------------------------------------------------- C Local variables: CHARACTER*20 FINP, FOUT, FMIN, POIFIL, INDIR CHARACTER*10 CIT, CEST LOGICAL EXI1, EXI2, MEXI, CEXI, ERR DATA ICODE / 1, 1, 1, 1 / C --------------------------------------------------------------------- C Default values for the input parameters. DATA NOD / 2 / DATA NOF / 0 / DATA ITTHR / 10 / DATA HEAL / 5 / DATA NEARN / 2 / DATA NSAMPL / 100 / DATA SIGMA / 4. / DATA COVTHR / 0.995 / DATA XDEPS / 1.D-2 / DATA FDEPS / 1.D-4 / DATA GDEPS / 1.D-3 / DATA IREST / 0 / DATA NDUMP / 100 / DATA FINP / 'in.dat' / DATA FOUT / 'out.dat' / DATA FMIN / 'MINIMA' / DATA POIFIL / 'POIMAR.DAT' / DATA INDIR / 'input' / C Default values for bounds and fix statuses. DATA XLL / MXV * -1.D0 / DATA XRL / MXV * 1.D0 / DATA IXAT / MXV * 1 / C --------------------------------------------------------------------- C Read input data. CALL FROMFL(INFIL,NOD,NOF,ITTHR,HEAL,NEARN,NSAMPL, & SIGMA,COVTHR,XDEPS,FDEPS,GDEPS,IREST,NDUMP, & FINP,FOUT,FMIN,POIFIL,INDIR,POINT,XLL,XRL, & IXAT,MXV,MXNNN,MXNS) TNOS = 0 ITC = 0 NOM = 0 C Initialize Function, Gradient, Hessian and Jacobian call-counters. NFTOT = 0 NGTOT = 0 NHTOT = 0 NJTOT = 0 MXSP = 0 C Number of local optimizations. NLOC = 0 C C If restore was specified, read the dump file. IF (IREST.NE.0) & CALL RESTOR(NOD,NOF,ITTHR,HEAL,NEARN,NSAMPL,SIGMA,COVTHR, & XDEPS,FDEPS,GDEPS,NDUMP,FINP,FOUT,FMIN,POIFIL, & INDIR,ITC,NOM,NFTOT,NGTOT,NHTOT,NJTOT,TNOS,XMINZ, & FMINZ,XLL,XRL,IXAT,MXV,MXNM,MXSP,NLOC) C C Informative printout. WRITE (*,*) WRITE (*,*) 'TML running with ...' WRITE (*,*) 'Iteration threshold: ', ITTHR WRITE (*,*) 'Healing parameter: ', HEAL WRITE (*,*) 'Number of nearest neighbors: ', NEARN WRITE (*,*) 'Sample size: ', NSAMPL WRITE (*,*) 'Sigma parameter: ', SIGMA WRITE (*,*) 'Coverage threshold: ', COVTHR WRITE (*,*) 'X distance criterion: ', XDEPS WRITE (*,*) 'F distance criterion: ', FDEPS WRITE (*,*) 'Gradient convergence criterion: ', GDEPS WRITE (*,*) 'Saving every ', NDUMP, ' iterations' WRITE (*,*) 'Merlin input file: ', FINP WRITE (*,*) 'Merlin output file: ', FOUT WRITE (*,*) C NOE = NSAMPL NNB = NEARN C C Set the appropriate Merlin flags for function evaluation. MODEV = 1 IF (NOF.EQ.0) THEN IFFORM = 0 ELSE IFFORM = 1 END IF C NOFV = 0 ALVOL = 0. BOXDI = 0. DO 1,I=1,NOD IF (IXAT(I) .EQ. 1) THEN NOFV = NOFV + 1 ALVOL = ALVOL + LOG(XRL(I)-XLL(I)) BOXDI = BOXDI + (XRL(I)-XLL(I))**2 END IF 1 CONTINUE BOXDI = 2*SQRT(BOXDI) C C /\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\ C Main loop starts below the line: ( 30 CONTINUE ) C /\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\ C 30 CONTINUE ITC = ITC + 1 C C Calculate the critical radious RC. CALL CRITIC(ITC,NOE,NOFV,ALVOL,SIGMA,RC) IF (LP) WRITE (*,*) 'Calculated critical radious: ', RC C C Sample NOE points in array S(,). DO 2,I=1,NOE DO 3,J=1,NOD SS(J) = RANM()*(XRL(J)-XLL(J))+XLL(J) IF (IXAT(J).EQ.0) SS(J) = POINT(J) S(J,I) = SS(J) 3 CONTINUE FVAL(I) = ACSQ(SS,NOD) NFTOT = NFTOT+1 2 CONTINUE C C Add the NOM minima to the set S(,) DO 4,I=1,NOM NI = NOE+I DO 5,J=1,NOD S(J,NI) = XMINZ(J,I) 5 CONTINUE FVAL(NI) = FMINZ(I) 4 CONTINUE NOV = NOE + NOM C C Calculate the Graph minima, with respect to NNB nearest neighbours. C IGMIN[I] = 0 means that the I_th point in S, is not a graph minimum C IGMIN[I] = 1 means that the I_th point in S, is a graph minimum C CALL GRAPHM(S,NOD,NOV,NNB,FVAL,IGMIN,IREGIO,NOGM,DIS,MXV,BOXDI) IF (LP) WRITE (*,*) 'Calculated graph minima ', NOGM C C Determine Start points. C ISTART(I) is a Flag-Array. C If the I_th point in S is a Start Point then ISTART(I) = 1 C Otherwise ISTART(I) is set to 0 (ZERO). C Start Points are chosen among the Graph minima. C A previously found minimum cannot be a Start Point. C NOSP = 0 DO 6,I = 1,NOE ISTART(I) = IGMIN(I) C Consider only points that correspond to graph minima. IF (ISTART(I) .EQ. 0) GOTO 70 C Calculate the distances dis(j) = |r(i)-r(j)| for j=1,NOV C (r(i) = S(.,i)) CALL DIST(S,NOD,NOV,I,DIS,MXV) C C Check if there is a point within the critical radious RC, C that corresponds to a lower function value. C In that case r(i) cannot be a start point. C DO 7,J=1,NOV IF (DIS(J).LT.RC) THEN IF (FVAL(J) .LT. FVAL(I)) THEN ISTART(I) = 0 GOTO 70 END IF END IF 7 CONTINUE 70 CONTINUE NOSP = NOSP + ISTART(I) 6 CONTINUE C Start points are determined. (A total of NOSP start points are found) IF (LP) WRITE (*,*) 'Start points: ', NOSP C EXX = EXP(-ITC/HEAL) TNOS = TNOS + NOSP + (NOGM-NOSP)*(1-EXX)/(1+EXX) C DO 200,I = 1, NOE IF (NOM .EQ. MXNM) THEN WRITE (*,*) 'TML: MAX number of allowed MINIMA reached.' STOP END IF IF (ISTART(I) .EQ. 0) GOTO 200 DO 11 J=1,NOD POINT(J) = S(J,I) 11 CONTINUE CALL OPTIMA(NOD,NOF,POINT,VAL,XLL,XRL,IXAT,ICODE, & FINP,FOUT,GRMS,NF,NG,NH,NJ) IF (LP) WRITE (*,*) 'Completed local search.' NLOC = NLOC+1 NFTOT = NFTOT + NF NGTOT = NGTOT + NG NHTOT = NHTOT + NH NJTOT = NJTOT + NJ CALL ADDMIN(LP,XDEPS,FDEPS,GDEPS,POINT,NOD,VAL,XLL,XRL, & DIS,GRAD,NOM,FMIN,XMINZ,FMINZ,MXV,MXNM,MXNT, & NFTOT,NGTOT) 200 CONTINUE MXSP = MAX(MXSP,NOSP) C ESTNOM = NOM*(TNOS-1)/(TNOS-NOM-2) IF (ESTNOM.LT.0.) ESTNOM = NOM+1 COVER = (TNOS-NOM-1)*(TNOS+NOM)/TNOS/(TNOS-1) IF (COVER.LT.0.) COVER = 0. C C Informative printout. CALL I2STR(ITC,CIT,LECIT) WRITE (CEST,'(F10.1)') ESTNOM DO 40,IEST=1,10 IF (CEST(IEST:IEST).NE.' ') GOTO 41 40 CONTINUE 41 CONTINUE WRITE (*,130) CIT(1:LECIT), COVER, CEST(IEST:) C C Stop the program if the command file exists. INQUIRE (FILE=CMDFIL,EXIST=CEXI) IF (CEXI) THEN CALL DUMPER(NOD,NOF,ITTHR,HEAL,NEARN,NSAMPL,SIGMA, & COVTHR,XDEPS,FDEPS,GDEPS,NDUMP,FINP,FOUT,FMIN, & POIFIL,INDIR,ITC,NOM,NFTOT,NGTOT,NHTOT,NJTOT, & TNOS,XMINZ,FMINZ,XLL,XRL,IXAT,MXV,MXNM, & MXSP,NLOC) OPEN (10,FILE=CMDFIL) CLOSE (10,STATUS='DELETE') STOP END IF C C Save the current state every NDUMP iterations. IF (NDUMP.NE.0) THEN IF (MOD(ITC,NDUMP).EQ.0) & CALL DUMPER(NOD,NOF,ITTHR,HEAL,NEARN,NSAMPL,SIGMA, & COVTHR,XDEPS,FDEPS,GDEPS,NDUMP,FINP,FOUT, & FMIN,POIFIL,INDIR,ITC,NOM,NFTOT,NGTOT, & NHTOT,NJTOT,TNOS,XMINZ,FMINZ,XLL,XRL,IXAT, & MXV,MXNM,MXSP,NLOC) END IF C IF (ITC .LT. ITTHR) GOTO 30 IF ( TNOS .LE. NOM + 2) GOTO 30 IF ( ESTNOM .LT. NOM + 0.5 .AND. COVER .GT. COVTHR ) THEN WRITE (*,*) WRITE (*,*) 'TML run completed.' IF (NOM.EQ.0) THEN WRITE (*,*) 'No local minima were found.' ELSE CLOSE (1) WRITE (*,*) 'Number of local minimizers found: ', NOM WRITE (*,*) 'The minimizers are disposed to file: ', FMIN END IF WRITE (*,*) 'Total number of Function calls: ', NFTOT WRITE (*,*) 'Total number of Gradient calls: ', NGTOT WRITE (*,*) 'Total number of Jacobian calls: ', NJTOT WRITE (*,*) 'Total number of Hessian calls : ', NHTOT WRITE (*,*) 'Total number of local optimizations : ', NLOC WRITE (*,*) 'Maximum number of starting points: ', MXSP CALL DUMPER(NOD,NOF,ITTHR,HEAL,NEARN,NSAMPL,SIGMA, & COVTHR,XDEPS,FDEPS,GDEPS,NDUMP,FINP,FOUT,FMIN, & POIFIL,INDIR,ITC,NOM,NFTOT,NGTOT,NHTOT,NJTOT, & TNOS,XMINZ,FMINZ,XLL,XRL,IXAT,MXV,MXNM, & MXSP,NLOC) STOP END IF GOTO 30 C /\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\ C End of main loop. C /\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\ 130 FORMAT ('Iterations: ',A,5X,'Coverage: ',F7.5,5X, & 'Est. minima: ',A) C END C --------------------------------------------------------------------- SUBROUTINE CHEMIN ( GN, EPS, POINT, NOD, VAL, XLL, XRL, GRAD, & NOC, NGR ) C --------------------------------------------------------------------- C C Description: C Calculates the root mean square gradient, for the C components that do not "coincide" with one of their margins. C If all variables coincide with one of their margins then C the rms Gradient (GN) is set to -1. C C --------------------------------------------------------------------- IMPLICIT DOUBLE PRECISION (A-H,O-Z) C --------------------------------------------------------------------- C Arguments: DIMENSION POINT(NOD), GRAD(NOD), XLL(NOD), XRL(NOD) C --------------------------------------------------------------------- C CALL DERMER(POINT,VAL,GRAD,NOC,NGR,1) GN = 0. DO 10,I = 1,NOD RGRAD = GRAD(I) YL = ABS( XLL(I) - POINT(I) ) YR = ABS( XRL(I) - POINT(I) ) YLX = MAX( ABS( XLL(I) ), ABS( POINT(I) ) )*EPS YRX = MAX( ABS( XRL(I) ), ABS( POINT(I) ) )*EPS IF (YL .LT. YLX) RGRAD = MIN(0.D0,GRAD(I)) IF (YR .LT. YRX) RGRAD = MAX(0.D0,GRAD(I)) GN = GN + RGRAD**2 10 CONTINUE GN = SQRT(GN/NOD) END C --------------------------------------------------------------------- SUBROUTINE CRITIC ( IT, NOE, NOFV, ALVOL, SIGMA, RC ) C --------------------------------------------------------------------- C C Description: C Calculates the critical distance RC for the MLSL algorithm. C Uses logarithms to avoid intermediate overflows. C C --------------------------------------------------------------------- IMPLICIT DOUBLE PRECISION (A-H,O-Z) C --------------------------------------------------------------------- C PI = ACOS(-1.D0) C C The following IF structure, calculates the C Logarithm of the Gamma function: C algam = Gamma(1+nofv/2) C for two cases: NOFV is even, and NOFV is odd. C IF ( (NOFV/2)*2 .EQ. NOFV ) THEN ALGAM = 0. DO 10,I=2,NOFV/2 ALGAM = ALGAM + LOG(DBLE(I)) 10 CONTINUE ELSE ALGAM = LOG(SQRT(PI)) - (NOFV+1)/2*LOG(2.D0) DO 20,I=1,NOFV,2 ALGAM = ALGAM + LOG(DBLE(I)) 20 CONTINUE END IF C DIN = DBLE(IT*NOE) ALRC = -LOG(SQRT(PI))+LOG( SIGMA*LOG(DIN)/DIN )/NOFV ALRC = ALRC + ALVOL/NOFV + ALGAM/NOFV RC = EXP(ALRC) END C --------------------------------------------------------------------- SUBROUTINE MINEL ( DIS, NOE, JN, VMIN ) C --------------------------------------------------------------------- C C Description: C Given an array (DIS) with NOE elements, subroutine MINEL C calculates the minimum element (VMIN) as well as its index (JN). C C --------------------------------------------------------------------- IMPLICIT DOUBLE PRECISION (A-H,O-Z) C --------------------------------------------------------------------- C Arguments: DIMENSION DIS(*) C --------------------------------------------------------------------- C JN = 1 VMIN = DIS(JN) DO 10,I=2,NOE IF (DIS(I) .LT. VMIN) THEN VMIN = DIS(I) JN = I END IF 10 CONTINUE END C --------------------------------------------------------------------- SUBROUTINE MAXEL ( DIS, NOE, JX, VMAX ) C --------------------------------------------------------------------- C C Description: C Given an array (DIS) with NOE elements, subroutine MAXEL C calculates the maximum element (VMAX) as well as its index (JX). C C --------------------------------------------------------------------- IMPLICIT DOUBLE PRECISION (A-H,O-Z) C --------------------------------------------------------------------- C Arguments: DIMENSION DIS(*) C --------------------------------------------------------------------- C JX = 1 VMAX = DIS(JX) DO 10,I=2,NOE IF (DIS(I) .GT. VMAX) THEN VMAX = DIS(I) JX = I END IF 10 CONTINUE END C --------------------------------------------------------------------- SUBROUTINE GRAPHM ( S, NOD, NOV, NNB, FVAL, IGMIN, IREGIO, NOGM, & DIS, MXV, BOXDI ) C --------------------------------------------------------------------- C C Description: C Calculates the graph minima, for NNB number of nearest neighbours. C The graph is the set of points contained in array S(MXV,*) C It returns the array IGMIN(*) C IGMIN(I) = 0 means that the I_th point in S, is not a graph minimum C IGMIN(I) = 1 means that the I_th point in S, is a graph minimum C NOGM is output. Corresponds to the number of graph minima C FVAL(I) is the value corresponding to the I_th point in S C C --------------------------------------------------------------------- IMPLICIT DOUBLE PRECISION (A-H,O-Z) C --------------------------------------------------------------------- C Arguments: DIMENSION S(MXV,*), IGMIN(*), IREGIO(*), FVAL(*), DIS(*) C --------------------------------------------------------------------- C NOGM = 0 DO 10,I=1,NOV CALL DIST(S,NOD,NOV,I,DIS,MXV) CALL NEARNE(DIS,NOV,I,IREGIO,NNB,BOXDI) DO 20,J=1,NNB DIS(J) = FVAL( IREGIO(J) ) 20 CONTINUE CALL MINEL(DIS,NNB,JM,FM) IF ( FVAL(I) .LE. FM ) THEN IGMIN(I) = 1 NOGM = NOGM + 1 ELSE IGMIN(I) = 0 END IF 10 CONTINUE END C --------------------------------------------------------------------- SUBROUTINE DIST ( S, NOD, NOV, I, DIS, MXV ) C --------------------------------------------------------------------- C C Description: C Given: C 1) An array S containing NOV vectors, with NOD components each, C 2) An Integer I, in [1,NOV], C It calculates the Euclidean distances C between the Ith and the Jth vector in S for all J=1,2,...,NOV C and returns them in array DIS C C --------------------------------------------------------------------- IMPLICIT DOUBLE PRECISION (A-H,O-Z) C --------------------------------------------------------------------- C Arguments: DIMENSION S(MXV,*), DIS(*) C --------------------------------------------------------------------- C DO 10,J=1,NOV D = 0.D0 DO 20,K=1,NOD D = D + (S(K,J)-S(K,I))**2 20 CONTINUE DIS(J) = SQRT(D) 10 CONTINUE END C --------------------------------------------------------------------- SUBROUTINE NEARNE ( DIS, NOV, I, IREGIO, NNB, BOXDI ) C --------------------------------------------------------------------- C C Description: C Calculates the NNB nearest neighbours of point with index I C and returns their indices in IREGIO(1:NNB) C DIS(1:NOV) contains the distances |X(I)-X(J)| for J=1,2,...,NOV C The array DIS is destroyed. C C --------------------------------------------------------------------- IMPLICIT DOUBLE PRECISION (A-H,O-Z) C --------------------------------------------------------------------- C Arguments: DIMENSION DIS(*), IREGIO(*) C --------------------------------------------------------------------- C DIS(I) = BOXDI DO 10,J=1,NNB CALL MINEL(DIS,NOV,JN,VN) IREGIO(J) = JN DIS(JN) = BOXDI 10 CONTINUE END C --------------------------------------------------------------------- SUBROUTINE ADDMIN ( LP, XDEPS, FDEPS, GDEPS, POINT, NOD, VAL, & XLL, XRL, DIS, GRAD, NOM, FMIN, XMINZ, & FMINZ, MXV, MXNM, MXNT, NFTOT, NGTOT ) C --------------------------------------------------------------------- IMPLICIT DOUBLE PRECISION (A-H,O-Z) DIMENSION POINT(NOD), XLL(NOD), XRL(NOD), GRAD(NOD) DIMENSION XMINZ(MXV,MXNM), FMINZ(MXNM), DIS(MXNT) LOGICAL LP CHARACTER*(*) FMIN c C Calculate the effective root mean square of the gradient. C Components "touching" the margins are not taken into account. CALL CHEMIN(GN,XDEPS,POINT,NOD,VAL,XLL,XRL,GRAD,NFC,NGC) NFTOT = NFTOT+NFC NGTOT = NGTOT+NGC IF (LP) WRITE (*,*) 'Effective RMS Gradient: ', GN IF (GN.GT.GDEPS) return NOM = NOM + 1 DO 13,J=1,NOD XMINZ(J,NOM) = POINT(J) 13 CONTINUE FMINZ(NOM) = VAL IF (NOM .EQ. 1) THEN OPEN (1,FILE=FMIN) CLOSE (1,STATUS='DELETE') OPEN (1,FILE=FMIN) WRITE (1,*) NOD C Always accept the first minimum. WRITE (*,*) '... The first minimum is found ...', & ' F: ', NFTOT, ' G: ', NGTOT ELSE C Check if this is a new minimum. Follow the steps: C I) At first accept it as a new minimum. C II) Then check if it is true and if not discard it. CALL DIST(XMINZ,NOD,NOM,NOM,DIS,MXV) CALL MINEL(DIS,NOM-1,JN,VMIN) IF ( ABS(FMINZ(JN)-VAL) .LE. & FDEPS*MAX(1.,ABS(VAL),ABS(FMINZ(JN))) .AND. & VMIN .LE. XDEPS ) THEN NOM = NOM - 1 IF (LP) WRITE (*,*) 'NOT a new MINIMUM' RETURN END IF WRITE (*,*) '... Number of minima found: ',NOM, & ' F: ', NFTOT, ' G: ', NGTOT END IF CALL WRITER(POINT,VAL,1,1,0) C END C --------------------------------------------------------------------- SUBROUTINE DUMPER ( NOD, NOF, ITTHR, HEAL, NEARN, NSAMPL, SIGMA, & COVTHR, XDEPS, FDEPS, GDEPS, NDUMP, FINP, & FOUT, FMIN, POIFIL, INDIR, ITC, NOM, NFTOT, & NGTOT, NHTOT, NJTOT, TNOS, XMINZ, FMINZ, XLL, & XRL, IXAT, MXV, MXNM, MXSP, NLOC ) C --------------------------------------------------------------------- IMPLICIT DOUBLE PRECISION (A-H,O-Z) CHARACTER*(*) FINP, FOUT, FMIN, POIFIL, INDIR DIMENSION XMINZ(MXV,MXNM), FMINZ(MXNM) DIMENSION XLL(NOD), XRL(NOD), IXAT(NOD) C DIMENSION IRLX(25) C OPEN (10,FILE='DUMP',FORM='UNFORMATTED') WRITE (10) NOD WRITE (10) NOF WRITE (10) ITTHR WRITE (10) HEAL WRITE (10) NEARN WRITE (10) NSAMPL WRITE (10) SIGMA WRITE (10) COVTHR WRITE (10) XDEPS WRITE (10) FDEPS WRITE (10) GDEPS WRITE (10) NDUMP WRITE (10) FINP WRITE (10) FOUT WRITE (10) FMIN WRITE (10) POIFIL WRITE (10) INDIR C WRITE (10) ITC WRITE (10) NOM WRITE (10) NFTOT WRITE (10) NGTOT WRITE (10) NHTOT WRITE (10) NJTOT WRITE (10) TNOS WRITE (10) MXSP WRITE (10) NLOC C DO 10,J=1,NOM DO 20,I=1,NOD WRITE (10) XMINZ(I,J) 20 CONTINUE 10 CONTINUE DO 30,I=1,NOM WRITE (10) FMINZ(I) 30 CONTINUE DO 40,I=1,NOD WRITE (10) XLL(I) 40 CONTINUE DO 50,I=1,NOD WRITE (10) XRL(I) 50 CONTINUE DO 60,I=1,NOD WRITE (10) IXAT(I) 60 CONTINUE C CALL RLUXUT(IRLX) DO 70,I=1,25 WRITE (10) IRLX(I) 70 CONTINUE CLOSE (10) WRITE (*,*) 'Current state has been saved.' END C --------------------------------------------------------------------- SUBROUTINE RESTOR ( NOD, NOF, ITTHR, HEAL, NEARN, NSAMPL, SIGMA, & COVTHR, XDEPS, FDEPS, GDEPS, NDUMP, FINP, & FOUT, FMIN, POIFIL, INDIR, ITC, NOM, NFTOT, & NGTOT, NHTOT, NJTOT, TNOS, XMINZ, FMINZ, XLL, & XRL, IXAT, MXV, MXNM, MXSP, NLOC ) IMPLICIT DOUBLE PRECISION (A-H,O-Z) CHARACTER*(*) FINP, FOUT, FMIN, POIFIL, INDIR DIMENSION XMINZ(MXV,MXNM), FMINZ(MXNM) DIMENSION XLL(MXV), XRL(MXV), IXAT(MXV) C DIMENSION IRLX(25) C OPEN (10,FILE='DUMP',FORM='UNFORMATTED') READ (10) NOD READ (10) NOF READ (10) ITTHR READ (10) HEAL READ (10) NEARN READ (10) NSAMPL READ (10) SIGMA READ (10) COVTHR READ (10) XDEPS READ (10) FDEPS READ (10) GDEPS READ (10) NDUMP READ (10) FINP READ (10) FOUT READ (10) FMIN READ (10) POIFIL READ (10) INDIR C READ (10) ITC READ (10) NOM READ (10) NFTOT READ (10) NGTOT READ (10) NHTOT READ (10) NJTOT READ (10) TNOS READ (10) MXSP READ (10) NLOC C DO 10,J=1,NOM DO 20,I=1,NOD READ (10) XMINZ(I,J) 20 CONTINUE 10 CONTINUE DO 30,I=1,NOM READ (10) FMINZ(I) 30 CONTINUE DO 40,I=1,NOD READ (10) XLL(I) 40 CONTINUE DO 50,I=1,NOD READ (10) XRL(I) 50 CONTINUE DO 60,I=1,NOD READ (10) IXAT(I) 60 CONTINUE C DO 70,I=1,25 READ (10) IRLX(I) 70 CONTINUE CALL RLUXIN(IRLX) CLOSE (10) WRITE (*,*) 'Previous state has been restored.' END C --------------------------------------------------------------------- SUBROUTINE FROMFL ( INFIL, NOD, NOF, ITTHR, HEAL, NEARN, NSAMPL, & SIGMA, COVTHR, XDEPS, FDEPS, GDEPS, IREST, & NDUMP, FINP, FOUT, FMIN, POIFIL, INDIR, POINT, & XLL, XRL, IXAT, MXV, MXNNN, MXNS ) C --------------------------------------------------------------------- IMPLICIT DOUBLE PRECISION (A-H,O-Z) CHARACTER*(*) INFIL, FINP, FOUT, FMIN, POIFIL, INDIR DIMENSION POINT(MXV), XLL(MXV), XRL(MXV), IXAT(MXV) C LOGICAL EXI1, EXI2, MEXI, ERR C C Check if the input file exists. C If it does, read the parameters from there. C Otherwise use the defaults & create a sample input file, C containing the default values, as an aid to the user. C INQUIRE (FILE=INFIL,EXIST=EXI1) IF (EXI1) THEN C Open and read the data file. OPEN (1,FILE=INFIL) READ (1,*) NOD READ (1,*) NOF READ (1,*) ITTHR READ (1,*) HEAL READ (1,*) NEARN READ (1,*) NSAMPL READ (1,*) SIGMA READ (1,*) COVTHR READ (1,*) XDEPS READ (1,*) FDEPS READ (1,*) GDEPS READ (1,*) IREST READ (1,*) NDUMP READ (1,'(A)') FINP READ (1,'(A)') FOUT READ (1,'(A)') FMIN READ (1,'(A)') POIFIL READ (1,'(A)') INDIR CLOSE (1) ELSE C Create a sample input file with default values. OPEN (1,FILE=INFIL,STATUS='NEW') WRITE (1,'(1X,I8,T25,A)') NOD, 'NOD' WRITE (1,'(1X,I8,T25,A)') NOF, 'NOF' WRITE (1,'(1X,I8,T25,A)') ITTHR, 'ITTHR' WRITE (1,'(1X,G21.14,T25,A)') HEAL, 'HEAL' WRITE (1,'(1X,I8,T25,A)') NEARN, 'NEARN' WRITE (1,'(1X,I8,T25,A)') NSAMPL, 'NSAMPL' WRITE (1,'(1X,G21.14,T25,A)') SIGMA, 'SIGMA' WRITE (1,'(1X,G21.14,T25,A)') COVTHR, 'COVTHR' WRITE (1,'(1X,G21.14,T25,A)') XDEPS, 'XDEPS' WRITE (1,'(1X,G21.14,T25,A)') FDEPS, 'FDEPS' WRITE (1,'(1X,G21.14,T25,A)') GDEPS, 'GDEPS' WRITE (1,'(1X,I8,T25,A)') IREST, 'IREST' WRITE (1,'(1X,I8,T25,A)') NDUMP, 'NDUMP' WRITE (1,'(A,T25,A)') FINP , 'FINP' WRITE (1,'(A,T25,A)') FOUT , 'FOUT' WRITE (1,'(A,T25,A)') FMIN , 'FMIN' WRITE (1,'(A,T25,A)') POIFIL, 'POIFIL' WRITE (1,'(A,T25,A)') INDIR, 'INDIR' CLOSE (1) WRITE (*,*) 'TML: File ',INFIL,' is created with ', & 'default values.' END IF C C Make sure NOD is correct. IF (NOD.LT.1) THEN WRITE (*,*) 'TML: Incorrect number of parameters ', & '(must be positive).' STOP END IF IF (NOD.GT.MXV) THEN WRITE (*,*) 'TML: The number of parameters is set to: ', & NOD,' > ', MXV,' which is the maximum allowed.' STOP END IF C C Perform initialization of required arrays. (POINT,XLL,XRL & IXAT) C INQUIRE (FILE=POIFIL,EXIST=EXI2) IF (EXI2) THEN C Read parameter bounds. OPEN (1,FILE=POIFIL,STATUS='OLD') DO 20,I=1,NOD READ (1,*) POINT(I), XLL(I), XRL(I), IXAT(I) 20 CONTINUE CLOSE (1) ELSE C Create a sample parameter file with default values. OPEN (1,FILE=POIFIL,STATUS='NEW') DO 9,I=1,NOD POINT(I) = 2*RANM()-1. WRITE (1,120) POINT(I), XLL(I), XRL(I), IXAT(I) 9 CONTINUE CLOSE (1) WRITE (*,*) 'TML: File ',POIFIL,' is created with ', & 'default values.' END IF C C Check if the input is meaningful. ERR = .FALSE. IF (NOF.LT.0) THEN WRITE (*,*) 'TML: Incorrect number of function terms ', & '(must be zero or positive).' STOP END IF C IF (ITTHR.LE.0) THEN WRITE (*,*) 'TML: Incorrect iteration threshold ',ITTHR, & ' (must be positive).' ERR = .TRUE. END IF C IF (NEARN.LE.0) THEN WRITE (*,*) 'TML: Incorrect NEARN ', NEARN, & ' (must be positive).' ERR = .TRUE. END IF IF (NEARN.GT.MXNNN) THEN WRITE (*,*) 'TML: NEARN is set to: ',NEARN,' > ', & MXNNN,' which is the maximum allowed.' ERR = .TRUE. END IF C IF (NSAMPL.LE.0) THEN WRITE (*,*) 'TML: Incorrect NSAMPL ', NSAMPL, & ' (must be positive).' ERR = .TRUE. END IF IF (NSAMPL.GT.MXNS) THEN WRITE (*,*) 'TML: NSAMPL is set to: ', NSAMPL,' > ', & MXNS,' which is the maximum allowed.' ERR = .TRUE. END IF C IF (SIGMA.LE.0.) THEN WRITE (*,*) 'TML: Incorrect SIGMA ', SIGMA ERR = .TRUE. END IF C IF (COVTHR.LE.0. .OR. COVTHR.GT.1.) THEN WRITE (*,*) 'TML: Incorrect COVTHR ', COVTHR ERR = .TRUE. END IF C IF (XDEPS.LE.0. .OR. XDEPS.GE.1.) THEN WRITE (*,*) 'TML: Incorrect XDEPS ', XDEPS ERR = .TRUE. END IF C IF (FDEPS.LE.0. .OR. FDEPS.GE.1.) THEN WRITE (*,*) 'TML: Incorrect FDEPS ', FDEPS ERR = .TRUE. END IF C IF (GDEPS.LE.0. .OR. GDEPS.GE.1.) THEN WRITE (*,*) 'TML: Incorrect GDEPS ', GDEPS ERR = .TRUE. END IF C IF (IREST.NE.0 .AND. IREST.NE.1) THEN WRITE (*,*) 'TML: Incorrect IREST ', IREST ERR = .TRUE. END IF C IF (NDUMP.LT.0) THEN WRITE (*,*) 'TML: Incorrect NDUMP ', NDUMP ERR = .TRUE. END IF C IF (FINP.EQ.' ') THEN WRITE (*,*) 'TML: Merlin input file not specified.' ERR = .TRUE. END IF C IF (FOUT.EQ.' ') THEN WRITE (*,*) 'TML: Merlin output file not specified.' ERR = .TRUE. END IF C INQUIRE (FILE=FINP,EXIST=MEXI) IF (.NOT.MEXI) THEN LEFI = LENGTH(FINP) WRITE (*,*) 'TML: The Merlin input file ',FINP(1:LEFI), & ' does not exist.' END IF C C Missing input files ? IF (.NOT. (EXI1.AND.EXI2)) THEN WRITE (*,*) 'TML: Please edit the iput files ', & 'to provide the desired input.' STOP END IF C C Check if parameters, bounds and fix statuses are meaningful. DO 100,I=1,NOD IF (XLL(I).GE.XRL(I)) THEN WRITE (*,*) 'TML: Incorrect bounds for parameter ', I ERR = .TRUE. ELSE IF (POINT(I).LT.XLL(I) .OR. POINT(I).GT.XRL(I)) THEN WRITE (*,*) 'TML: Out of bounds parameter ', I ERR = .TRUE. END IF IF (IXAT(I).NE.0 .AND. IXAT(I).NE.1) THEN WRITE (*,*) 'TML: Incorrect fix status for variable ', I ERR = .TRUE. END IF 100 CONTINUE C C Any error ? IF (ERR) STOP C 120 FORMAT (F8.5,3X,2(F4.1,3X),I1) C END