1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
|
SUBROUTINE PATHS
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
INCLUDE 'SIZES'
COMMON /PATH / LATOM,LPARAM,REACT(200)
COMMON /GEOVAR/ NVAR, LOC(2,MAXPAR), IDUMY, XPARAM(MAXPAR)
COMMON /KEYWRD/ KEYWRD
C ***** Modified by Jiro Toyoda at 1994-05-25 *****
C COMMON /TIME / TIME0
COMMON /TIMEC / TIME0
C ***************************** at 1994-05-25 *****
COMMON /GEOM / GEO(3,NUMATM), XCOORD(3,NUMATM)
COMMON /ALPARM/ ALPARM(3,MAXPAR),X0, X1, X2, ILOOP
************************************************************************
*
* PATH FOLLOWS A REACTION COORDINATE. THE REACTION COORDINATE IS ON
* ATOM LATOM, AND IS A DISTANCE IF LPARAM=1,
* AN ANGLE IF LPARAM=2,
* AN DIHEDRALIF LPARAM=3.
*
************************************************************************
DIMENSION GD(MAXPAR),XLAST(MAXPAR),MDFP(20),XDFP(20)
CHARACTER*241 KEYWRD
CHARACTER*10 TYPE(3)
SAVE TYPE, GD, XLAST, MDFP, XDFP
DATA TYPE / 'ANGSTROMS ','DEGREES ','DEGREES '/
ILOOP=1
IF(INDEX(KEYWRD,'RESTAR') .NE. 0) THEN
MDFP(9)=0
CALL DFPSAV(TOTIME,XPARAM,GD,XLAST,FUNCT1,MDFP,XDFP)
WRITE(6,'(//10X,'' RESTARTING AT POINT'',I3)')ILOOP
ENDIF
IF(ILOOP.GT.1) GOTO 10
WRITE(6,'('' ABOUT TO ENTER FLEPO FROM PATH'')')
TIME0=SECOND()
CALL FLEPO(XPARAM,NVAR,FUNCT)
WRITE(6,'('' OPTIMIZED VALUES OF PARAMETERS, INITIAL POINT'')')
CALL WRITMO(TIME0,FUNCT)
TIME0=SECOND()
10 CONTINUE
IF(ILOOP.GT.2) GOTO 40
GEO(LPARAM,LATOM)=REACT(2)
IF(ILOOP.EQ.1) THEN
X0=REACT(1)
X1=X0
X2=REACT(2)
IF(X2.LT. -100.D0) STOP
DO 20 I=1,NVAR
ALPARM(2,I)=XPARAM(I)
20 ALPARM(1,I)=XPARAM(I)
ILOOP=2
ENDIF
CALL FLEPO(XPARAM,NVAR,FUNCT)
RNORD=REACT(2)
IF(LPARAM.GT.1) RNORD=RNORD*57.29577951D0
WRITE(6,'(1X,16(''*****'')//17X,''REACTION COORDINATE = ''
1,F12.4,2X,A10,19X//1X,16(''*****''))')RNORD,TYPE(LPARAM)
CALL WRITMO(TIME0,FUNCT)
TIME0=SECOND()
DO 30 I=1,NVAR
30 ALPARM(3,I)=XPARAM(I)
C
C NOW FOR THE MAIN INTERPOLATION ROUTE
C
IF(ILOOP.EQ.2)ILOOP=3
40 CONTINUE
LPR=ILOOP
DO 110 ILOOP = LPR,100
C
IF(REACT(ILOOP).LT. -100.D0) RETURN
C
RNORD=REACT(ILOOP)
IF(LPARAM.GT.1) RNORD=RNORD*57.29577951D0
WRITE(6,'(1X,16(''*****'')//19X,''REACTION COORDINATE = ''
1,F12.4,2X,A10,19X//1X,16(''*****''))')RNORD,TYPE(LPARAM)
C
X3=REACT(ILOOP)
C3=(X0**2-X1**2)*(X1-X2)-(X1**2-X2**2)*(X0-X1)
C WRITE(6,'('' C3:'',F13.7)')C3
IF (ABS(C3) .LT. 1.D-8) THEN
C
C WE USE A LINEAR INTERPOLATION
C
CC1=0.D0
CC2=0.D0
ELSE
C WE DO A QUADRATIC INTERPOLATION
C
CC1=(X1-X2)/C3
CC2=(X0-X1)/C3
ENDIF
CB1=1.D0/(X1-X2)
CB2=(X1**2-X2**2)*CB1
C
C NOW TO CALCULATE THE INTERPOLATED COORDINATES
C
DO 50 I=1,NVAR
DELF0=ALPARM(1,I)-ALPARM(2,I)
DELF1=ALPARM(2,I)-ALPARM(3,I)
ACONST = CC1*DELF0-CC2*DELF1
BCONST = CB1*DELF1-ACONST*CB2
CCONST = ALPARM(3,I) - BCONST*X2 - ACONST*X2**2
XPARAM(I)=CCONST+BCONST*X3+ACONST*X3**2
ALPARM(1,I)=ALPARM(2,I)
50 ALPARM(2,I)=ALPARM(3,I)
C
C NOW TO CHECK THAT THE GUESSED GEOMETRY IS NOT TOO ABSURD
C
DO 60 I=1,NVAR
60 IF(ABS(XPARAM(I)-ALPARM(3,I)) .GT. 0.2) GOTO 70
GOTO 90
70 WRITE(6,'('' GEOMETRY TOO UNSTABLE FOR EXTRAPOLATION TO BE USED
1''/ ,'' - THE LAST GEOMETRY IS BEING USED TO START THE NEXT''
2,'' CALCULATION'')')
DO 80 I=1,NVAR
80 XPARAM(I)=ALPARM(3,I)
90 CONTINUE
X0=X1
X1=X2
X2=X3
GEO(LPARAM,LATOM)=REACT(ILOOP)
CALL FLEPO(XPARAM,NVAR,FUNCT)
CALL WRITMO(TIME0,FUNCT)
TIME0=SECOND()
DO 100 I=1,NVAR
100 ALPARM(3,I)=XPARAM(I)
110 CONTINUE
END
|
