ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/freemol/trunk/src/mengine/src/ktorsion.c
Revision: 111
Committed: Thu Mar 12 14:12:33 2009 UTC (11 years, 4 months ago) by gilbertke
File size: 16408 byte(s)
Log Message:
more atom localization
Line File contents
1 #define EXTERN extern
2
3 #include "pcwin.h"
4 #include "torsions.h"
5 #include "atom_k.h"
6
7 int isbond(int,int);
8 int is_delocalbond(int,int);
9 void numeral(int,char *,int);
10 void four(char *,char *,char *,char *, char *);
11 int is_ring42(int,int);
12 int is_ring54(int, int, int, int);
13 void message_alert(char *, char *);
14 int have_ring3(void);
15 int have_ring4(void);
16 int have_ring5(void);
17 int have_ring6(void);
18 int get_field(void);
19 void ktorsion(int *type,int *tclass,int **iat,int **bo);
20
21 EXTERN struct t_allene {
22 int nallene, ntor[10];
23 } allene;
24
25
26 EXTERN struct t_torkn1 {
27 int use_tor4, use_tor5;
28 int ntor, ntor4, ntor5, ntordel, torindex[MAXTORDEL];
29 char kv[MAXTORCONST][13], kv4[MAXTOR4CONST][13], kv5[MAXTOR5CONST][13];
30 char kvdel[MAXTORDEL][13];
31 float tv1[MAXTORCONST], tv2[MAXTORCONST], tv3[MAXTORCONST];
32 float tv4[MAXTORCONST], tv5[MAXTORCONST], tv6[MAXTORCONST];
33 int phase1[MAXTORCONST], phase2[MAXTORCONST], phase3[MAXTORCONST];
34 int phase4[MAXTORCONST], phase5[MAXTORCONST], phase6[MAXTORCONST];
35 float tv41[MAXTOR4CONST], tv42[MAXTOR4CONST], tv43[MAXTOR4CONST];
36 int phase41[MAXTOR4CONST], phase42[MAXTOR4CONST], phase43[MAXTOR4CONST];
37 float tv51[MAXTOR5CONST], tv52[MAXTOR5CONST], tv53[MAXTOR5CONST];
38 int phase51[MAXTOR5CONST], phase52[MAXTOR5CONST], phase53[MAXTOR5CONST];
39 float tvdel1[MAXTORDEL], tvdel2[MAXTORDEL], tvdel3[MAXTORDEL];
40 int phasedel1[MAXTORDEL], phasedel2[MAXTORDEL], phasedel3[MAXTORDEL];
41 } torkn1;
42
43 EXTERN int Missing_constants;
44
45 void ktorsion(int *type,int *tclass,int **iat,int **bo)
46 {
47 int i, j, ierr, itor, k,field;
48 int ia, ib, ic, id;
49 int ita, itb, itc, itd;
50 int cl_a, cl_b, cl_c, cl_d;
51 int use_ring4, use_ring5;
52 char izero[4];
53 char pa[4],pb[4],pc[4],pd[4],pt[13], kv1[13];
54 char k1[13], k2[13], k3[13], k4[13], k5[13], k6[13], k7[13], k8[13], k9[13], k10[13],
55 k11[13], k12[13], k13[13], k14[13], k15[13];
56
57 itor = 0;
58 use_ring4 = FALSE;
59 if ( (have_ring4()) && torkn1.ntor4 > 0)
60 use_ring4 = TRUE;
61
62 use_ring5 = FALSE;
63 if ( (have_ring5()) && torkn1.ntor5 > 0)
64 use_ring5 = TRUE;
65 field = get_field();
66 for (i=0; i < torsions.ntor; i++)
67 {
68 ia = torsions.i14[i][0];
69 ib = torsions.i14[i][1];
70 ic = torsions.i14[i][2];
71 id = torsions.i14[i][3];
72
73 ita = type[ia];
74 itb = type[ib];
75 itc = type[ic];
76 itd = type[id];
77
78 cl_a = tclass[ia];
79 cl_b = tclass[ib];
80 cl_c = tclass[ic];
81 cl_d = tclass[id];
82
83 if (field == MMX)
84 {
85 if (ita == 40)
86 ita = 2;
87 if (itb == 40)
88 itb = 2;
89 if (itc == 40)
90 itc = 2;
91 if (itd == 40)
92 itd = 2;
93
94 if (ita == 56 && itb == 56 && itc == 56 && itd == 56)
95 {
96 ita = 1; itb = 1; itc = 1; itd = 1;
97 }
98 }
99 numeral(ita,pa,3);
100 numeral(itb,pb,3);
101 numeral(itc,pc,3);
102 numeral(itd,pd,3);
103
104 if (itb < itc )
105 four(pa,pb, pc, pd, pt);
106 else if (itb > itc)
107 four(pd,pc, pb, pa, pt);
108 else if (ita < itd)
109 four(pa,pb, pc, pd, pt);
110 else
111 four(pd,pc, pb, pa, pt);
112 strcpy(izero," 0");
113
114 if (field == MMFF94)
115 {
116
117 numeral(atom_k.tclass1[ita],pa,3);
118 numeral(atom_k.tclass[itb],pb,3);
119 numeral(atom_k.tclass[itc],pc,3);
120 numeral(atom_k.tclass1[itd],pd,3);
121
122
123 if (atom_k.tclass[itb] < atom_k.tclass[itc] )
124 {
125 four(pa, pb, pc, izero,k2);
126 four(izero ,pb, pc, pd,k3);
127 } else if (atom_k.tclass[itb] > atom_k.tclass[itc])
128 {
129 four(izero,pc, pb, pa, k2);
130 four(pd,pc, pb, izero, k3);
131 } else if (atom_k.tclass1[ita] < atom_k.tclass1[itd])
132 {
133 four(pa,pb, pc, izero, k2);
134 four(izero,pb, pc, pd, k3);
135 } else
136 {
137 four(izero,pc, pb, pa, k2);
138 four(pd,pc, pb, izero, k3);
139 }
140 } else
141 {
142 four( izero, pc, pb, pa, k1 );
143 four( pd, pc, pb, izero, k2 );
144 four( izero, pb, pc, pd, k3 );
145 }
146
147 numeral(ita,pa,3);
148 numeral(itb,pb,3);
149 numeral(itc,pc,3);
150 numeral(itd,pd,3);
151
152 four( pa, pb, pc, izero, k4 );
153 four( izero, pc, pb, izero, k5 );
154 four( izero, pb, pc, izero, k6 );
155 four( pd, pc, izero, izero, k7 );
156 four( izero, izero, pc,pd, k8 );
157
158 four( izero, izero, pb, pa, k9 );
159 four( pa, pb, izero, izero, k10 );
160 four( izero, izero, pc, izero, k11 );
161 four( izero, pc, izero, izero, k12 );
162 four( izero, izero, pb, izero, k13 );
163 four( izero, pb, izero, izero, k14 );
164
165 numeral(cl_a,pa,3);
166 numeral(cl_b,pb,3);
167 numeral(cl_c,pc,3);
168 numeral(cl_d,pd,3);
169 if (itb < itc )
170 four(pa,pb, pc, pd, k15);
171 else if (itb > itc)
172 four(pd,pc, pb, pa, k15);
173 else if (ita < itd)
174 four(pa,pb, pc, pd, k15);
175 else
176 four(pd,pc, pb, pa, k15);
177
178 ierr = FALSE;
179 // check for four membered rings
180 if ( isbond(ia,id) )
181 {
182 for(j=0; j < torkn1.ntor4; j++)
183 {
184 strcpy(kv1,torkn1.kv4[j]);
185 if (strcmp(pt,kv1) == 0 || strcmp(k1,kv1) == 0 || strcmp(k2,kv1) == 0 || strcmp(k3,kv1) == 0 ||
186 strcmp(k4,kv1) == 0 || strcmp(k5,kv1) == 0 || strcmp(k6,kv1) == 0 || strcmp(k7,kv1) == 0 ||
187 strcmp(k8,kv1) == 0 || strcmp(k9,kv1) == 0 || strcmp(k10,kv1) == 0 || strcmp(k11,kv1) == 0 ||
188 strcmp(k12,kv1) == 0 || strcmp(k13,kv1) == 0 || strcmp(k14,kv1) == 0 )
189 {
190 torsions.v1[i] = torkn1.tv41[j];
191 torsions.ph1[i] = torkn1.phase41[j];
192 torsions.v2[i] = torkn1.tv42[j];
193 torsions.ph2[i] = torkn1.phase42[j];
194 torsions.v3[i] = torkn1.tv43[j];
195 torsions.ph3[i] = torkn1.phase43[j];
196 ierr = TRUE;
197 goto L_10;
198 }
199 }
200 }
201 // delocalized torsions
202 if ( is_delocalbond(ib,ic) )
203 {
204 for(j=0; j < torkn1.ntordel; j++)
205 {
206 strcpy(kv1,torkn1.kvdel[j]);
207 if (strcmp(pt,kv1) == 0 || strcmp(k1,kv1) == 0 || strcmp(k2,kv1) == 0 || strcmp(k3,kv1) == 0 ||
208 strcmp(k4,kv1) == 0 || strcmp(k5,kv1) == 0 || strcmp(k6,kv1) == 0 || strcmp(k7,kv1) == 0 ||
209 strcmp(k8,kv1) == 0 || strcmp(k9,kv1) == 0 || strcmp(k10,kv1) == 0 || strcmp(k11,kv1) == 0 ||
210 strcmp(k12,kv1) == 0 || strcmp(k13,kv1) == 0 || strcmp(k14,kv1) == 0 )
211 {
212 torsions.v1[i] = torkn1.tvdel1[j];
213 torsions.ph1[i] = torkn1.phasedel1[j];
214 torsions.v2[i] = torkn1.tvdel2[j];
215 torsions.ph2[i] = torkn1.phasedel2[j];
216 torsions.v3[i] = torkn1.tvdel3[j];
217 torsions.ph3[i] = torkn1.phasedel3[j];
218 ierr = TRUE;
219 goto L_10;
220 }
221 }
222 }
223 // check for five membered rings
224 if (is_ring54(ia,ib,ic,id) )
225 {
226 for(j=0; j < torkn1.ntor5; j++)
227 {
228 strcpy(kv1,torkn1.kv5[j]);
229 if (strcmp(pt,kv1) == 0 || strcmp(k1,kv1) == 0 || strcmp(k2,kv1) == 0 || strcmp(k3,kv1) == 0 ||
230 strcmp(k4,kv1) == 0 || strcmp(k5,kv1) == 0 || strcmp(k6,kv1) == 0 || strcmp(k7,kv1) == 0 ||
231 strcmp(k8,kv1) == 0 || strcmp(k9,kv1) == 0 || strcmp(k10,kv1) == 0 || strcmp(k11,kv1) == 0 ||
232 strcmp(k12,kv1) == 0 || strcmp(k13,kv1) == 0 || strcmp(k14,kv1) == 0 )
233 {
234 torsions.v1[i] = torkn1.tv51[j];
235 torsions.ph1[i] = torkn1.phase51[j];
236 torsions.v2[i] = torkn1.tv52[j];
237 torsions.ph2[i] = torkn1.phase52[j];
238 torsions.v3[i] = torkn1.tv53[j];
239 torsions.ph3[i] = torkn1.phase53[j];
240 ierr = TRUE;
241 goto L_10;
242 }
243 }
244 }
245 // check regular parameters
246 if (field == MMFF94)
247 {
248 for(j=0; j < torkn1.ntor; j++) // specific constants
249 {
250 strcpy(kv1,torkn1.kv[j]);
251 if (strcmp(pt,kv1) == 0)
252 {
253 torsions.v1[i] = torkn1.tv1[j];
254 torsions.ph1[i] = torkn1.phase1[j];
255 torsions.v2[i] = torkn1.tv2[j];
256 torsions.ph2[i] = torkn1.phase2[j];
257 torsions.v3[i] = torkn1.tv3[j];
258 torsions.ph3[i] = torkn1.phase3[j];
259 torsions.v4[i] = torkn1.tv4[j];
260 torsions.ph4[i] = torkn1.phase4[j];
261 torsions.v5[i] = torkn1.tv5[j];
262 torsions.ph5[i] = torkn1.phase5[j];
263 torsions.v6[i] = torkn1.tv6[j];
264 torsions.ph6[i] = torkn1.phase6[j];
265 goto L_10;
266 }
267 }
268 for(j=0; j < torkn1.ntor; j++) // class 2
269 {
270 strcpy(kv1,torkn1.kv[j]);
271 if (strcmp(k2,kv1) == 0 || strcmp(k3,kv1) == 0 )
272 {
273 torsions.v1[i] = torkn1.tv1[j];
274 torsions.ph1[i] = torkn1.phase1[j];
275 torsions.v2[i] = torkn1.tv2[j];
276 torsions.ph2[i] = torkn1.phase2[j];
277 torsions.v3[i] = torkn1.tv3[j];
278 torsions.ph3[i] = torkn1.phase3[j];
279 torsions.v4[i] = torkn1.tv4[j];
280 torsions.ph4[i] = torkn1.phase4[j];
281 torsions.v5[i] = torkn1.tv5[j];
282 torsions.ph5[i] = torkn1.phase5[j];
283 torsions.v6[i] = torkn1.tv6[j];
284 torsions.ph6[i] = torkn1.phase6[j];
285 goto L_10;
286 }
287 }
288 for(j=0; j < torkn1.ntor; j++) // wild cards
289 {
290 strcpy(kv1,torkn1.kv[j]);
291 if (strcmp(k5,kv1) == 0 || strcmp(k6,kv1) == 0 || strcmp(k7,kv1) == 0 || strcmp(k8,kv1) == 0 )
292 {
293 torsions.v1[i] = torkn1.tv1[j];
294 torsions.ph1[i] = torkn1.phase1[j];
295 torsions.v2[i] = torkn1.tv2[j];
296 torsions.ph2[i] = torkn1.phase2[j];
297 torsions.v3[i] = torkn1.tv3[j];
298 torsions.ph3[i] = torkn1.phase3[j];
299 torsions.v4[i] = torkn1.tv4[j];
300 torsions.ph4[i] = torkn1.phase4[j];
301 torsions.v5[i] = torkn1.tv5[j];
302 torsions.ph5[i] = torkn1.phase5[j];
303 torsions.v6[i] = torkn1.tv6[j];
304 torsions.ph6[i] = torkn1.phase6[j];
305 goto L_10;
306 }
307 }
308 }
309 // check regular specific parameters
310 for(j=0; j < torkn1.ntor; j++)
311 {
312 strcpy(kv1,torkn1.kv[j]);
313 if (strcmp(pt,kv1) == 0)
314 {
315 torsions.v1[i] = torkn1.tv1[j];
316 torsions.ph1[i] = torkn1.phase1[j];
317 torsions.v2[i] = torkn1.tv2[j];
318 torsions.ph2[i] = torkn1.phase2[j];
319 torsions.v3[i] = torkn1.tv3[j];
320 torsions.ph3[i] = torkn1.phase3[j];
321 torsions.v4[i] = torkn1.tv4[j];
322 torsions.ph4[i] = torkn1.phase4[j];
323 torsions.v5[i] = torkn1.tv5[j];
324 torsions.ph5[i] = torkn1.phase5[j];
325 torsions.v6[i] = torkn1.tv6[j];
326 torsions.ph6[i] = torkn1.phase6[j];
327 // check for pi torsions = mm3 does not define pitorsion parameter so we need to pick it up here
328 if (itb == 2 && itc == 2 && field == MMX)
329 {
330 for (k=0; k < MAXIAT; k++)
331 {
332 if (iat[ib][k] == ic)
333 {
334 if (bo[ib][k] == 1 && is_delocalbond(ib,ic) )
335 {
336 torsions.v1[i] = 0.5;
337 torsions.v2[i] = 1.0;
338 torsions.v3[i] = 0.0;
339 ierr = TRUE;
340 goto L_10;
341 break;
342 }
343 }
344 }
345 }
346 goto L_10;
347 break;
348 }
349 }
350 // check regular generalized parameters
351 for(j=0; j < torkn1.ntor; j++)
352 {
353 strcpy(kv1,torkn1.kv[j]);
354 if (strcmp(pt,kv1) == 0 || strcmp(k1,kv1) == 0 || strcmp(k2,kv1) == 0 || strcmp(k3,kv1) == 0 ||
355 strcmp(k4,kv1) == 0 || strcmp(k5,kv1) == 0 || strcmp(k6,kv1) == 0 || strcmp(k7,kv1) == 0 ||
356 strcmp(k8,kv1) == 0 || strcmp(k9,kv1) == 0 || strcmp(k10,kv1) == 0 || strcmp(k11,kv1) == 0 ||
357 strcmp(k12,kv1) == 0 || strcmp(k13,kv1) == 0 || strcmp(k14,kv1) == 0 || strcmp(k15,kv1) == 0)
358 {
359 torsions.v1[i] = torkn1.tv1[j];
360 torsions.ph1[i] = torkn1.phase1[j];
361 torsions.v2[i] = torkn1.tv2[j];
362 torsions.ph2[i] = torkn1.phase2[j];
363 torsions.v3[i] = torkn1.tv3[j];
364 torsions.ph3[i] = torkn1.phase3[j];
365 torsions.v4[i] = torkn1.tv4[j];
366 torsions.ph4[i] = torkn1.phase4[j];
367 torsions.v5[i] = torkn1.tv5[j];
368 torsions.ph5[i] = torkn1.phase5[j];
369 torsions.v6[i] = torkn1.tv6[j];
370 torsions.ph6[i] = torkn1.phase6[j];
371 // check for pi torsions = mm3 does not define pitorsion parameter so we need to pick it up here
372 if (itb == 2 && itc == 2 && field == MMX)
373 {
374 for (k=0; k < MAXIAT; k++)
375 {
376 if (iat[ib][k] == ic)
377 {
378 if (bo[ib][k] == 1 && is_delocalbond(ib,ic) )
379 {
380 torsions.v1[i] = 0.5;
381 torsions.v2[i] = 1.0;
382 torsions.v3[i] = 0.0;
383 ierr = TRUE;
384 goto L_10;
385 break;
386 }
387 }
388 }
389 }
390 goto L_10;
391 break;
392 }
393 }
394 // missing parameters
395 if (ierr == FALSE)
396 {
397 // Missing_constants = TRUE;
398 fprintf(pcmlogfile,"Torsion constants missing: angle: %d %d %d %d types: %d %d %d %d\n",ia,ib,ic,id,
399 ita, itb, itc, itd);
400 }
401 L_10:
402 continue;
403 }
404
405 if (allene.nallene > 0)
406 {
407 for (i = 0; i < allene.nallene; i++)
408 {
409 torsions.v1[allene.ntor[i]] = 0;
410 torsions.v2[allene.ntor[i]] = -11.5;
411 torsions.v3[allene.ntor[i]] = 0;
412 }
413 }
414 }
415 /* -------------------------------------------- */
416 void four(char *pa, char *pb, char *pc, char *pd, char *pt)
417 {
418 strcpy(pt,pa);
419 strcat(pt,pb);
420 strcat(pt,pc);
421 strcat(pt,pd);
422 }