ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/freemol/trunk/smi23d/src/mengine/kangle.c
Revision: 3
Committed: Mon Jun 9 21:38:26 2008 UTC (11 years, 4 months ago) by tjod
File size: 24249 byte(s)
Log Message:
test

Line File contents
1 #define EXTERN extern
2 #include "pcwin.h"
3 #include "pcmod.h"
4
5 #include "angles.h"
6 #include "rings.h"
7 #include "field.h"
8 #include "bonds_ff.h"
9 #include "atom_k.h"
10
11 int is_ring43(int, int, int);
12 int is_ring53(int, int, int);
13 int is_ring63(int, int, int);
14 int is_delocalbond(int, int);
15 int is_TS(int, int, int);
16 int isbond(int,int);
17 int find_bond(int,int);
18 void transtheta(int *, float *, float *, char *, char *, char *, char *);
19 void ts_four(int *, float *, float *,char *, char *,char *,char *);
20 long ipack(int ,int ,int ,int);
21 long kang(int,int,int);
22 void numeral(int,char *,int);
23 float get_general_angle(int);
24
25
26 EXTERN struct t_ts_bondorder {
27 float fbnd[15];
28 } ts_bondorder;
29
30 EXTERN struct t_angk1 {
31 int use_ang3, use_ang4, use_ang5;
32 int nang, nang3, nang4, nang5;
33 int ndel, ndel3, ndel4;
34 char ktype[MAXANGCONST][10],ktype3[MAXANG3CONST][10],ktype4[MAXANG4CONST][10],ktype5[MAXANG5CONST][10];
35 char kdel[MAXANGDEL][10],kdel3[MAXANG3DEL][10],kdel4[MAXANG4DEL][10];
36 float con[MAXANGCONST], ang[MAXANGCONST][3];
37 float con3[MAXANG3CONST], ang3[MAXANG3CONST][3];
38 float con4[MAXANG4CONST], ang4[MAXANG4CONST][3];
39 float con5[MAXANG5CONST], ang5[MAXANG5CONST][3];
40 float condel[MAXANGDEL], angdel[MAXANGDEL][3];
41 float condel3[MAXANG3DEL], angdel3[MAXANG3DEL][3];
42 float condel4[MAXANG4DEL], angdel4[MAXANG4DEL][3];
43 } angk1;
44
45 EXTERN struct t_angkp1 {
46 int npiang;
47 char kpa[MAXPIANGCONST][10];
48 float pacon[MAXPIANGCONST], panat[MAXPIANGCONST][3], piobp[MAXPIANGCONST];
49 } angkp1;
50
51 EXTERN struct t_angf {
52 int use_angf, nfang;
53 char kftype[MAXANGCONST][10];
54 float fcon[MAXANGCONST],fc0[MAXANGCONST],fc1[MAXANGCONST],fc2[MAXANGCONST];
55 float fc3[MAXANGCONST],fc4[MAXANGCONST],fc5[MAXANGCONST],fc6[MAXANGCONST];
56 } angf;
57
58 EXTERN int Missing_constants;
59
60 void kangle()
61 {
62 long int pi_mask, aromatic_mask;
63 int i, j, k, izero, ierr, nh, nRc, nligand, jji;
64 int ia, ib, ic, ie;
65 int iaa, ibb, icc;
66 int cl_a, cl_b, cl_c;
67 int nbo1, nb1, nb2;
68 char pa[4],pb[4],pc[4],pt[10], pt1[10], pt2[10],pt3[10], pt4[10], pt5[10], pt6[10], pt7[10];
69 char zero[4];
70 float bkan, at;
71
72 pi_mask = (1L << PI_MASK);
73 aromatic_mask = (1 << AROMATIC_MASK);
74
75 if( angles.nang != 0 )
76 {
77 for( i = 0; i < angles.nang; i++ )
78 {
79 ia = angles.i13[i][0];
80 ib = angles.i13[i][1];
81 ic = angles.i13[i][2];
82 iaa = atom[ia].type;
83 ibb = atom[ib].type;
84 icc = atom[ic].type;
85 cl_a = atom[ia].tclass;
86 cl_b = atom[ib].tclass;
87 cl_c = atom[ic].tclass;
88
89 /* get bond index each bond for MMFF94 */
90 nb1 = find_bond(ia,ib);
91 nb2 = find_bond(ib,ic);
92
93 /* special metal cases */
94 if( iaa >= 300)
95 cl_a = 300;
96
97 if( ibb >= 300)
98 cl_b = 300;
99
100 if( icc >= 300)
101 cl_c = 300;
102
103 /* special MMX cases */
104 if ( field.type == MMX)
105 {
106 if( iaa == 40 )
107 iaa = 2;
108
109 if( ibb == 40 )
110 ibb = 2;
111
112 if( icc == 40 )
113 icc = 2;
114
115 }
116 numeral(iaa,pa,3);
117 numeral(ibb,pb,3);
118 numeral(icc,pc,3);
119 if( iaa <= icc )
120 {
121 strcpy(pt,pa);
122 strcat(pt,pb);
123 strcat(pt,pc);
124 }
125 if( iaa > icc )
126 {
127 strcpy(pt,pc);
128 strcat(pt,pb);
129 strcat(pt,pa);
130 }
131 izero = 0;
132 strcpy(zero," 0");
133 strcpy(pt1,pa); strcat(pt1,pb); strcat(pt1,zero);
134 strcpy(pt2,pc); strcat(pt2,pb); strcat(pt2,zero);
135 strcpy(pt3,zero); strcat(pt3,pb); strcat(pt3,pc);
136 strcpy(pt4,zero); strcat(pt4,pb); strcat(pt4,pa);
137 strcpy(pt5,zero); strcat(pt5,pb); strcat(pt5,zero);
138
139 numeral(cl_a,pa,3);
140 numeral(cl_b,pb,3);
141 numeral(cl_c,pc,3);
142 if (cl_a < cl_c)
143 {
144 strcpy(pt6,pa); strcat(pt6,pb); strcat(pt6,pc);
145 }else
146 {
147 strcpy(pt6,pc); strcat(pt6,pb); strcat(pt6,pa);
148 }
149
150 if (field.type == MMFF94)
151 {
152 cl_a = atom_k.tclass1[atom[ia].type];
153 cl_b = atom_k.tclass1[atom[ib].type];
154 cl_c = atom_k.tclass1[atom[ic].type];
155 numeral(cl_a,pa,3);
156 numeral(cl_b,pb,3);
157 numeral(cl_c,pc,3);
158 if (cl_a < cl_c)
159 {
160 strcpy(pt7,pa); strcat(pt7,pb); strcat(pt7,pc);
161 }else
162 {
163 strcpy(pt7,pc); strcat(pt7,pb); strcat(pt7,pa);
164 }
165 }
166 angles.acon[i] = 0.0;
167 angles.anat[i] = 0.0;
168 angles.angtype[i] = HARMONIC;
169
170 ierr = FALSE;
171 /* check TS constants */
172 /* check delocalized angles in MMFF94 cyclopropanes */
173 if (field.type == MMFF94 && angk1.ndel3 > 0)
174 {
175 if (isbond(ia,ic) )
176 {
177 nbo1 = bonds_ff.index[nb1] + bonds_ff.index[nb2];
178 if (nbo1 == 1 || nbo1 == 2)
179 {
180 for(j=0; j < angk1.ndel3; j++)
181 {
182 if (strcmp(angk1.kdel3[j],pt) == 0)
183 {
184 angles.acon[i] = angk1.condel3[j];
185 angles.anat[i] = angk1.angdel3[j][0];
186 if (nbo1 == 1)
187 angles.index[i] = 5;
188 else if (nbo1 == 2)
189 angles.index[i] = 6;
190 ierr = TRUE;
191 break;
192 }
193 }
194 }
195
196 if (ierr == TRUE)
197 goto L_10;
198 }
199 }
200 /* check three membered rings */
201 if ( (rings.nring3 > 0) && (angk1.nang3 > 0) )
202 {
203 if (isbond(ia,ic))
204 {
205 /* search three membered ring parameters */
206 for (j = 0; j < angk1.nang3; j++)
207 {
208 if ( (strcmp(angk1.ktype3[j],pt) == 0) || (strcmp(angk1.ktype3[j],pt1) == 0) ||
209 (strcmp(angk1.ktype3[j],pt2) == 0) || (strcmp(angk1.ktype3[j],pt3) == 0) ||
210 (strcmp(angk1.ktype3[j],pt4) == 0) || (strcmp(angk1.ktype3[j],pt5) == 0))
211 {
212 if ( angk1.ang3[j][1] == 0.0 && angk1.ang3[j][2] == 0)
213 {
214 angles.acon[i] = angk1.con3[j];
215 angles.anat[i] = angk1.ang3[j][0];
216 angles.index[i] = 3;
217 ierr = TRUE;
218 break;
219 } else
220 {
221 nh = 0;
222 for (k=0; k < MAXIAT; k++)
223 {
224 ie = atom[ib].iat[k];
225 if ( (atom[ie].type == 5 || atom[ie].type == 36) && (ie != ia) && (ie != ic))
226 nh++;
227 }
228 angles.acon[i] = angk1.con3[j];
229 angles.anat[i] = angk1.ang3[j][nh];
230 angles.index[i] = 3;
231 ierr = TRUE;
232 break;
233 }
234 }
235 }
236 if (ierr == TRUE)
237 goto L_10;
238 else
239 {
240 Missing_constants = TRUE;
241 }
242 }
243 }
244 /* check delocalized angles in MMFF94 cyclobutanes */
245 if (field.type == MMFF94 && angk1.ndel4 > 0)
246 {
247 nRc = is_ring43(ia, ib, ic);
248 if (nRc == TRUE)
249 {
250 nbo1 = bonds_ff.index[nb1] + bonds_ff.index[nb2];
251 if (nbo1 == 1 || nbo1 == 2)
252 {
253 for(j=0; j < angk1.ndel4; j++)
254 {
255 if (strcmp(angk1.kdel4[j],pt) == 0)
256 {
257 angles.acon[i] = angk1.condel4[j];
258 angles.anat[i] = angk1.angdel4[j][0];
259 if (nbo1 == 1)
260 angles.index[i] = 7;
261 else if (nbo1 == 2)
262 angles.index[i] = 8;
263 ierr = TRUE;
264 break;
265 }
266 }
267 }
268 if (ierr == TRUE)
269 goto L_10;
270 }
271 }
272 /* check four memebered rings */
273 if ( (rings.nring4 > 0) && (angk1.nang4 > 0) )
274 {
275 nRc = is_ring43(ia, ib, ic);
276 if (nRc == TRUE)
277 {
278 /* search four membered ring parameters */
279 for (j = 0; j < angk1.nang4; j++)
280 {
281 if ( (strcmp(angk1.ktype4[j],pt) == 0) || (strcmp(angk1.ktype4[j],pt1) == 0) ||
282 (strcmp(angk1.ktype4[j],pt2) == 0) || (strcmp(angk1.ktype4[j],pt3) == 0) ||
283 (strcmp(angk1.ktype4[j],pt4) == 0) || (strcmp(angk1.ktype4[j],pt5) == 0))
284 {
285 if ( angk1.ang4[j][1] == 0.0 && angk1.ang4[j][2] == 0)
286 {
287 angles.acon[i] = angk1.con4[j];
288 angles.anat[i] = angk1.ang4[j][0];
289 angles.index[i] = 4;
290 ierr = TRUE;
291 break;
292 } else
293 {
294 nh = 0;
295 for (k=0; k < MAXIAT; k++)
296 {
297 ie = atom[ib].iat[k];
298 if ( (atom[ie].type == 5 || atom[ie].type == 36) && (ie != ia) && (ie != ic))
299 nh++;
300 }
301 angles.acon[i] = angk1.con4[j];
302 angles.anat[i] = angk1.ang4[j][nh];
303 angles.index[i] = 4;
304 ierr = TRUE;
305 break;
306 }
307 }
308 }
309 if (ierr == TRUE)
310 goto L_10;
311 else
312 {
313 Missing_constants = TRUE;
314 }
315 }
316 }
317 /* check five membered rings */
318 if ( (rings.nring5 > 0) && (angk1.nang5 > 0) )
319 {
320 nRc = is_ring53(ia,ib, ic);
321 if (nRc == TRUE)
322 {
323 /* search five memebered ring parameters */
324 for (j = 0; j < angk1.nang5; j++)
325 {
326 if ( (strcmp(angk1.ktype5[j],pt) == 0) || (strcmp(angk1.ktype5[j],pt1) == 0) ||
327 (strcmp(angk1.ktype5[j],pt2) == 0) || (strcmp(angk1.ktype5[j],pt3) == 0) ||
328 (strcmp(angk1.ktype5[j],pt4) == 0) || (strcmp(angk1.ktype5[j],pt5) == 0) )
329 {
330 if ( angk1.ang5[j][1] == 0.0 && angk1.ang5[j][2] == 0)
331 {
332 angles.acon[i] = angk1.con5[j];
333 angles.anat[i] = angk1.ang5[j][0];
334 angles.index[i] = 0;
335 ierr = TRUE;
336 break;
337 } else
338 {
339 nh = 0;
340 for (k=0; k < MAXIAT; k++)
341 {
342 ie = atom[ib].iat[k];
343 if ( (atom[ie].type == 5 || atom[ie].type == 36) && (ie != ia) && (ie != ic))
344 nh++;
345 }
346 angles.acon[i] = angk1.con5[j];
347 angles.anat[i] = angk1.ang5[j][nh];
348 angles.index[i] = 0;
349 ierr = TRUE;
350 break;
351 }
352 }
353 }
354 if (ierr == TRUE)
355 goto L_10;
356 // don't fail on missing 5 ring parameters
357 }
358 }
359 /* check delocalized angles in MMFF94 */
360 if (field.type == MMFF94 && angk1.ndel > 0)
361 {
362 if ( bonds_ff.index[nb1] == 1 || bonds_ff.index[nb2] == 1 )
363 {
364 nbo1 = bonds_ff.index[nb1] + bonds_ff.index[nb2];
365 for(j=0; j < angk1.ndel; j++)
366 {
367 if (strcmp(angk1.kdel[j],pt) == 0)
368 {
369 angles.acon[i] = angk1.condel[j];
370 angles.anat[i] = angk1.angdel[j][0];
371 if (nbo1 == 1)
372 angles.index[i] = 1;
373 else if (nbo1 == 2)
374 angles.index[i] = 2;
375 ierr = TRUE;
376 break;
377 }
378 }
379 if (ierr == TRUE)
380 goto L_10;
381 }
382 }
383 // look for fourier angles
384 if (angf.nfang > 0)
385 {
386 for (j=0; j < angf.nfang; j++)
387 {
388 if (strcmp(angf.kftype[j],pt) == 0)
389 {
390 angles.fcon[i] = angf.fcon[j];
391 angles.c0[i] = angf.fc0[j];
392 angles.c1[i] = angf.fc1[j];
393 angles.c2[i] = angf.fc2[j];
394 angles.c3[i] = angf.fc3[j];
395 angles.c4[i] = angf.fc4[j];
396 angles.c5[i] = angf.fc5[j];
397 angles.c6[i] = angf.fc6[j];
398 angles.angtype[i] = FOURIER;
399 ierr = TRUE;
400 break;
401 }
402 }
403 if (ierr == TRUE)
404 goto L_10;
405 }
406
407 /* check regular parameters for specific angles*/
408 for (j = 0; j < angk1.nang; j++)
409 {
410 if ( (strcmp(angk1.ktype[j],pt) == 0) )
411 {
412 if ( angk1.ang[j][1] == 0.0 && angk1.ang[j][2] == 0)
413 {
414 angles.acon[i] = angk1.con[j];
415 angles.anat[i] = angk1.ang[j][0];
416 angles.index[i] = 0;
417 ierr = TRUE;
418 break;
419 } else
420 {
421 nh = 0;
422 for (k=0; k < MAXIAT; k++)
423 {
424 ie = atom[ib].iat[k];
425 if ( (atom[ie].atomnum == 1) && (ie != ia) && (ie != ic))
426 nh++;
427 }
428 angles.acon[i] = angk1.con[j];
429 angles.anat[i] = angk1.ang[j][nh];
430 angles.index[i] = 0;
431 ierr = TRUE;
432 break;
433 }
434 }
435 }
436 if (ierr == TRUE)
437 goto L_10;
438 // did not find any specific look for generalized
439 for (j = 0; j < angk1.nang; j++)
440 {
441 if ( (strcmp(angk1.ktype[j],pt) == 0) || (strcmp(angk1.ktype[j],pt1) == 0) ||
442 (strcmp(angk1.ktype[j],pt2) == 0) || (strcmp(angk1.ktype[j],pt3) == 0) ||
443 (strcmp(angk1.ktype[j],pt4) == 0) || (strcmp(angk1.ktype[j],pt5) == 0) ||
444 (strcmp(angk1.ktype[j],pt6) == 0))
445 {
446 if ( angk1.ang[j][1] == 0.0 && angk1.ang[j][2] == 0)
447 {
448 angles.acon[i] = angk1.con[j];
449 angles.anat[i] = angk1.ang[j][0];
450 angles.index[i] = 0;
451 ierr = TRUE;
452 break;
453 } else
454 {
455 nh = 0;
456 for (k=0; k < MAXIAT; k++)
457 {
458 ie = atom[ib].iat[k];
459 if ( (atom[ie].atomnum == 1) && (ie != ia) && (ie != ic))
460 nh++;
461 }
462 angles.acon[i] = angk1.con[j];
463 angles.anat[i] = angk1.ang[j][nh];
464 angles.index[i] = 0;
465 ierr = TRUE;
466 break;
467 }
468 }
469 }
470 if (ierr == TRUE)
471 goto L_10;
472 /* check MMFF for class1 parameters */
473 if (field.type == MMFF94)
474 {
475 for (j = 0; j < angk1.nang; j++)
476 {
477 if ( (strcmp(angk1.ktype[j],pt7) == 0) )
478 {
479 if ( angk1.ang[j][1] == 0.0 && angk1.ang[j][2] == 0)
480 {
481 angles.acon[i] = angk1.con[j];
482 angles.anat[i] = angk1.ang[j][0];
483 angles.index[i] = 0;
484 ierr = TRUE;
485 break;
486 } else
487 {
488 nh = 0;
489 for (k=0; k < MAXIAT; k++)
490 {
491 ie = atom[ib].iat[k];
492 if ( (atom[ie].atomnum == 1) && (ie != ia) && (ie != ic))
493 nh++;
494 }
495 angles.acon[i] = angk1.con[j];
496 angles.anat[i] = angk1.ang[j][nh];
497 angles.index[i] = 0;
498 ierr = TRUE;
499 break;
500 }
501 }
502 }
503 if (ierr == TRUE)
504 goto L_10;
505 }
506 if (ierr != TRUE)
507 {
508 bkan = get_general_angle(ibb);
509 angles.acon[i] = 0.80;
510 angles.anat[i] = bkan;
511 angles.index[i] = 0;
512 // ==== done with angle lookup =============
513 L_10:
514 if (field.type == MMX && ierr != TRUE)
515 {
516 if (ibb >= 300)
517 {
518 nligand = 0;
519 jji = 0;
520 for ( k=0; k < MAXIAT; k++)
521 {
522 if (atom[ib].iat[k] != 0 && atom[ib].bo[k] != 9)
523 jji++;
524 if (atom[ib].bo[k] == 9)
525 nligand++;
526 }
527 if (jji == 2 && nligand == 0)
528 {
529 angles.acon[i] = .1;
530 angles.anat[i] = 180.;
531 }
532 if ((jji == 3 && nligand == 0) || (jji == 2 && nligand == 1) )
533 {
534 angles.acon[i] = .1;
535 angles.anat[i] = 120.;
536 }
537 if ((jji == 4 && nligand == 0) || (jji == 3 && nligand == 1) || (jji == 2 && nligand == 2))
538 {
539 if (jji == 4 && (atom[ib].type >= 300 && (atom[ib].flags & (1L << SATMET_MASK) &&
540 atom[ib].flags & (1L << GT18e_MASK) && atom[ib].flags & (1L << SQPLAN_MASK) &&
541 !(atom[ib].flags & (1L << LOWSPIN_MASK)))) )
542 {
543 angles.acon[i] = .0001;
544 angles.anat[i] = 90.;
545 } else
546 {
547 angles.acon[i] = .1;
548 angles.anat[i] = 109.5;
549 }
550 }
551
552 }
553 }
554 }
555 }
556 }
557 }
558 // ===================================
559 // generalized angles for MMFF with no hydrogen minimization
560 float get_general_angle(int ia)
561 {
562 float angs[80] = {
563 109.0, 120.0, 120.0, 180.0, 0.000, 109.0, 0.000, 109.0, 120.0, 118.0,
564 0.000, 0.000, 0.000, 0.000, 109.0, 109.0, 120.0, 109.0, 109.0, 90.00,
565 0.000, 60.00, 0.000, 0.000, 109.0, 109.0, 0.000, 0.000, 0.000, 90.00,
566 0.000, 0.000, 0.000, 109.0, 0.000, 0.000, 120.0, 120.0, 120.0, 120.0,
567 120.0, 180.0, 120.0, 120.0, 120.0, 120.0, 0.000, 0.000, 109.0, 0.000,
568 109.0, 0.000, 180.0, 120.0, 120.0, 120.0, 120.0, 120.0, 120.0, 180.0,
569 180.0, 120.0, 120.0, 120.0, 120.0, 120.0, 120.0, 109.0, 120.0, 109.0,
570 0.000, 109.0, 109.0, 120.0, 120.0, 109.0, 0.000, 109.0, 109.0, 109.0 };
571
572 if (ia < 80)
573 return (angs[ia-1]);
574 else
575 return 109.0;
576 }
577 /* --------------------------------- */
578 long kang(int i,int j,int k)
579 {
580 long int kang_v;
581 static long i10000 = 10000;
582 static long i100 = 100;
583
584 kang_v = i10000*i + i100*j + k;
585 return( kang_v );
586 }
587
588 /* -------------------- */
589 long ipack(int i,int j,int k,int l)
590 {
591 long int ipack_v;
592 static long itbig = 200000000;
593 static long itmid = 2000000;
594 static long ithou = 1000;
595
596 ipack_v = itbig*i + itmid*j + ithou*k + l;
597 return( ipack_v );
598 }