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

Line User Rev File contents
1 tjod 3 #define EXTERN extern
2    
3     #include "pcwin.h"
4     #include "pcmod.h"
5     #include "bonds_ff.h"
6     #include "field.h"
7     #include "atom_k.h"
8     #include "utility.h"
9     #include "pdb.h"
10    
11     #define MAXOPLS 600
12    
13     void numeral(int, char *, int);
14     void compute_fcharge(void);
15     void assign_amber_chrg(int);
16     void assign_gast_charge();
17     void compute_opls(void);
18     int find_rsize(int,int);
19     void search_rings(int);
20     void message_alert(char *,char *);
21     void find_residues(void);
22     void bbchk(void);
23     void pireset(void);
24     int get_oplsconst(int,int,char *);
25     void ksortr(int, int *,int *,int *);
26     void get_rsize(int,int,int, int *);
27     void ksort_ring(int,int,int *);
28     void get_firstlevel(int,int *,int *,int *,int *);
29     void get_secondlevel(int,int,int *,int *,int *,int *);
30     void get_residue_type(int,int*);
31    
32     EXTERN struct ElementType {
33     char symbol[3];
34     int atomnum;
35     float weight, covradius, vdwradius;
36     int s,p,d,f, type;
37     } Elements[] ;
38    
39     EXTERN struct BioType {
40     char name_3[4], name_1[2];
41     int type, natom;
42     } BioGroups[];
43    
44     typedef struct Gastype {
45     int type;
46     float a,b,c,d;
47     } Gastype;
48    
49     EXTERN Gastype Gastchrg[];
50    
51     EXTERN struct t_bondpk {
52     int npibond,npibond44i,npibond45i,npibond46i,npibond55i,npibond56i,npibond66i;
53     int npibond44o,npibond45o,npibond46o,npibond55o,npibond56o,npibond66o;
54     int npibond40,npibond50,npibond60;
55     int use_pibond,use_pibond44i,use_pibond45i,use_pibond46i,use_pibond55i,use_pibond56i,use_pibond66i;
56     int use_pibond44o,use_pibond45o,use_pibond46o,use_pibond55o,use_pibond56o,use_pibond66o;
57     int use_pibond40,use_pibond50,use_pibond60;
58     char kb[MAXPIBONDCONST][7],kb44i[20][7],kb45i[20][7],kb46i[20][7],kb55i[20][7],kb56i[20][7],kb66i[20][7];
59     char kb44o[20][7],kb45o[20][7],kb46o[20][7],kb55o[20][7],kb56o[20][7],kb66o[20][7];
60     char kb40[20][7],kb50[20][7],kb60[20][7];
61     float bk[MAXPIBONDCONST], bl[MAXPIBONDCONST], bmom[MAXPIBONDCONST],sslop[MAXPIBONDCONST], tslop[MAXPIBONDCONST], tslop2[MAXPIBONDCONST];
62     float bk44i[20], bl44i[20], bmom44i[20],sslop44i[20], tslop44i[20], tslop244i[20];
63     float bk45i[20], bl45i[20], bmom45i[20],sslop45i[20], tslop45i[20], tslop245i[20];
64     float bk46i[20], bl46i[20], bmom46i[20],sslop46i[20], tslop46i[20], tslop246i[20];
65     float bk55i[20], bl55i[20], bmom55i[20],sslop55i[20], tslop55i[20], tslop255i[20];
66     float bk56i[20], bl56i[20], bmom56i[20],sslop56i[20], tslop56i[20], tslop256i[20];
67     float bk66i[20], bl66i[20], bmom66i[20],sslop66i[20], tslop66i[20], tslop266i[20];
68     float bk44o[20], bl44o[20], bmom44o[20],sslop44o[20], tslop44o[20], tslop244o[20];
69     float bk45o[20], bl45o[20], bmom45o[20],sslop45o[20], tslop45o[20], tslop245o[20];
70     float bk46o[20], bl46o[20], bmom46o[20],sslop46o[20], tslop46o[20], tslop246o[20];
71     float bk55o[20], bl55o[20], bmom55o[20],sslop55o[20], tslop55o[20], tslop255o[20];
72     float bk56o[20], bl56o[20], bmom56o[20],sslop56o[20], tslop56o[20], tslop256o[20];
73     float bk66o[20], bl66o[20], bmom66o[20],sslop66o[20], tslop66o[20], tslop266o[20];
74     float bk40[20], bl40[20], bmom40[20],sslop40[20], tslop40[20], tslop240[20];
75     float bk50[20], bl50[20], bmom50[20],sslop50[20], tslop50[20], tslop250[20];
76     float bk60[20], bl60[20], bmom60[20],sslop60[20], tslop60[20], tslop260[20];
77     } bondpk;
78    
79     EXTERN struct t_charge_k {
80     int ncharge, nbndchrg, nbndchrgdel;
81     int type[MAXATOMTYPE], btype[MAXBONDCONST], btypedel[MAXBONDCONST];
82     float charge[MAXATOMTYPE], bcharge[MAXBONDCONST], formchrg[MAXATOMTYPE], bchargedel[MAXBONDCONST];
83     float typechrg[MAXATOMTYPE];
84     } charge_k;
85    
86     EXTERN struct t_amberchrg_k {
87     int ncharge, type[MAXATOMTYPE], res_type[MAXATOMTYPE];
88     char symbol[MAXATOMTYPE][4];
89     float chrg[MAXATOMTYPE];
90     } amberchrg_k;
91    
92     EXTERN struct t_oplschrg_k {
93     int ncharge, type[MAXOPLS];
94     char chrgstring[MAXOPLS][25];
95     float chrg[MAXOPLS];
96     } oplschrg_k;
97    
98     EXTERN struct t_dipole_k {
99     int ndipole,ndipole3,ndipole4,ndipole5;
100     char kb[MAXBONDCONST][7], kb3[MAXBOND3CONST][7], kb4[MAXBOND4CONST][7], kb5[MAXBOND5CONST][7];
101     float bmom[MAXBONDCONST],bmom3[MAXBOND3CONST],bmom4[MAXBOND4CONST],bmom5[MAXBOND5CONST];
102     } dipole_k;
103    
104     EXTERN struct t_residues {
105     int nchain, ichainstart[10];
106     int ngroup, iresnum[200], irestype[200], istartatom[200];
107     } residues;
108     EXTERN struct t_ringdata {
109     int nring3, nring4, nring5, nring6;
110     int tot3, tot4, tot5, tot6;
111     int **iring3,**iring4,**iring5,**iring6;
112     } ringdata;
113    
114    
115     void kcharge()
116     {
117     int i, j, ia, ib, it, iit, kit, jji;
118     long int mask;
119     char pa[4],pb[4],pt[7];
120     float bmomj;
121    
122     mask = 1L << 0;
123    
124     for (i=1; i <= natom; i++)
125     {
126     atom[i].sigma_charge = 0.0;
127     if (field.type == MMX && (atom[i].type < 300))
128     atom[i].charge = 0.0;
129     else if (field.type != MMX && atom[i].type < 300)
130     atom[i].charge = 0.0;
131     }
132    
133     if (field.type == MMFF94)
134     {
135     compute_fcharge();
136     for(i=0; i < bonds_ff.nbnd; i++)
137     {
138     ia = bonds_ff.i12[i][0];
139     ib = bonds_ff.i12[i][1];
140     iit = atom[ia].type;
141     kit = atom[ib].type;
142     if (iit < kit)
143     it = iit*100+kit;
144     else
145     it = kit*100+iit;
146     if (bonds_ff.index[i] == 1)
147     {
148     for(j=0; j < charge_k.nbndchrgdel; j++)
149     {
150     if (it == charge_k.btypedel[j])
151     {
152     if (iit < kit)
153     {
154     atom[ia].charge -= charge_k.bchargedel[j];
155     atom[ib].charge += charge_k.bchargedel[j];
156     }else
157     {
158     atom[ib].charge -= charge_k.bchargedel[j];
159     atom[ia].charge += charge_k.bchargedel[j];
160     }
161     break;
162     }
163     }
164     } else
165     {
166     for(j=0; j < charge_k.nbndchrg; j++)
167     {
168     if (it == charge_k.btype[j])
169     {
170     if (iit < kit)
171     {
172     atom[ia].charge -= charge_k.bcharge[j];
173     atom[ib].charge += charge_k.bcharge[j];
174     }else
175     {
176     atom[ib].charge -= charge_k.bcharge[j];
177     atom[ia].charge += charge_k.bcharge[j];
178     }
179     break;
180     }
181     }
182     }
183     }
184     // add in formal charges
185     for (i=1; i <= natom; i++)
186     {
187     if (atom[i].atomnum == 15)
188     atom[i].formal_charge = 0.0;
189     if (atom[i].atomnum == 16 && atom[i].type != 72)
190     atom[i].formal_charge = 0.0;
191     }
192     for (i=1; i <= natom; i++)
193     {
194     atom[i].charge += (1.0 - atom_k.ligands[atom[i].type]*charge_k.formchrg[atom[i].type])*atom[i].formal_charge ;
195     for(j=0; j < MAXIAT; j++)
196     {
197     if (atom[i].iat[j] != 0)
198     {
199     if (atom[atom[i].iat[j]].formal_charge < 0.0)
200     atom[i].charge += charge_k.formchrg[atom[atom[i].iat[j]].type]*atom[atom[i].iat[j]].formal_charge;
201     }
202     }
203     }
204     }
205     }
206     /* -------------------------------------------------------- */
207     void compute_fcharge()
208     {
209     int i, j, jjbo;
210     int jatm, k, jji, attached[4];
211     int *used;
212     float charge;
213    
214     used = ivector(1,natom+1);
215     for (i=1; i <= natom; i++)
216     used[i] = FALSE;
217    
218     for (i=1; i <= natom; i++)
219     {
220     atom[i].formal_charge = charge_k.typechrg[atom[i].type];
221     }
222     // find variable charges: types 32, 72, 76, 81
223     for (i=1; i <= natom; i++)
224     {
225     if (used[i] == FALSE)
226     {
227     if (atom[i].formal_charge != 0.0)
228     {
229     if (atom[i].type == 32)
230     {
231     jatm = atom[i].iat[0];
232     jji = 0;
233     jjbo = 0;
234     for (k=0; k < MAXIAT; k++)
235     {
236     if (atom[jatm].iat[k] != 0 && atom[jatm].bo[k] != 9)
237     {
238     if (atom[atom[jatm].iat[k]].type == 32)
239     {
240     if (atom[jatm].bo[k] == 1)
241     {
242     jjbo++;
243     attached[jji] = atom[jatm].iat[k];
244     jji++;
245     } else
246     {
247     attached[jji] = atom[jatm].iat[k];
248     jji++;
249     }
250     }
251     }
252     }
253     charge = -(double)jjbo/jji;
254     if (jji == 1) // only 1 type 32
255     atom[i].formal_charge = 0.0;
256     else if (jji == 2)
257     {
258     if (atom[jatm].atomnum == 15) // p-o treat as single negative charge;
259     charge = -1.00/jji;
260     if (atom[jatm].type == 73)
261     {
262     charge = -0.500;
263     atom[attached[0]].formal_charge = charge;
264     atom[attached[1]].formal_charge = charge;
265     used[attached[0]] = TRUE;
266     used[attached[1]] = TRUE;
267     } else if (atom[jatm].atomnum != 7)
268     {
269     atom[attached[0]].formal_charge = charge;
270     atom[attached[1]].formal_charge = charge;
271     used[attached[0]] = TRUE;
272     used[attached[1]] = TRUE;
273     }else
274     {
275     atom[attached[0]].formal_charge = 0.0;
276     atom[attached[1]].formal_charge = 0.0;
277     used[attached[0]] = TRUE;
278     used[attached[1]] = TRUE;
279     }
280     } else if (jji == 3)
281     {
282     atom[attached[0]].formal_charge = charge;
283     atom[attached[1]].formal_charge = charge;
284     atom[attached[2]].formal_charge = charge;
285     used[attached[0]] = TRUE;
286     used[attached[1]] = TRUE;
287     used[attached[2]] = TRUE;
288     } else if (jji == 4)
289     {
290     atom[attached[0]].formal_charge = charge;
291     atom[attached[1]].formal_charge = charge;
292     atom[attached[2]].formal_charge = charge;
293     atom[attached[3]].formal_charge = charge;
294     used[attached[0]] = TRUE;
295     used[attached[1]] = TRUE;
296     used[attached[2]] = TRUE;
297     used[attached[3]] = TRUE;
298     }
299     }else if (atom[i].type == 72)
300     {
301     jatm = atom[i].iat[0];
302     jji = 0;
303     jjbo = 0;
304     for (k=0; k < MAXIAT; k++)
305     {
306     if (atom[jatm].iat[k] != 0 && atom[jatm].bo[k] != 9)
307     {
308     if (atom[atom[jatm].iat[k]].type == 72)
309     {
310     if (atom[jatm].bo[k] == 1)
311     {
312     jjbo++;
313     attached[jji] = atom[jatm].iat[k];
314     jji++;
315     } else
316     {
317     attached[jji] = atom[jatm].iat[k];
318     jji++;
319     }
320     }
321     }
322     }
323     if (jji == 1)
324     {
325     if (atom[jatm].atomnum == 15)
326     atom[i].formal_charge = 0.0;
327     else
328     atom[i].formal_charge = -1.00;
329     }else
330     {
331     atom[attached[0]].formal_charge = -0.5;
332     atom[attached[1]].formal_charge = -0.5;
333     used[attached[0]] = TRUE;
334     used[attached[1]] = TRUE;
335     }
336     } else if (atom[i].type == 81)
337     {
338     jji = 0;
339     for (k=0; k < MAXIAT ; k++)
340     {
341     if (atom[i].iat[k] != 0 && atom[atom[i].iat[k]].type == 80)
342     {
343     jatm = atom[i].iat[k];
344     jji = 0;
345     for (j=0; j < MAXIAT; j++)
346     {
347     if (atom[jatm].iat[j] != 0 && atom[jatm].bo[j] != 9)
348     {
349     if (atom[atom[jatm].iat[j]].atomnum == 7)
350     {
351     attached[jji] = atom[jatm].iat[j];
352     jji++;
353     }
354     }
355     }
356     }
357     }
358     if (jji == 0)
359     {
360     charge = 1.00;
361     atom[i].formal_charge = charge;
362     used[i] = TRUE;
363     } else
364     {
365     charge = 1.00/jji;
366     for (j=0; j < jji; j++)
367     {
368     atom[attached[j]].formal_charge = charge;
369     used[attached[j]] = TRUE;
370     }
371     }
372     }
373     }
374     }
375     }
376     free_ivector(used,1, natom+1);
377     }