ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/freemol/trunk/src/mengine/src/get_mem.c
Revision: 110
Committed: Thu Mar 12 01:43:50 2009 UTC (12 years, 4 months ago) by gilbertke
File size: 9760 byte(s)
Log Message:
further cleanup and localization of atom data
Line User Rev File contents
1 wdelano 58 #define EXTERN extern
2    
3     #include "pcwin.h"
4     #include "pcmod.h"
5 gilbertke 63 #include "attached.h"
6 wdelano 58 #include "derivs.h"
7     #include "hess.h"
8     #include "utility.h"
9     #include "nonbond.h"
10     #include "torsions.h"
11    
12 gilbertke 103 void get_molecule_memory(int);
13     int use_solv(void);
14 gilbertke 105 void allocate_rings(int);
15 gilbertke 110 void max_torsions(int natom,int *type,int **iat,int **bo);
16 wdelano 58
17     #define STILL 1
18     #define HCT 2
19    
20     EXTERN struct t_solvent {
21     int type;
22 gilbertke 103 double EPSin, EPSsolv;
23 wdelano 58 double doffset, p1,p2,p3,p4,p5;
24     double *shct,*asolv,*rsolv,*vsolv,*gpol,*rborn;
25     } solvent;
26    
27 gilbertke 103 // =======================================
28 gilbertke 104 void get_molecule_memory(int niatom)
29 gilbertke 103 {
30     int i,j;
31 wdelano 58
32 gilbertke 105 MAXATOM = niatom+10;
33     MAXBND = 6*MAXATOM/5;
34     MAXANG = 12*MAXATOM/5;
35     MAXTOR = 4*MAXATOM;
36    
37 gilbertke 103 atom.type = ivector(0,niatom+10);
38     atom.tclass = ivector(0,niatom+10);
39     atom.mmx_type = ivector(0,niatom+10);
40     atom.mm3_type = ivector(0,niatom+10);
41     atom.mmff_type = ivector(0,niatom+10);
42     atom.atomnum = ivector(0,niatom+10);
43     atom.use = ivector(0,niatom+10);
44     atom.flags = ilvector(0,niatom+10);
45     atom.iat = imatrix(0, niatom+10,0,MAXIAT);
46     atom.bo = imatrix(0,niatom+10, 0, MAXIAT);
47     atom.x = dvector(0,niatom+10);
48     atom.y = dvector(0,niatom+10);
49     atom.z = dvector(0,niatom+10);
50     atom.charge = dvector(0,niatom+10);
51     atom.formal_charge = dvector(0,niatom+10);
52     atom.sigma_charge = dvector(0,niatom+10);
53     atom.atomwt = dvector(0,niatom+10);
54     atom.radius = dvector(0,niatom+10);
55 gilbertke 104 atom.name = malloc( (niatom+10)*sizeof(LABEL));
56 gilbertke 103
57     for (i=0; i <= niatom; i++)
58     {
59     atom.type[i] = 0;
60     atom.tclass[i] = 0;
61     atom.mmx_type[i] = 0;
62     atom.mm3_type[i] = 0;
63     atom.mmff_type[i] = 0;
64     atom.atomnum[i] = 0;
65     atom.use[i] = 0;
66     atom.flags[i] = 0;
67     atom.x[i] = 0.0;
68     atom.y[i] = 0.0;
69     atom.z[i] = 0.0;
70     atom.charge[i] = 0.0;
71     atom.formal_charge[i] = 0.0;
72     atom.sigma_charge[i] = 0.0;
73     atom.atomwt[i] = 0.0;
74     for (j=0; j < MAXIAT; j++)
75     {
76     atom.iat[i][j] = 0;
77     atom.bo[i][j] = 0;
78     }
79 gilbertke 104 strcpy(atom.name[i],"");
80 gilbertke 103 }
81 gilbertke 105 allocate_rings(niatom); // allocate space for ring data
82 gilbertke 104 }
83 wdelano 58 /* ================================================================ */
84     void get_memory()
85     {
86     int i, j;
87     int ntor;
88     int ntypes, found, itype[MAXATOMTYPE];
89    
90 gilbertke 89 skip = imatrix(0,natom+1,0,natom+1);
91     for (i=1; i <= natom; i++)
92     {
93     for (j=1; j <=natom; j++)
94     skip[i][j] = 0;
95     }
96 wdelano 58
97     deriv.d1 = dmatrix(0,natom+1, 0,3);
98     deriv.deb = dmatrix(0,natom+1, 0,3);
99     deriv.dea = dmatrix(0,natom+1, 0,3);
100     deriv.destb = dmatrix(0,natom+1, 0,3);
101     deriv.deopb = dmatrix(0,natom+1, 0,3);
102     deriv.detor = dmatrix(0,natom+1, 0,3);
103     deriv.de14 = dmatrix(0,natom+1, 0,3);
104     deriv.devdw = dmatrix(0,natom+1, 0,3);
105     deriv.deqq = dmatrix(0,natom+1, 0,3);
106     deriv.deaa = dmatrix(0,natom+1, 0,3);
107     deriv.destor = dmatrix(0,natom+1, 0,3);
108     deriv.dehb = dmatrix(0,natom+1, 0,3);
109     deriv.deimprop = dmatrix(0,natom+1, 0,3);
110     deriv.deub = dmatrix(0,natom+1, 0,3);
111     deriv.desolv = dmatrix(0,natom+1, 0,3);
112     deriv.degeom = dmatrix(0,natom+1, 0,3);
113     deriv.drb = dvector(0,natom+1);
114    
115     ntypes = 0;
116     nonbond.maxnbtype = 0;
117    
118     for (i=1; i <= natom; i++)
119     {
120     found = FALSE;
121 gilbertke 103 if (atom.type[i] > nonbond.maxnbtype)
122     nonbond.maxnbtype = atom.type[i];
123 wdelano 58 for (j=0; j < ntypes; j++)
124     {
125 gilbertke 103 if (atom.type[i] == itype[j])
126 wdelano 58 {
127     found = TRUE;
128     break;
129     }
130     }
131     if (found == FALSE)
132     {
133 gilbertke 103 itype[ntypes] = atom.type[i];
134 wdelano 58 ntypes++;
135     }
136     }
137    
138     j=0;
139     for (i=0; i <= ntypes; i++)
140     j+= i;
141    
142     nonbond.nonbond = j;
143     nonbond.iNBtype = imatrix(0,nonbond.maxnbtype+1, 0,nonbond.maxnbtype+1);
144     nonbond.ipif = imatrix(0,nonbond.maxnbtype+1, 0,nonbond.maxnbtype+1);
145 gilbertke 103 nonbond.vrad = dmatrix(0,nonbond.maxnbtype+1, 0,nonbond.maxnbtype+1);
146     nonbond.veps = dmatrix(0,nonbond.maxnbtype+1, 0,nonbond.maxnbtype+1);
147     nonbond.vrad14 = dmatrix(0,nonbond.maxnbtype+1, 0,nonbond.maxnbtype+1);
148     nonbond.veps14 = dmatrix(0,nonbond.maxnbtype+1, 0,nonbond.maxnbtype+1);
149 wdelano 58
150 gilbertke 110 max_torsions(natom,atom.type,atom.iat,atom.bo);
151 wdelano 58 ntor = torsions.ntor;
152     if (ntor > 0)
153     {
154     torsions.i14 = imatrix(0,ntor,0,4);
155     torsions.v1 = vector(0,ntor);
156     torsions.v2 = vector(0,ntor);
157     torsions.v3 = vector(0,ntor);
158     torsions.v4 = vector(0,ntor);
159     torsions.v5 = vector(0,ntor);
160     torsions.v6 = vector(0,ntor);
161     torsions.ph1 = ivector(0,ntor);
162     torsions.ph2 = ivector(0,ntor);
163     torsions.ph3 = ivector(0,ntor);
164     torsions.ph4 = ivector(0,ntor);
165     torsions.ph5 = ivector(0,ntor);
166     torsions.ph6 = ivector(0,ntor);
167     for (i=0; i < ntor; i++)
168     {
169     torsions.v1[i] = 0.0;
170     torsions.v2[i] = 0.0;
171     torsions.v3[i] = 0.0;
172     torsions.v4[i] = 0.0;
173     torsions.v5[i] = 0.0;
174     torsions.v6[i] = 0.0;
175     torsions.ph1[i] = 0;
176     torsions.ph2[i] = 0;
177     torsions.ph3[i] = 0;
178     torsions.ph4[i] = 0;
179     torsions.ph5[i] = 0;
180     torsions.ph6[i] = 0;
181     }
182     }
183     hess.hessx = matrix(0,natom+1, 0,3);
184     hess.hessy = matrix(0,natom+1, 0,3);
185     hess.hessz = matrix(0,natom+1, 0,3);
186    
187     attached.n13 = ivector(0,natom+1);
188     attached.n14 = ivector(0,natom+1);
189     attached.i13 = imatrix(0, 20, 0,natom+1);
190     attached.i14 = imatrix(0,144, 0,natom+1);
191    
192     for(i=0; i <= natom; i++)
193     {
194     for (j=0; j < 3; j++)
195     {
196     deriv.d1[i][j] = 0.0;
197     deriv.deb[i][j] = 0.0;
198     deriv.dea[i][j] = 0.0;
199     deriv.destb[i][j] = 0.0;
200     deriv.deopb[i][j] = 0.0;
201     deriv.detor[i][j] = 0.0;
202     deriv.de14[i][j] = 0.0;
203     deriv.devdw[i][j] = 0.0;
204     deriv.deqq[i][j] = 0.0;
205     deriv.deaa[i][j] = 0.0;
206     deriv.destor[i][j] = 0.0;
207     deriv.dehb[i][j] = 0.0;
208     deriv.deimprop[i][j] = 0.0;
209     deriv.deub[i][j] = 0.0;
210    
211     hess.hessx[i][j] = 0.0;
212     hess.hessy[i][j] = 0.0;
213     hess.hessz[i][j] = 0.0;
214     }
215     }
216 gilbertke 103 if (use_solv())
217 wdelano 58 {
218     solvent.asolv = dvector(0,natom+1);
219     solvent.rsolv = dvector(0,natom+1);
220     solvent.rborn = dvector(0,natom+1);
221     if (solvent.type == STILL)
222     {
223     solvent.vsolv = dvector(0,natom+1);
224     solvent.gpol = dvector(0,natom+1);
225     } else if (solvent.type == HCT)
226     {
227     solvent.shct = dvector(0,natom+1);
228     }
229     }
230     }
231    
232     void free_memory()
233     {
234     int ntor;
235 gilbertke 89 free_imatrix(skip, 0, natom+1,0,natom+1);
236 wdelano 58 free_dmatrix(deriv.d1, 0,natom+1, 0,3);
237     free_dmatrix(deriv.deb, 0,natom+1, 0,3);
238     free_dmatrix(deriv.dea, 0,natom+1, 0,3);
239     free_dmatrix(deriv.destb, 0,natom+1, 0,3);
240     free_dmatrix(deriv.deopb, 0,natom+1, 0,3);
241     free_dmatrix(deriv.detor, 0,natom+1, 0,3);
242     free_dmatrix(deriv.de14, 0,natom+1, 0,3);
243     free_dmatrix(deriv.devdw, 0,natom+1, 0,3);
244     free_dmatrix(deriv.deqq, 0,natom+1, 0,3);
245     free_dmatrix(deriv.deaa, 0,natom+1, 0,3);
246     free_dmatrix(deriv.destor, 0,natom+1, 0,3);
247     free_dmatrix(deriv.dehb, 0,natom+1, 0,3);
248     free_dmatrix(deriv.deimprop,0,natom+1, 0,3);
249     free_dmatrix(deriv.deub, 0,natom+1, 0,3);
250     free_dmatrix(deriv.desolv ,0,natom+1, 0,3);
251     free_dmatrix(deriv.degeom ,0,natom+1, 0,3);
252     free_dvector(deriv.drb ,0,natom+1);
253    
254     free_imatrix(nonbond.iNBtype ,0,nonbond.maxnbtype+1, 0,nonbond.maxnbtype+1);
255     free_imatrix(nonbond.ipif ,0,nonbond.maxnbtype+1, 0,nonbond.maxnbtype+1);
256 gilbertke 103 free_dmatrix(nonbond.vrad ,0,nonbond.maxnbtype+1, 0,nonbond.maxnbtype+1);
257     free_dmatrix(nonbond.veps ,0,nonbond.maxnbtype+1, 0,nonbond.maxnbtype+1);
258     free_dmatrix(nonbond.vrad14 ,0,nonbond.maxnbtype+1, 0,nonbond.maxnbtype+1);
259     free_dmatrix(nonbond.veps14 ,0,nonbond.maxnbtype+1, 0,nonbond.maxnbtype+1);
260 wdelano 58
261     // ntor = 9*natom;
262     ntor = torsions.ntor;
263     if (ntor > 0)
264     {
265     free_imatrix(torsions.i14 ,0,ntor,0,4);
266     free_vector(torsions.v1 ,0,ntor);
267     free_vector(torsions.v2 ,0,ntor);
268     free_vector(torsions.v3 ,0,ntor);
269     free_vector(torsions.v4 ,0,ntor);
270     free_vector(torsions.v5 ,0,ntor);
271     free_vector(torsions.v6 ,0,ntor);
272     free_ivector(torsions.ph1 ,0,ntor);
273     free_ivector(torsions.ph2 ,0,ntor);
274     free_ivector(torsions.ph3 ,0,ntor);
275     free_ivector(torsions.ph4 ,0,ntor);
276     free_ivector(torsions.ph5 ,0,ntor);
277     free_ivector(torsions.ph6 ,0,ntor);
278     }
279    
280     free_matrix(hess.hessx, 0,natom+1, 0,3);
281     free_matrix(hess.hessy, 0,natom+1, 0,3);
282     free_matrix(hess.hessz, 0,natom+1, 0,3);
283    
284     free_ivector(attached.n13 ,0,natom+1);
285     free_ivector(attached.n14 ,0,natom+1);
286     free_imatrix(attached.i13 ,0, 20, 0,natom+1);
287     free_imatrix(attached.i14 ,0,144, 0,natom+1);
288    
289 gilbertke 103 if (use_solv())
290 wdelano 58 {
291     free_dvector(solvent.asolv ,0,natom+1);
292     free_dvector(solvent.rsolv ,0,natom+1);
293     free_dvector(solvent.rborn ,0,natom+1);
294     if (solvent.type == STILL)
295     {
296     free_dvector(solvent.vsolv ,0,natom+1);
297     free_dvector(solvent.gpol ,0,natom+1);
298     } else if (solvent.type == HCT)
299     {
300     free_dvector(solvent.shct ,0,natom+1);
301     }
302     }
303     }
304