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