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

Line File contents
1 #define EXTERN extern
2 #include "pcwin.h"
3 #include "pcmod.h"
4 #include "field.h"
5
6 EXTERN struct t_vdw1 {
7 int nvdw;
8 float rad[MAXVDWCONST], eps[MAXVDWCONST];
9 int lpd[MAXVDWCONST], ihtyp[MAXVDWCONST], ihdon[MAXVDWCONST];
10 float alpha[MAXVDWCONST],n[MAXVDWCONST],a[MAXVDWCONST],g[MAXVDWCONST];
11 char da[MAXVDWCONST][2];
12 } vdw1;
13
14 EXTERN struct t_cbonds {
15 int nbnd, icbonds[MAXCBND][3], lpd[MAXCBND], ia1[MAXCBND],ia2[MAXCBND];
16 double vrad[MAXCBND],veps[MAXCBND];
17 } cbonds;
18
19 EXTERN int Missing_constants;
20 EXTERN FILE *errfile;
21
22 void kcoord_bonds()
23 {
24 int i, ia, ib,ierr, ita,itb;
25 int jj, nusi;
26 long int mask5,mask6,mask7,mask8;
27 float rad,eps;
28
29 mask5 = 1L << SATMET_MASK;
30 mask6 = 1L << GT18e_MASK;
31 mask7 = 1L << LOWSPIN_MASK;
32 mask8 = 1L << SQPLAN_MASK;
33
34 for (i=0; i < cbonds.nbnd; i++)
35 {
36 ia = cbonds.icbonds[i][0];
37 ib = cbonds.icbonds[i][1];
38 ita = atom[ia].type;
39 itb = atom[ib].type;
40 nusi = 0;
41 rad = vdw1.rad[ita];
42 eps = vdw1.eps[ita];
43 ierr = FALSE;
44 if (strcmp(field.radiusrule,"ARITHMETIC") == 0)
45 {
46 cbonds.vrad[i] = rad + vdw1.rad[atom[ib].type];
47 } else if (strcmp(field.radiusrule,"GEOMETRIC") == 0)
48 {
49 cbonds.vrad[i] = 2.0*( sqrt(rad*vdw1.rad[atom[ib].type]));
50 } else
51 cbonds.vrad[i] = rad + vdw1.rad[atom[ib].type];
52
53 if (strcmp(field.epsrule,"ARITHMETIC") == 0)
54 {
55 cbonds.veps[i] = 0.5*(eps + vdw1.eps[atom[ib].type]);
56 } else if (strcmp(field.epsrule,"GEOMETRIC") == 0)
57 {
58 cbonds.veps[i] = sqrt( eps*vdw1.eps[atom[ib].type] );
59 }else
60 cbonds.veps[i] = sqrt( eps*vdw1.eps[atom[ib].type] );
61
62 if ( (atom[ia].flags & mask6) && !( atom[ia].flags & mask5 ))
63 nusi = 2;
64 else if( ( atom[ia].flags &mask5 ) && ( atom[ia].flags & mask6 ) )
65 {
66 if( ( atom[ia].flags & mask7 ) && !( atom[ia].flags & mask8 ) )
67 nusi = 1;
68 else if(( atom[ia].flags & mask8 ) && !( atom[ia].flags & mask7 ) )
69 nusi = 3;
70 }
71
72 // metal to pi atom
73 if (atom[ib].type == 2 || atom[ib].type == 4 || atom[ib].type == 29 ||
74 atom[ib].type == 30 || atom[ib].type == 40 || atom[ib].type == 48 )
75 {
76 cbonds.ia1[i] = 0;
77 cbonds.ia2[i] = 0;
78 for (jj = 0 ; jj < MAXIAT; jj++)
79 {
80 if (atom[ib].iat[jj] != 0 && atom[ib].iat[jj] != ia && cbonds.ia1[i] == 0)
81 cbonds.ia1[i] = atom[ib].iat[jj];
82 else if (atom[ib].iat[jj] != 0 && atom[ib].iat[jj] != ia && cbonds.ia2[i] == 0)
83 cbonds.ia2[i] = atom[ib].iat[jj];
84 }
85 cbonds.vrad[i] = rad + 0.66;
86 if (nusi > 0)
87 cbonds.vrad[i] += 0.15;
88 cbonds.lpd[i] = vdw1.lpd[atom[ib].type];
89 if (nusi > 0)
90 cbonds.lpd[i] *= 0.85;
91 } else
92 { // metal to nitrogen and oxygen lone pairs
93 cbonds.vrad[i] = rad + vdw1.rad[atom[ib].type]; // use radius of lone pair
94 cbonds.vrad[i] *= 0.6; // cut radius in half
95
96 cbonds.lpd[i] = vdw1.lpd[atom[ib].type];
97 if (atom[ib].type == 7)
98 cbonds.lpd[i] *= 3.0;
99 else if (atom[ib].type == 10)
100 cbonds.lpd[i] *= 4.0;
101
102 if( ( atom[ia].flags & mask6 ) && !( atom[ia].flags & mask5 ) )
103 {
104 cbonds.vrad[i] += .15;
105 cbonds.lpd[i] *= .85;
106 }
107
108 if( nusi == 1 || (nusi == 3 && (atom[atom[ib].iat[0]].type == 6 || atom[atom[ib].iat[0]].type == 15)) )
109 cbonds.vrad[i] -= .1;
110 }
111 }
112
113 }