hdiff output

r33230/Dijinit.f90 2017-08-27 10:30:16.494218340 +0100 r33229/Dijinit.f90 2017-08-27 10:30:17.686234232 +0100
 22: !  Dijkstra connection algorithm for pathsample. 22: !  Dijkstra connection algorithm for pathsample.
 23: ! 23: !
 24: SUBROUTINE DIJINIT(NWORST) 24: SUBROUTINE DIJINIT(NWORST)
 25: USE PORFUNCS 25: USE PORFUNCS
 26: USE COMMONS 26: USE COMMONS
 27: USE UTILS,ONLY : GETUNIT 27: USE UTILS,ONLY : GETUNIT
 28: IMPLICIT NONE 28: IMPLICIT NONE
 29:  29: 
 30: INTEGER J1, J2, J4, PARENT(NMIN), JMINW, NPERM, J5, LJ1, LJ2, NWORST, NSTEPS, NMINSTART, NMINEND, J6, MUNIT, J7 30: INTEGER J1, J2, J4, PARENT(NMIN), JMINW, NPERM, J5, LJ1, LJ2, NWORST, NSTEPS, NMINSTART, NMINEND, J6, MUNIT, J7
 31: INTEGER JN, JM, NPOSITION, LUNIT 31: INTEGER JN, JM, NPOSITION, LUNIT
 32: INTEGER NMINGAP, NPRUNEDONE, NPRUNEPAIRS, NPRUNEMIN,NPRUNEPAIRSOLD, RELAXED(NMIN), NNEIGH 32: INTEGER NMINGAP, NPRUNEDONE, NPRUNEPAIRS, NPRUNEMIN,NPRUNEPAIRSOLD
 33: INTEGER, ALLOCATABLE :: LOCATIONSTART(:), LOCATIONEND(:), PRUNEPAIRS(:,:) 33: INTEGER, ALLOCATABLE :: LOCATIONSTART(:), LOCATIONEND(:), PRUNEPAIRS(:,:)
 34: LOGICAL PERMANENT(NMIN), ISA(NMIN), ISB(NMIN), ISSTART(NMIN), NOTDONE, PRUNEMIN(NMIN), REDODIJKSTRA 34: LOGICAL PERMANENT(NMIN), ISA(NMIN), ISB(NMIN), ISSTART(NMIN), NOTDONE, PRUNEMIN(NMIN), REDODIJKSTRA
 35: DOUBLE PRECISION MINWEIGHT, DUMMY, TNEW, ELAPSED, PFTOTALSTART, HUGESAVE, THRESH, LPOINTS1(NOPT), LPOINTS2(NOPT) 35: DOUBLE PRECISION MINWEIGHT, DUMMY, TNEW, ELAPSED, PFTOTALSTART, HUGESAVE, THRESH, LPOINTS1(NOPT), LPOINTS2(NOPT)
 36: DOUBLE PRECISION MAXWEIGHT, SCALEFAC, PDMAX, PD, MINGAPTHRESH, DIST2, DISTANCE, RMAT(3,3), MINNWEIGHT, MINNDIST 36: DOUBLE PRECISION MAXWEIGHT, SCALEFAC, PDMAX, PD, MINGAPTHRESH, DIST2, DISTANCE, RMAT(3,3)
 37: ! 37: !
 38: ! KIND=16 is not supported by Portland. If you want extra precision, uncomment the following line 38: ! KIND=16 is not supported by Portland. If you want extra precision, uncomment the following line
 39: ! and use NAG. 39: ! and use NAG.
 40: ! 40: !
 41: ! REAL(KIND=16) :: TMPWEIGHT, WEIGHT(NMIN) 41: ! REAL(KIND=16) :: TMPWEIGHT, WEIGHT(NMIN)
 42: REAL(KIND=8) :: TMPWEIGHT, WEIGHT(NMIN) 42: REAL(KIND=8) :: TMPWEIGHT, WEIGHT(NMIN)
 43:  43: 
 44: IF (DIJPRUNET) THEN 44: IF (DIJPRUNET) THEN
 45:    MUNIT=GETUNIT() 45:    MUNIT=GETUNIT()
 46:    IF (.NOT.(PRUNECYCLET)) NPRUNE=1 46:    IF (.NOT.(PRUNECYCLET)) NPRUNE=1
 95:          PDMAX=ABS(ALLPAIRS(J2)) 95:          PDMAX=ABS(ALLPAIRS(J2))
 96: !        IF (DEBUG) PRINT '(A,G20.10)','Dijinit> maximum neighbour metric value increased to',PDMAX  96: !        IF (DEBUG) PRINT '(A,G20.10)','Dijinit> maximum neighbour metric value increased to',PDMAX 
 97: !        IF (DEBUG) PRINT '(A,I8,G20.10)','Dijinit> J2,ALLPAIRS=',J2,ALLPAIRS(J2) 97: !        IF (DEBUG) PRINT '(A,I8,G20.10)','Dijinit> J2,ALLPAIRS=',J2,ALLPAIRS(J2)
 98:       ENDIF 98:       ENDIF
 99:    ENDDO 99:    ENDDO
100: ELSE100: ELSE
101:    DO J2=1,NMIN101:    DO J2=1,NMIN
102:       DO J5=1,PAIRDISTMAX102:       DO J5=1,PAIRDISTMAX
103:          IF (PAIRDIST(J2,J5).GT.PDMAX) THEN103:          IF (PAIRDIST(J2,J5).GT.PDMAX) THEN
104:             PDMAX=PAIRDIST(J2,J5)104:             PDMAX=PAIRDIST(J2,J5)
105:             IF (DEBUG) PRINT '(A,G20.10)','Dijinit> maximum neighbour metric value increased to',PDMAX 105: !           IF (DEBUG) PRINT '(A,G20.10)','Dijinit> maximum neighbour metric value increased to',PDMAX 
106:             IF (DEBUG) PRINT '(A,2I8,G20.10,I8)','Dijinit> J2,J5,PAIRDIST,PAIRLIST=',J2,J5,PAIRDIST(J2,J5),PAIRLIST(J2,J5)106: !           IF (DEBUG) PRINT '(A,2I8,G20.10)','Dijinit> J2,J5,PAIRDIST=',J2,J5,PAIRDIST(J2,J5)
107:          ENDIF107:          ENDIF
108:       ENDDO108:       ENDDO
109:    ENDDO109:    ENDDO
110:    PRINT '(A,G20.10)','Dijinit> maximum neighbour metric value=',PDMAX110:    PRINT '(A,G20.10)','Dijinit> maximum neighbour metric value=',PDMAX
111: ENDIF111: ENDIF
112: 112: 
113: IF (MINGAPT) THEN113: IF (MINGAPT) THEN
114:    IF (MINGAPRATIOT) THEN114:    IF (MINGAPRATIOT) THEN
115:       MINGAPTHRESH=PDMAX*MINGAPINP115:       MINGAPTHRESH=PDMAX*MINGAPINP
116:    ELSE116:    ELSE
128: MAXWEIGHT=HUGE(1.0D0)/1.0D1128: MAXWEIGHT=HUGE(1.0D0)/1.0D1
129: ! MAXWEIGHT=1.0D6129: ! MAXWEIGHT=1.0D6
130: loopstart: DO J1=1,NMINSTART ! cycle over all minima in the starting state130: loopstart: DO J1=1,NMINSTART ! cycle over all minima in the starting state
131:    SCALEFAC=1.0D0131:    SCALEFAC=1.0D0
132: 222   LJ1=LOCATIONSTART(J1)132: 222   LJ1=LOCATIONSTART(J1)
133:    WEIGHT(1:NMIN)=HUGE(1.0D0)133:    WEIGHT(1:NMIN)=HUGE(1.0D0)
134:    HUGESAVE=WEIGHT(1)134:    HUGESAVE=WEIGHT(1)
135:    WEIGHT(LJ1)=0.0D0135:    WEIGHT(LJ1)=0.0D0
136:    PERMANENT(1:NMIN)=.FALSE.136:    PERMANENT(1:NMIN)=.FALSE.
137:    PERMANENT(LJ1)=.TRUE.137:    PERMANENT(LJ1)=.TRUE.
138:    RELAXED(1:NMIN)=0 
139:    NPERM=1138:    NPERM=1
140:    PARENT(1:NMIN)=0 ! parent is initially undefined139:    PARENT(1:NMIN)=0 ! parent is initially undefined
141:    J4=LJ1140:    J4=LJ1
142:    dijkstraloop: DO141:    dijkstraloop: DO
143:       NNEIGH=0 
144:       MINNWEIGHT=HUGE(1.0D0) 
145:       MINNDIST=HUGE(1.0D0) 
146:       DO J2=1,NMIN142:       DO J2=1,NMIN
147:          IF (J2.EQ.J4) CYCLE143:          IF (J2.EQ.J4) CYCLE
148:          IF (PERMANENT(J2)) CYCLE144:          IF (PERMANENT(J2)) CYCLE
149:          PD=1.0D4*PDMAX145:          PD=1.0D4*PDMAX
150:          JM=MIN(J4,J2)146:          JM=MIN(J4,J2)
151:          JN=MAX(J4,J2)147:          JN=MAX(J4,J2)
152:          NPOSITION=((JN-2)*(JN-1))/2+JM148:          NPOSITION=((JN-2)*(JN-1))/2+JM
153: !        IF (INITIALDIST) PRINT '(A,I8,A,I8,A,I10,A,G20.10)','Dijinit> minima ',J2,' and ',J4,' position ',NPOSITION,' distance ',ALLPAIRS(NPOSITION)149: !        IF (INITIALDIST) PRINT '(A,I8,A,I8,A,I10,A,G20.10)','Dijinit> minima ',J2,' and ',J4,' position ',NPOSITION,' distance ',ALLPAIRS(NPOSITION)
154:          IF (.NOT.PAIRSIGNORET) THEN !for pruning the database all minima count not just the ones not searched yet150:          IF (.NOT.PAIRSIGNORET) THEN !for pruning the database all minima count not just the ones not searched yet
155:             DO J5=1,NPAIRDONE ! skip151:             DO J5=1,NPAIRDONE ! skip
221:                   GOTO 973217:                   GOTO 973
222:                ENDIF218:                ENDIF
223:             ENDDO219:             ENDDO
224:          ENDIF 220:          ENDIF 
225:          IF (INITIALDIST) THEN221:          IF (INITIALDIST) THEN
226:             PD=ABS(ALLPAIRS(NPOSITION))222:             PD=ABS(ALLPAIRS(NPOSITION))
227:          ELSE223:          ELSE
228:             DO J5=1,PAIRDISTMAX224:             DO J5=1,PAIRDISTMAX
229:                IF (PAIRLIST(J4,J5).EQ.J2) THEN225:                IF (PAIRLIST(J4,J5).EQ.J2) THEN
230:                   PD=PAIRDIST(J4,J5)226:                   PD=PAIRDIST(J4,J5)
231:                   NNEIGH=NNEIGH+1 
232:                   IF (WEIGHT(J2).LT.MINNWEIGHT) MINNWEIGHT=WEIGHT(J2) 
233:                   IF (PD.LT.MINNDIST) MINNDIST=PD 
234:                   GOTO 973227:                   GOTO 973
235:                ENDIF228:                ENDIF
236:             ENDDO229:             ENDDO
237:             DO J5=1,PAIRDISTMAX230:             DO J5=1,PAIRDISTMAX
238:                IF (PAIRLIST(J2,J5).EQ.J4) THEN231:                IF (PAIRLIST(J2,J5).EQ.J4) THEN
239:                   PD=PAIRDIST(J2,J5)232:                   PD=PAIRDIST(J2,J5)
240:                   NNEIGH=NNEIGH+1 
241:                   IF (WEIGHT(J2).LT.MINNWEIGHT) MINNWEIGHT=WEIGHT(J2) 
242:                   IF (PD.LT.MINNDIST) MINNDIST=PD 
243:                   GOTO 973233:                   GOTO 973
244:                ENDIF234:                ENDIF
245:             ENDDO235:             ENDDO
246:          ENDIF236:          ENDIF
247: 973      CONTINUE237: 973      CONTINUE
248:          TMPWEIGHT=PD*SCALEFAC 238:          TMPWEIGHT=PD*SCALEFAC 
249:          IF (TMPWEIGHT.LT.HUGE(1.0D0)/10.0D0) THEN ! don;t raise a huge number to any power!239:          IF (TMPWEIGHT.LT.HUGE(1.0D0)/10.0D0) THEN ! don;t raise a huge number to any power!
250:             IF (INDEXCOSTFUNCTION) THEN 240:             IF (INDEXCOSTFUNCTION) THEN 
251:                IF (TMPWEIGHT.EQ.0.0D0) THEN ! minima are connected!241:                IF (TMPWEIGHT.EQ.0.0D0) THEN ! minima are connected!
252:                ELSE242:                ELSE
274:                ELSEIF (COSTFUNCTIONPOWER.EQ.-1) THEN264:                ELSEIF (COSTFUNCTIONPOWER.EQ.-1) THEN
275:                   TMPWEIGHT=1.0D0/TMPWEIGHT265:                   TMPWEIGHT=1.0D0/TMPWEIGHT
276:                ELSE266:                ELSE
277:                   TMPWEIGHT=TMPWEIGHT**COSTFUNCTIONPOWER 267:                   TMPWEIGHT=TMPWEIGHT**COSTFUNCTIONPOWER 
278:                ENDIF268:                ENDIF
279:             ENDIF269:             ENDIF
280:          ENDIF270:          ENDIF
281:          271:          
282: !        PRINT '(A,2I10,3G20.10)','J2,J4,TMPWEIGHT,WEIGHT(J4),WEIGHT(J2)=',J2,J4,TMPWEIGHT,WEIGHT(J4),WEIGHT(J2)272: !        PRINT '(A,2I10,3G20.10)','J2,J4,TMPWEIGHT,WEIGHT(J4),WEIGHT(J2)=',J2,J4,TMPWEIGHT,WEIGHT(J4),WEIGHT(J2)
283:          IF (TMPWEIGHT+WEIGHT(J4).LT.WEIGHT(J2)) THEN ! relax J2273:          IF (TMPWEIGHT+WEIGHT(J4).LT.WEIGHT(J2)) THEN ! relax J2
284:             RELAXED(J2)=RELAXED(J2)+1 
285:             WEIGHT(J2)=WEIGHT(J4)+TMPWEIGHT274:             WEIGHT(J2)=WEIGHT(J4)+TMPWEIGHT
286:             PARENT(J2)=J4275:             PARENT(J2)=J4
287: !           PRINT '(A,2I10)','J2,PARENT=',J2,PARENT(J2)276: !           PRINT '(A,2I10)','J2,PARENT=',J2,PARENT(J2)
288:          ENDIF277:          ENDIF
289:       ENDDO278:       ENDDO
290: 279: 
291:       MINWEIGHT=HUGE(1.0D0)280:       MINWEIGHT=HUGE(1.0D0)
292:       NOTDONE=.TRUE.281:       NOTDONE=.TRUE.
293:       DO J2=1,NMIN282:       DO J2=1,NMIN
294:          IF (.NOT.PERMANENT(J2)) THEN283:          IF (.NOT.PERMANENT(J2)) THEN
295:             IF (WEIGHT(J2).LT.MINWEIGHT) THEN284:             IF (WEIGHT(J2).LT.MINWEIGHT) THEN
296:                MINWEIGHT=WEIGHT(J2)285:                MINWEIGHT=WEIGHT(J2)
297:                JMINW=J2286:                JMINW=J2
298:                NOTDONE=.FALSE.287:                NOTDONE=.FALSE.
299:             ENDIF288:             ENDIF
300:          ENDIF289:          ENDIF
301: !        PRINT '(A,I10,L5,2G20.10,I10)','J2,PERMANENT,WEIGHT,MINWEIGHT,JMINW=',J2,PERMANENT(J2),WEIGHT(J2),MINWEIGHT,JMINW290: !        PRINT '(A,I10,L5,2G20.10,I10)','J2,PERMANENT,WEIGHT,MINWEIGHT,JMINW=',J2,PERMANENT(J2),WEIGHT(J2),MINWEIGHT,JMINW
302:       ENDDO291:       ENDDO
303:       IF (NOTDONE) THEN292:       IF (NOTDONE) THEN
304:          PRINT '(A,I8,A,I8)','dijinit> WARNING - JMINW not set - value=',JMINW,' J4=',J4293:          PRINT '(A,I8,A,I8)','dijinit> WARNING - JMINW not set - value=',JMINW,' J4=',J4
305:          PRINT '(A,I8)','dijinit> NPERM=',NPERM294: !        PRINT '(A,I8)','dijinit> NPERM=',NPERM
306:          DO J2=1,NMIN295:          DO J2=1,NMIN
307:             PRINT '(A,I8,L5,2G20.10)','J2,PERMANENT,WEIGHT,MINWEIGHT=',J2,PERMANENT(J2),WEIGHT(J2),MINWEIGHT296: !           PRINT '(A,I8,L5,2G20.10)','J2,PERMANENT,WEIGHT,MINWEIGHT=',J2,PERMANENT(J2),WEIGHT(J2),MINWEIGHT
308:             IF (.NOT.PERMANENT(J2)) THEN297:             IF (.NOT.PERMANENT(J2)) THEN
309:                IF (WEIGHT(J2).LT.MINWEIGHT) THEN298:                IF (WEIGHT(J2).LT.MINWEIGHT) THEN
310:                   MINWEIGHT=WEIGHT(J2)299:                   MINWEIGHT=WEIGHT(J2)
311:                   JMINW=J2300:                   JMINW=J2
312:                   NOTDONE=.FALSE.301:                   NOTDONE=.FALSE.
313:                ENDIF302:                ENDIF
314:             ENDIF303:             ENDIF
315:          ENDDO304:          ENDDO
316:          STOP !!! DJW305:          STOP !!! DJW
317:       ENDIF306:       ENDIF
318: 307: 
319:       IF (.NOT.INITIALDIST) PRINT '(A,I10,A,I10,A,2G20.10)','Dijinit> Number of non-permanent nodes in neighbour list for ', & 
320:   &                       J4,' is ',NNEIGH,' min weight and dist=',MINNWEIGHT,MINNDIST 
321:       J4=JMINW308:       J4=JMINW
322:       PERMANENT(J4)=.TRUE.309:       PERMANENT(J4)=.TRUE.
323:       NPERM=NPERM+1310:       NPERM=NPERM+1
324: !     PRINT '(A,2I8,G20.10,I8)','permanent minimum J4,NPERM,WEIGHT,relaxations=',J4,NPERM,WEIGHT(J4),RELAXED(J4)311: !     PRINT '(A,2I8,G20.10)','permanent minimum J4,NPERM,WEIGHT=',J4,NPERM,WEIGHT(J4) 
325:       IF (WEIGHT(J4).GT.MAXWEIGHT) THEN312:       IF (WEIGHT(J4).GT.MAXWEIGHT) THEN
326:          SCALEFAC=SCALEFAC/10.0D0313:          SCALEFAC=SCALEFAC/10.0D0
327:          PRINT '(A,G20.10)','dijinit> Maximum weight is too large - scaling by ',SCALEFAC314:          PRINT '(A,G20.10)','dijinit> Maximum weight is too large - scaling by ',SCALEFAC
328:          GOTO 222315:          GOTO 222
329:       ENDIF316:       ENDIF
330: 317: 
331:       IF (NPERM.EQ.NMIN) EXIT dijkstraloop318:       IF (NPERM.EQ.NMIN) EXIT dijkstraloop
332: 319: 
333:    ENDDO dijkstraloop320:    ENDDO dijkstraloop
334: 321: 
345: IF (MINGAPT) NMINGAP=0332: IF (MINGAPT) NMINGAP=0
346: PRINT '(A)','Dijinit> Summary of best path based on missing connection metric - note distance scaling is removed'333: PRINT '(A)','Dijinit> Summary of best path based on missing connection metric - note distance scaling is removed'
347: PRINT '(A)','    min1          energy        min2          energy             metric          edge weight            weight'334: PRINT '(A)','    min1          energy        min2          energy             metric          edge weight            weight'
348: DO 335: DO 
349:    IF (PARENT(J5).EQ.0) THEN336:    IF (PARENT(J5).EQ.0) THEN
350:       PRINT '(A,I6,A)','Dijinit> ERROR - parent for J5=',J5,' is zero'337:       PRINT '(A,I6,A)','Dijinit> ERROR - parent for J5=',J5,' is zero'
351:       PRINT '(A)',     'Dijinit> Suggests all possible pairs have been tried!'338:       PRINT '(A)',     'Dijinit> Suggests all possible pairs have been tried!'
352:       STOP339:       STOP
353:    ENDIF340:    ENDIF
354:    DUMMY=1.0D4*PDMAX*SCALEFAC341:    DUMMY=1.0D4*PDMAX*SCALEFAC
355: !  PRINT '(A,4G20.10)','PDMAX,SCALEFAC,DUMMY,DUMMY/SCALEFAC=',PDMAX,SCALEFAC,DUMMY,DUMMY/SCALEFAC 
356:    IF (.NOT.PAIRSIGNORET) THEN342:    IF (.NOT.PAIRSIGNORET) THEN
357:      DO J2=1,NPAIRDONE ! skip343:      DO J2=1,NPAIRDONE ! skip
358:        IF ((PAIR1(J2).EQ.J5).AND.(PAIR2(J2).EQ.PARENT(J5))) THEN344:        IF ((PAIR1(J2).EQ.J5).AND.(PAIR2(J2).EQ.PARENT(J5))) THEN
359:           IF (INITIALDIST) THEN345:           IF (INITIALDIST) THEN
360:              JM=MIN(J5,PARENT(J5))346:              JM=MIN(J5,PARENT(J5))
361:              JN=MAX(J5,PARENT(J5))347:              JN=MAX(J5,PARENT(J5))
362:              NPOSITION=((JN-2)*(JN-1))/2+JM348:              NPOSITION=((JN-2)*(JN-1))/2+JM
363:              DUMMY=ABS(ALLPAIRS(NPOSITION))*SCALEFAC349:              DUMMY=ABS(ALLPAIRS(NPOSITION))*SCALEFAC
364:              IF (DUMMY.EQ.0.0D0) THEN350:              IF (DUMMY.EQ.0.0D0) THEN
365:                 GOTO 864351:                 GOTO 864


r33230/mindouble.f90 2017-08-27 10:30:16.738221591 +0100 r33229/mindouble.f90 2017-08-27 10:30:17.966237965 +0100
 94:       DEALLOCATE(MINFRQ2) 94:       DEALLOCATE(MINFRQ2)
 95:       ALLOCATE(MINFRQ2(2*MAXMIN)) 95:       ALLOCATE(MINFRQ2(2*MAXMIN))
 96:       MINFRQ2(1:MAXMIN)=VDP(1:MAXMIN) 96:       MINFRQ2(1:MAXMIN)=VDP(1:MAXMIN)
 97:    ENDIF 97:    ENDIF
 98:    98:   
 99: ! 99: !
100: ! If the PAIRDIST vector has not been zeroed in a MERGEDBT run executing the100: ! If the PAIRDIST vector has not been zeroed in a MERGEDBT run executing the
101: ! next block can give a SIGFPE101: ! next block can give a SIGFPE
102: !102: !
103:    IF (DIJINITT.AND.INITIALDIST) THEN103:    IF (DIJINITT.AND.INITIALDIST) THEN
104:       PRINT '(A,I20,A)','mindouble> Increasing pair distance array dimension to ',2*MAXMIN,' minima '104:       PRINT '(A,I20)','mindouble> Increasing pair distance array dimension to ',2*MAXMIN,' minima '
105:       IF (ALLOCATED(VDPVEC)) DEALLOCATE(VDPVEC)105:       IF (ALLOCATED(VDPVEC)) DEALLOCATE(VDPVEC)
106:       OLDSIZE=(MAXMIN*(MAXMIN-1))/2106:       OLDSIZE=(MAXMIN*(MAXMIN-1))/2
107:       NEWSIZE=(2*MAXMIN*(2*MAXMIN-1))/2107:       NEWSIZE=(2*MAXMIN*(2*MAXMIN-1))/2
108:       ALLOCATE(VDPVEC(OLDSIZE))108:       ALLOCATE(VDPVEC(OLDSIZE))
109:       VDPVEC(1:OLDSIZE)=ALLPAIRS(1:OLDSIZE)109:       VDPVEC(1:OLDSIZE)=ALLPAIRS(1:OLDSIZE)
110:       DEALLOCATE(ALLPAIRS)110:       DEALLOCATE(ALLPAIRS)
111:       ALLOCATE(ALLPAIRS(NEWSIZE))111:       ALLOCATE(ALLPAIRS(NEWSIZE))
112:       ALLPAIRS(1:NEWSIZE)=1.0D100112:       ALLPAIRS(1:NEWSIZE)=1.0D100
113:       ALLPAIRS(1:OLDSIZE)=VDPVEC(1:OLDSIZE)113:       ALLPAIRS(1:OLDSIZE)=VDPVEC(1:OLDSIZE)
114:       DEALLOCATE(VDPVEC)114:       DEALLOCATE(VDPVEC)


r33230/minpermdist.f90 2017-08-27 10:30:17.030225484 +0100 r33229/minpermdist.f90 2017-08-27 10:30:18.250241752 +0100
206: ENDIF206: ENDIF
207: !207: !
208: ! It is possible for the standard orientation to result in a distance that is worse than208: ! It is possible for the standard orientation to result in a distance that is worse than
209: ! the starting distance. Hence we need to set XBEST here.209: ! the starting distance. Hence we need to set XBEST here.
210: !210: !
211: DUMMYA(1:3*NATOMS)=COORDSA(1:3*NATOMS)211: DUMMYA(1:3*NATOMS)=COORDSA(1:3*NATOMS)
212: DUMMYB(1:3*NATOMS)=COORDSB(1:3*NATOMS)212: DUMMYB(1:3*NATOMS)=COORDSB(1:3*NATOMS)
213: DBEST=1.0D100213: DBEST=1.0D100
214: CALL NEWMINDIST(DUMMYB,DUMMYA,NATOMS,DISTANCE,BULKT,TWOD,'AX    ',.FALSE.,RIGID,DEBUG,RMAT)214: CALL NEWMINDIST(DUMMYB,DUMMYA,NATOMS,DISTANCE,BULKT,TWOD,'AX    ',.FALSE.,RIGID,DEBUG,RMAT)
215: DBEST=DISTANCE**2215: DBEST=DISTANCE**2
216: ! IF (DEBUG) PRINT '(A,G20.10)',' minpermdist> Initial distance before standard orientation=',DISTANCE216: IF (DEBUG) PRINT '(A,G20.10)',' minpermdist> Initial distance before standard orientation=',DISTANCE
217: 217: 
218: ! If global translation is possible (i.e. no frozen atoms and not periodic boundary conditions) then218: ! If global translation is possible (i.e. no frozen atoms and not periodic boundary conditions) then
219: ! translate the origin of structure A to the CoM of structure B.219: ! translate the origin of structure A to the CoM of structure B.
220: IF ((NFREEZE.GT.0).OR.BULKT.OR.MKTRAPT.OR.NOTRANSROTT) THEN220: IF ((NFREEZE.GT.0).OR.BULKT.OR.MKTRAPT.OR.NOTRANSROTT) THEN
221:    XBEST(1:3*NATOMS)=DUMMYA(1:3*NATOMS)221:    XBEST(1:3*NATOMS)=DUMMYA(1:3*NATOMS)
222: ELSE222: ELSE
223:    DO J1=1,NATOMS223:    DO J1=1,NATOMS
224:       XBEST(3*(J1-1)+1)=DUMMYA(3*(J1-1)+1)-CMBX224:       XBEST(3*(J1-1)+1)=DUMMYA(3*(J1-1)+1)-CMBX
225:       XBEST(3*(J1-1)+2)=DUMMYA(3*(J1-1)+2)-CMBY225:       XBEST(3*(J1-1)+2)=DUMMYA(3*(J1-1)+2)-CMBY
226:       XBEST(3*(J1-1)+3)=DUMMYA(3*(J1-1)+3)-CMBZ226:       XBEST(3*(J1-1)+3)=DUMMYA(3*(J1-1)+3)-CMBZ


r33230/setup.f 2017-08-27 10:30:17.426230757 +0100 r33229/setup.f 2017-08-27 10:30:18.518245324 +0100
1503:       IF (DIJINITT) THEN1503:       IF (DIJINITT) THEN
1504:          DO J1=1,NTS1504:          DO J1=1,NTS
1505: !1505: !
1506: ! JMC n.b. don't apply the nconnmin criteria at this point, hence the huge(1) 's 1506: ! JMC n.b. don't apply the nconnmin criteria at this point, hence the huge(1) 's 
1507: ! in place of NCONN() for the plus and minus minima.1507: ! in place of NCONN() for the plus and minus minima.
1508: !1508: !
1509:             CALL CHECKTS(ETS(J1),EMIN(PLUS(J1)),EMIN(MINUS(J1)),KPLUS(J1),KMINUS(J1),HUGE(1),HUGE(1), 1509:             CALL CHECKTS(ETS(J1),EMIN(PLUS(J1)),EMIN(MINUS(J1)),KPLUS(J1),KMINUS(J1),HUGE(1),HUGE(1), 
1510:      &                   PLUS(J1),MINUS(J1),.TRUE.,CUT_UNDERFLOW,DEADTS)1510:      &                   PLUS(J1),MINUS(J1),.TRUE.,CUT_UNDERFLOW,DEADTS)
1511:          ENDDO1511:          ENDDO
1512:          IF (INITIALDIST) THEN1512:          IF (INITIALDIST) THEN
1513:             ALLPAIRS(1:(NMIN*(NMIN-1))/2)=-DISBOUND1513:             ALLPAIRS(1:(NMIN*(NMIN-1)/2))=-DISBOUND
1514:             PRINT '(A,I12)','nmin*(nmin-1)/2=',(NMIN*(NMIN-1))/2 
1515:             INQUIRE(FILE='allpairs',EXIST=YESNO)1514:             INQUIRE(FILE='allpairs',EXIST=YESNO)
1516:             IF (YESNO) THEN1515:             IF (YESNO) THEN
1517:                LUNIT=GETUNIT()1516:                LUNIT=GETUNIT()
1518:                OPEN(UNIT=LUNIT,FILE='allpairs',STATUS='OLD')1517:                OPEN(UNIT=LUNIT,FILE='allpairs',STATUS='OLD')
1519:                DO J1=1,(NMIN*(NMIN-1))/21518:                READ(LUNIT,*) ALLPAIRS(1:(NMIN*(NMIN-1)/2))
1520:                   READ(LUNIT,*) ALLPAIRS(J1) 
1521:                ENDDO 
1522:                CLOSE(LUNIT)1519:                CLOSE(LUNIT)
1523:                PRINT '(A,I8)','setup> Pair distance values read'1520:                PRINT '(A,I8)','setup> Pair distance values read'
1524:             ELSE1521:             ELSE
1525:                CALL GETMETRIC(1,NMIN)1522:                CALL GETMETRIC(1,NMIN)
1526:                LUNIT=GETUNIT()1523:                LUNIT=GETUNIT()
1527:                OPEN(UNIT=LUNIT,FILE='allpairs',STATUS='UNKNOWN')1524:                OPEN(UNIT=LUNIT,FILE='allpairs',STATUS='UNKNOWN')
1528:                WRITE(LUNIT,'(G20.10)') ALLPAIRS(1:(NMIN*(NMIN-1))/2)1525:                WRITE(LUNIT,'(G20.10)') ALLPAIRS(1:(NMIN*(NMIN-1))/2)
1529:                CLOSE(LUNIT)1526:                CLOSE(LUNIT)
1530:             ENDIF1527:             ENDIF
1531:          ELSE1528:          ELSE


legend
Lines Added 
Lines changed
 Lines Removed

hdiff - version: 2.1.0