ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/freemol/trunk/src/mengine/src/gmmx_run.c
Revision: 110
Committed: Thu Mar 12 01:43:50 2009 UTC (12 years, 6 months ago) by gilbertke
File size: 5770 byte(s)
Log Message:
further cleanup and localization of atom data
Line File contents
1 /* NOTICE: this source code file has been modified for use with FreeMOL */
2 #define notKEG_DEBUG
3
4
5 #define EXTERN extern
6 #include "pcwin.h"
7 #include "pcmod.h"
8 #include "utility.h"
9 #include "vibrate.h"
10
11 #include <errno.h>
12
13 struct t_logp {
14 float logp;
15 } logp_calc;
16
17 EXTERN struct t_files {
18 int nfiles, append, batch, icurrent, ibatno;
19 } files;
20
21 EXTERN struct t_minim_values {
22 int iprint, ndc, nconst;
23 float dielc;
24 } minim_values;
25 EXTERN struct t_minim_control {
26 int type, method, field, added_const;
27 char added_path[256],added_name[256];
28 } minim_control;
29
30 int rd_sdf(FILE *);
31 int setup_calculation(void);
32 void end_calculation(void);
33 void minimize(int natom,int *use,double *x,double *y,double *z);
34 static int close_contact(int);
35 void initialize(void);
36 void type(void);
37 void read_gmmxinp(char*,FILE*,int);
38 void hdel(int);
39 void hadd(void);
40 int FetchRecord(FILE *, char *);
41 void write_sdf(int);
42 void zero_data(void);
43 void type_mmx(void);
44 void read_datafiles(char *);
45 float xlogp(int natom,int *atomnum,int **iat,int **bo,long int *flags,int *type);
46 void print_energy_totals(void);
47 char *get_structure_title(void);
48 void vibrate(int natom,double *atomwt,double *x,double *y,double *z,double *charge);
49 // =====================================
50 void read_gmmxinp(char *paramfile,
51 FILE *infile,
52 int flags)
53 {
54 int i,ncount,nret;
55 int icount;
56 int ngood, nbparam, nbatom,nbcontact;
57 int batom = FALSE;
58 float logp;
59 //char *std_file;
60 char bad_atom[80],bad_contact[80],bad_param[80];
61 FILE *logfile;
62
63 // check filetype
64 ncount = 0;
65
66 /*
67 std_file = strdup(out_file);
68 strcpy(Savebox.fname, out_file);
69 */
70
71 // build filenames
72 strcpy(bad_atom,"_batom.sdf");
73 strcpy(bad_contact,"_bcontact.sdf");
74 strcpy(bad_param,"_bparam.sdf");
75
76 // open logfile
77 logfile = pcmlogfile;
78
79 files.nfiles = 0;
80 //
81 minim_control.field = MMFF94;
82 //strcpy(line,paramfile);
83 zero_data();
84 read_datafiles(paramfile);
85 minim_control.method = 3; // just do first deriv minimization
86 icount = 0;
87 //
88 ngood = nbparam = nbatom = nbcontact = 0;
89 i = 0;
90 while (1)
91 // for (i=1; i <= files.nfiles; i++)
92 {
93 ++i;
94 initialize();
95 nret = rd_sdf(infile);
96 if (-1 == nret) return ;
97 if (FALSE == nret)
98 {
99 fprintf(logfile,"Strnum: %d Cannot handle aromatic bonds. Try Kekule version.\n",i);
100 goto L_10;
101 natom = 0;
102 }
103
104 if (natom == 1)
105 {
106 fprintf(logfile,"Strnum: %d has only one atom\n",i);
107 goto L_10;
108 }
109 if (nret == TRUE)
110 {
111 // TJO
112 if (flags & DO_ADDH)
113 {
114 hdel(0);
115 hadd();
116 }
117 type();
118
119 #ifdef KEG_DEBUG
120 minim_values.iprint = TRUE;
121 // for(i=1; i<= natom; i++)
122 // printf("Atom: %d %d %x\n",i,atom.type[i],atom.flags[i]);
123 #endif
124 nret = setup_calculation();
125
126 #ifdef KEG_DEBUG
127 print_energy_totals();
128 // for(i=1; i<= natom; i++)
129 //printf("Atom: %d %d %x\n",i,atom[i].type,atom[i].flags);
130 //if (minim_values.iprint == TRUE) exit(0);
131 minim_values.iprint = FALSE;
132 #endif
133 if (nret == TRUE)
134 {
135 minimize(natom,atom.use,atom.x,atom.y,atom.z);
136 // compute xlogp
137 if (flags & DO_XLOGP) {
138 if (VERBOSE) fprintf(pcmlogfile, " Evaluating XlogP\n");
139 logp = xlogp(natom,atom.atomnum,atom.iat,atom.bo,atom.flags,atom.type);
140 logp_calc.logp = logp;
141 }
142 // compute dipole moment
143 if (flags & DO_DIPOLE) {
144 if (VERBOSE) fprintf(pcmlogfile, " Evaluating dipole moment\n");
145 charge_dipole();
146 }
147 // compute vib
148 if (flags & DO_VIBRATION) {
149 if (VERBOSE) fprintf(pcmlogfile, " Evaluating vibrational parameters\n");
150 vibrate(natom,atom.atomwt,atom.x,atom.y,atom.z,atom.charge);
151 }
152 write_sdf(flags);
153 ++ngood;
154 ++files.nfiles;
155 #ifdef KEG_DEBUG
156 print_energy_totals();
157 #endif
158 } else
159 {
160 strcpy(Savebox.fname,bad_param);
161 files.append = TRUE;
162 //write_sdf(flags);
163 nbparam++;
164 //strcpy(Savebox.fname,std_file);
165 fprintf(logfile,"Strnum: %d CID: %s fails on missing parameters\n",files.nfiles,get_structure_title());
166 if (VERBOSE)
167 fprintf(stdout, "Missing parameters - no calc - %s\n",get_structure_title());
168 }
169 end_calculation();
170 if (VERBOSE) {
171 fprintf(stdout,
172 "Strnum: %d %s natom %d \n",files.nfiles,get_structure_title(), natom);
173 }
174 } else // read_sdf failed
175 {
176 fprintf(logfile,"Strnum: %d CID: %s fails on unsupported atom\n",files.nfiles,get_structure_title());
177 strcpy(Savebox.fname,bad_atom);
178 files.append = TRUE;
179 if (batom == FALSE)
180 {
181 batom = TRUE;
182 files.append = FALSE;
183 }
184 //write_sdf(flags);
185 nbatom++;
186 //strcpy(Savebox.fname,std_file);
187 }
188 L_10:
189 continue;
190 }
191 // fclose(infile);
192 fprintf(logfile,"Ngood: %d Nbparam: %d Nbatom: %d Nbcontact: %d \n",ngood,nbparam,nbatom,nbcontact);
193 if (VERBOSE) {
194 fprintf(stdout,
195 "Ngood: %d Nbparam: %d Nbatom: %d Nbcontact: %d \n",
196 ngood,nbparam,nbatom,nbcontact);
197 }
198
199 // fclose(logfile);
200 }
201 // ================================================
202 int close_contact(int mode)
203 {
204 int i, j;
205 double dx,dy,dz;
206
207 for (i=1; i < natom; i++)
208 {
209 for (j=i+1; j <= natom; j++)
210 {
211 dx = fabs(atom.x[i] - atom.x[j]);
212 dy = fabs(atom.y[i] - atom.y[j]);
213 dz = fabs(atom.z[i] - atom.z[j]);
214 if ( (dx + dy + dz) < 0.1)
215 {
216 return FALSE;
217 }
218 }
219 }
220 return TRUE;
221 }