ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/freemol/trunk/src/mengine/src/kvdw.c
(Generate patch)
# Line 2 | Line 2
2   #include "pcwin.h"
3   #include "pcmod.h"
4   #include "nonbond.h"
5 #include "pot.h"
6 #include "field.h"
5   #include "atom_k.h"
6  
7   void numeral(int,char *,int);
8 + int get_field(void);
9 + char * get_radiustype(void);
10 + char * get_radiussize(void);
11 + char * get_radiusrule(void);
12 + char * get_epsrule(void);
13        
14   EXTERN struct t_vdw1 {
15          int  nvdw;
# Line 33 | Line 36
36          
37   void kvdw()
38   {
39 <   int i, j, k;
39 >  int i, j, k,field;
40     int it, ita, itb, found, ianum,ibnum;
41     int dapair;
42     char pa[4],pb[4],pt[7];
40   long int hbmask;
43     float rad_ia, rad_ib;
44     float eps1, eps2, seps1, seps2;
45     float rii,rjj,gammaij, rij, epsij;
46    
47 +   field = get_field();
48   // mm3 setup
49 <   if (field.type == MM3)
49 >   if (field == MM3)
50     {
51          mm3hbond.npair = 0;
52          for (i=0; i < vdwpr_k.nvdwpr; i++)
# Line 66 | Line 69
69     }
70   //                
71     nonbond.npair = 0;
69   hbmask = 1 << HBOND_MASK;
72     for (i=0; i <= nonbond.maxnbtype; i++)
73     {
74         for (j=0; j <= nonbond.maxnbtype; j++)
# Line 75 | Line 77
77    
78     for (i=1; i < natom; i++)
79     {
80 <       ita = atom[i].type;
81 <       ianum = atom[i].atomnum;
80 >       ita = atom.type[i];
81 >       ianum = atom.atomnum[i];
82         for (j= i+1; j <= natom; j++)
83         {
84 <           itb = atom[j].type;
85 <           ibnum = atom[j].atomnum;
84 >           itb = atom.type[j];
85 >           ibnum = atom.atomnum[j];
86             if (ita < itb)
87                it = 1000*ita + itb;
88             else
# Line 92 | Line 94
94  
95             if (found == FALSE)
96             {
97 <               if (field.type == MMFF94)
97 >               if (field == MMFF94)
98                 {
99                          dapair = FALSE;
100                          rii = vdw1.a[ita]*pow(vdw1.alpha[ita],0.25);
101                          rjj = vdw1.a[itb]*pow(vdw1.alpha[itb],0.25);
102 <                        if ((atom[i].atomnum == 1 && ita != 5) || (atom[j].atomnum ==1  && itb != 5) )  // polar hydrogens
102 >                        if ((atom.atomnum[i] == 1 && ita != 5) || (atom.atomnum[j] ==1  && itb != 5) )  // polar hydrogens
103                          {
104                                  rij = 0.5*(rii+rjj);
105                                  if ( (strcmp(vdw1.da[ita],"A") == 0) && (strcmp(vdw1.da[itb],"D") == 0) )
# Line 138 | Line 140
140                 {
141                    rad_ia = vdw1.rad[ita];
142                    rad_ib = vdw1.rad[itb];
143 <                  if (strcmp(field.radiustype,"SIGMA") == 0)
143 >                  if (strcmp(get_radiustype(),"SIGMA") == 0)
144                    {
145                        rad_ia *= 1.122462048;
146                        rad_ib *= 1.122462048;
147                    }  
148 <                  if (strcmp(field.radiussize,"DIAMETER") == 0)
148 >                  if (strcmp(get_radiussize(),"DIAMETER") == 0)
149                    {
150                        rad_ia /= 2.0;
151                        rad_ib /= 2.0;
# Line 151 | Line 153
153                    eps1 = vdw1.eps[ita];
154                    eps2 = vdw1.eps[itb];
155        
156 <                  if (strcmp(field.radiusrule,"ARITHMETIC") == 0)
156 >                  if (strcmp(get_radiusrule(),"ARITHMETIC") == 0)
157                    {
158                       nonbond.vrad[ita][itb] = rad_ia + rad_ib;
159                       nonbond.vrad[itb][ita] = rad_ia + rad_ib;
160 <                  } else if (strcmp(field.radiusrule,"GEOMETRIC") == 0)
160 >                  } else if (strcmp(get_radiusrule(),"GEOMETRIC") == 0)
161                    {
162                       nonbond.vrad[ita][itb] = 2.0*( sqrt(rad_ia)*sqrt(rad_ib));
163                       nonbond.vrad[itb][ita] = 2.0*( sqrt(rad_ia)*sqrt(rad_ib));
164 <                  } else if (strcmp(field.radiusrule,"CUBIC-MEAN") == 0)
164 >                  } else if (strcmp(get_radiusrule(),"CUBIC-MEAN") == 0)
165                    {
166 <                     if ((atom[i].atomnum == 1 && ita != 5) || (atom[j].atomnum ==1  && itb != 5) )  // polar hydrogens
166 >                     if ((atom.atomnum[i] == 1 && ita != 5) || (atom.atomnum[j] ==1  && itb != 5) )  // polar hydrogens
167                       {
168                          nonbond.vrad[ita][itb] = rad_ia + rad_ib;
169                          nonbond.vrad[itb][ita] = rad_ia + rad_ib;
# Line 176 | Line 178
178                       nonbond.vrad[itb][ita] = rad_ia + rad_ib;
179                    }
180        
181 <                  if (strcmp(field.epsrule,"ARITHMETIC") == 0)
181 >                  if (strcmp(get_epsrule(),"ARITHMETIC") == 0)
182                    {
183                        nonbond.veps[ita][itb] = 0.5*(eps1 + eps2);
184                        nonbond.veps[itb][ita] = 0.5*(eps1 + eps2);
185 <                  } else if (strcmp(field.epsrule,"GEOMETRIC") == 0)
185 >                  } else if (strcmp(get_epsrule(),"GEOMETRIC") == 0)
186                    {
187                        nonbond.veps[ita][itb] = sqrt( eps1*eps2 );
188                        nonbond.veps[itb][ita] = sqrt( eps1*eps2 );
189 <                  } else if (strcmp(field.epsrule,"HARMONIC") == 0)
189 >                  } else if (strcmp(get_epsrule(),"HARMONIC") == 0)
190                   {
191                        nonbond.veps[ita][itb] = 2.0* (eps1*eps2)/(eps1+eps2);
192                        nonbond.veps[itb][ita] = 2.0* (eps1*eps2)/(eps1+eps2);
193 <                  } else if (strcmp(field.epsrule,"HHG") == 0)
193 >                  } else if (strcmp(get_epsrule(),"HHG") == 0)
194                    {
195 <                     if ((atom[i].atomnum == 1 && ita != 5) || (atom[j].atomnum ==1  && itb != 5) )  // polar hydrogens
195 >                     if ((atom.atomnum[i] == 1 && ita != 5) || (atom.atomnum[j] ==1  && itb != 5) )  // polar hydrogens
196                       {
197                         nonbond.veps[ita][itb] = 0.5*(eps1 + eps2);
198                         nonbond.veps[itb][ita] = 0.5*(eps1 + eps2);
# Line 207 | Line 209
209                        nonbond.veps[itb][ita] = sqrt( eps1*eps2 );
210                    }
211  
212 <                  if (field.type == MM3)
212 >                  if (field == MM3)
213                    {
214                        nonbond.ipif[ita][itb] = 1;
215                        nonbond.ipif[itb][ita] = 1;
216 <                  } else if (field.type == MMX)
216 >                  } else if (field == MMX)
217                    {
218                        nonbond.ipif[ita][itb] = 1;
219                        nonbond.ipif[itb][ita] = 1;
# Line 228 | Line 230
230                     nonbond.veps14[itb][ita] = nonbond.veps[itb][ita];
231                 }
232    // look for any special vdw pairs
233 <               numeral(atom[i].type,pa,3);
234 <               numeral(atom[j].type,pb,3);
233 >               numeral(atom.type[i],pa,3);
234 >               numeral(atom.type[j],pb,3);
235        
236 <              if (atom[i].type < atom[j].type)
236 >              if (atom.type[i] < atom.type[j])
237                {
238                  strcpy(pt,pa);
239                  strcat(pt,pb);
# Line 269 | Line 271
271               nonbond.iNBtype[ita][itb] = it;
272               nonbond.iNBtype[itb][ita] = it;
273              
272             if (pot.use_hbond)
273             {
274                 if ( atom[i].flags & hbmask || atom[i].mmx_type == 20)
275                 {
276                     for (k=0; k < MAXIAT; k++)
277                     {
278                         if (atom[j].iat[k] != 0 && atom[j].bo[k] != 9)
279                         {
280                             if (atom[atom[j].iat[k]].flags & hbmask || atom[atom[j].iat[k]].mmx_type == 20)
281                             {
282                                nonbond.vrad[ita][itb] -= 0.25;
283                                nonbond.vrad[itb][ita] -= 0.25;
284                                break;
285                             }
286                         }
287                     }
288                 } else if (atom[j].flags & hbmask || atom[j].mmx_type == 20)
289                 {
290                     for (k=0; k < MAXIAT; k++)
291                     {
292                         if (atom[i].iat[k] != 0 && atom[i].bo[k] != 9)
293                         {
294                             if (atom[atom[i].iat[k]].flags & hbmask || atom[atom[i].iat[k]].mmx_type == 20)
295                             {
296                                nonbond.vrad[ita][itb] -= 0.25;
297                                nonbond.vrad[itb][ita] -= 0.25;
298                                break;
299                             }
300                         }
301                     }
302                 }
303             }
274               nonbond.npair++;
275             }
276         }
277     }
278     for (i=1; i <= natom ; i++)
279     {
280 <       ita = atom[i].type;
280 >       ita = atom.type[i];
281         if (vdw1.rad[ita] == 0)
282         {
283              rii = vdw1.a[ita]*pow(vdw1.alpha[ita],0.25);

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines