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