ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/freemol/trunk/src/mengine/src/kstrbnd.c
Revision: 64
Committed: Thu Dec 4 02:03:23 2008 UTC (12 years, 9 months ago) by gilbertke
File size: 9335 byte(s)
Log Message:
fixed bug in routine to recognize delocalized bonds in kbond.c
Line File contents
1 #define EXTERN extern
2 #include "pcwin.h"
3 #include "pcmod.h"
4
5 #include "angles.h"
6 #include "bonds_ff.h"
7 #include "field.h"
8
9 int find_bond(int, int);
10 void numeral(int, char *, int);
11
12 EXTERN struct t_crossterm_k {
13 int nangang, nstrbnd, nstrtor;
14 int ang_ang[MAXAA], stbnindex[MAXSTBN] ;
15 char str_tor[MAXSTRTOR][7], stbn[MAXSTBN][10];
16 float aacon[MAXAA][3], stbncon[MAXSTBN][3], str_torcon[MAXSTRTOR];
17 } crossterm_k;
18
19 struct t_strbnd {
20 int nstrbnd, isb[MAXANG][3];
21 float ksb1[MAXANG], ksb2[MAXANG];
22 } strbnd;
23 struct stbndat {
24 int irow, krow, lrow;
25 float con1, con2;
26 } stbn[] =
27 {
28 0, 1, 0, 0.15, 0.15,
29 0, 1, 1, 0.10, 0.30,
30 0, 1, 2, 0.05, 0.35,
31 0, 1, 3, 0.05, 0.35,
32 0, 1, 4, 0.05, 0.35,
33 0, 2, 0, 0.00, 0.00,
34 0, 2, 1, 0.00, 0.15,
35 0, 2, 2, 0.00, 0.15,
36 0, 2, 3, 0.00, 0.15,
37 0, 2, 4, 0.00, 0.15,
38 1, 1, 1, 0.30, 0.30,
39 1, 1, 2, 0.30, 0.50,
40 1, 1, 3, 0.30, 0.50,
41 1, 1, 4, 0.30, 0.50,
42 2, 1, 2, 0.50, 0.50,
43 2, 1, 3, 0.50, 0.50,
44 2, 1, 4, 0.50, 0.50,
45 3, 1, 3, 0.50, 0.50,
46 3, 1, 4, 0.50, 0.50,
47 4, 1, 4, 0.50, 0.50,
48 1, 2, 1, 0.30, 0.30,
49 1, 2, 2, 0.25, 0.25,
50 1, 2, 3, 0.25, 0.25,
51 1, 2, 4, 0.25, 0.25,
52 2, 2, 2, 0.25, 0.25,
53 2, 2, 3, 0.25, 0.25,
54 2, 2, 4, 0.25, 0.25,
55 3, 2, 3, 0.25, 0.25,
56 3, 2, 4, 0.25, 0.25,
57 4, 2, 4, 0.25, 0.25 };
58
59 EXTERN int Missing_constants;
60 //EXTERN FILE *errfile;
61
62 void kstrbnd()
63 {
64 int i, j, iz, nb1, nb2, irow, krow, lrow;
65 int ia, ib, ic, itb, ita, itc, ierr;
66 int k, nstbn, index;
67 char pa[4],pb[4],pc[4],pt[10],pt1[10];
68
69 strbnd.nstrbnd = 0;
70 nb1 = -1;
71 nb2 = -1;
72
73 for (i=0; i < angles.nang; i++)
74 {
75 ierr = FALSE;
76 ia = angles.i13[i][0];
77 ib = angles.i13[i][1];
78 ic = angles.i13[i][2];
79 ita = atom[ia].type;
80 itb = atom[ib].type;
81 itc = atom[ic].type;
82 if (field.type == MMX && (atom[ia].type >= 300 || atom[ib].type >= 300 || atom[ic].type >= 300))
83 goto L_10;
84 if (field.type == MMX && (atom[ia].type >= 20 || atom[ic].type >= 20))
85 goto L_10;
86 if (field.type == MM3 && (atom[ia].type >= 300 || atom[ib].type >= 300 || atom[ic].type >= 300) )
87 goto L_10;
88 nstbn = -1;
89 numeral(ita,pa,3);
90 numeral(itb,pb,3);
91 numeral(itc,pc,3);
92 if (ita <= itc)
93 {
94 strcpy(pt,pa); strcat(pt,pb); strcat(pt,pc);
95 }else
96 {
97 strcpy(pt,pc); strcat(pt,pb); strcat(pt,pa);
98 }
99 itb = atom[ib].tclass;
100 numeral(itb,pb,3);
101 strcpy(pt1," 0"); strcat(pt1,pb); strcat(pt1," 0");
102 for (j=0; j < crossterm_k.nstrbnd; j++)
103 {
104 if (strcmp(crossterm_k.stbn[j],pt) == 0 || strcmp(crossterm_k.stbn[j],pt1) == 0 )
105 {
106 nstbn = j;
107 ierr = TRUE;
108 break;
109 }
110 }
111 for (j=0; j < bonds_ff.nbnd; j++)
112 {
113 iz = find_bond(ia, ib);
114 if (iz != -1)
115 {
116 nb1 = iz;
117 break;
118 }
119 }
120 for (j=0; j < bonds_ff.nbnd; j++)
121 {
122 iz = find_bond(ib, ic);
123 if (iz != -1)
124 {
125 nb2 = iz;
126 break;
127 }
128 }
129 if (angles.index[i] == 0)
130 index = 0;
131 else if (angles.index[i] == 1)
132 {
133 if (bonds_ff.index[nb1] == 1)
134 index = 1;
135 else
136 index = 2;
137 } else if (angles.index[i] == 2)
138 index = 3;
139 else if (angles.index[i] == 3)
140 index = 5;
141 else if (angles.index[i] == 4)
142 index = 4;
143 else if (angles.index[i] == 5)
144 {
145 if (bonds_ff.index[nb1] == 1)
146 index = 6;
147 else
148 index = 7;
149 } else if (angles.index[i] == 6)
150 index = 8;
151 else if (angles.index[i] == 7)
152 {
153 if (bonds_ff.index[nb1] == 1)
154 index = 9;
155 else
156 index = 10;
157 } else if (angles.index[i] == 8)
158 index = 11;
159 if (crossterm_k.stbnindex[nstbn] != index)
160 nstbn = -1;
161
162 if (field.type == MMFF94)
163 {
164 ierr = TRUE;
165 if (angles.anat[i] < 175.0)
166 {
167 if (nstbn == -1) // constant not found assign constants base on periodic table
168 {
169 if (atom[ia].atomnum <= 1)
170 irow = 0;
171 else if (atom[ia].atomnum <= 10)
172 irow = 1;
173 else if (atom[ia].atomnum <= 18)
174 irow = 2;
175 else if (atom[ia].atomnum <= 36)
176 irow = 3;
177 else if (atom[ia].atomnum <= 54)
178 irow = 4;
179
180 if (atom[ib].atomnum <= 1)
181 krow = 0;
182 else if (atom[ib].atomnum <= 10)
183 krow = 1;
184 else if (atom[ib].atomnum <= 18)
185 krow = 2;
186 else if (atom[ib].atomnum <= 36)
187 krow = 3;
188 else if (atom[ib].atomnum <= 54)
189 krow = 4;
190
191 if (atom[ic].atomnum <= 1)
192 lrow = 0;
193 else if (atom[ic].atomnum <= 10)
194 lrow = 1;
195 else if (atom[ic].atomnum <= 18)
196 lrow = 2;
197 else if (atom[ic].atomnum <= 36)
198 lrow = 3;
199 else if (atom[ic].atomnum <= 54)
200 lrow = 4;
201
202 for (j=0; j < 30; j++)
203 {
204 if ( irow == stbn[j].irow && krow == stbn[j].krow && lrow == stbn[j].lrow )
205 {
206 strbnd.isb[strbnd.nstrbnd][0] = i;
207 strbnd.isb[strbnd.nstrbnd][1] = nb1;
208 strbnd.isb[strbnd.nstrbnd][2] = nb2;
209 strbnd.ksb1[strbnd.nstrbnd] = stbn[j].con1;
210 strbnd.ksb2[strbnd.nstrbnd] = stbn[j].con2;
211 strbnd.nstrbnd++;
212 break;
213 } else if (lrow == stbn[j].irow && krow == stbn[j].krow && irow == stbn[j].lrow)
214 {
215 strbnd.isb[strbnd.nstrbnd][0] = i;
216 strbnd.isb[strbnd.nstrbnd][1] = nb1;
217 strbnd.isb[strbnd.nstrbnd][2] = nb2;
218 strbnd.ksb1[strbnd.nstrbnd] = stbn[j].con2;
219 strbnd.ksb2[strbnd.nstrbnd] = stbn[j].con1;
220 strbnd.nstrbnd++;
221 break;
222 }
223 }
224 }else
225 {
226 strbnd.isb[strbnd.nstrbnd][0] = i;
227 strbnd.isb[strbnd.nstrbnd][1] = nb1;
228 strbnd.isb[strbnd.nstrbnd][2] = nb2;
229 if (angles.anat[i] > 175.0)
230 {
231 strbnd.ksb1[strbnd.nstrbnd] = 0.0;
232 strbnd.ksb2[strbnd.nstrbnd] = 0.0;
233 } else
234 {
235 if (ita < itc)
236 {
237 strbnd.ksb1[strbnd.nstrbnd] = crossterm_k.stbncon[nstbn][0];
238 strbnd.ksb2[strbnd.nstrbnd] = crossterm_k.stbncon[nstbn][1];
239 } else
240 {
241 strbnd.ksb1[strbnd.nstrbnd] = crossterm_k.stbncon[nstbn][1];
242 strbnd.ksb2[strbnd.nstrbnd] = crossterm_k.stbncon[nstbn][0];
243 }
244 }
245 strbnd.nstrbnd++;
246 }
247 }
248 } else
249 {
250 k = 0;
251 if ( atom[ia].atomnum <= 1)
252 {
253 k++;
254 if (crossterm_k.stbncon[nstbn][2] == 0.0)
255 nb1 = -1;
256 }
257 if (atom[ic].atomnum <= 1)
258 {
259 k++;
260 if (crossterm_k.stbncon[nstbn][2] == 0.0)
261 nb2 = -1;
262 }
263 if (crossterm_k.stbncon[nstbn][k] != 0.0)
264 {
265 strbnd.isb[strbnd.nstrbnd][0] = i;
266 strbnd.isb[strbnd.nstrbnd][1] = nb1;
267 strbnd.isb[strbnd.nstrbnd][2] = nb2;
268 strbnd.ksb1[strbnd.nstrbnd] = crossterm_k.stbncon[nstbn][k];
269 strbnd.ksb2[strbnd.nstrbnd] = crossterm_k.stbncon[nstbn][k];
270 strbnd.nstrbnd++;
271 }
272 }
273 // if (ierr == FALSE)
274 // {
275 // Missing_constants = TRUE;
276 // fprintf(errfile,"Missing strbnd constants for %d %d %d \n",ia,ib,ic);
277 // }
278 L_10:
279 continue;
280 }
281 }
282
283