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 |
void charge_dipole(int natom,double *x,double *y,double *z,double *atomwt,double *charge); |
50 |
void set_connectivity(void); |
51 |
// ===================================== |
52 |
void read_gmmxinp(char *paramfile, |
53 |
FILE *infile, |
54 |
int flags) |
55 |
{ |
56 |
int i,ncount,nret; |
57 |
int icount; |
58 |
int ngood, nbparam, nbatom,nbcontact; |
59 |
int batom = FALSE; |
60 |
float logp; |
61 |
//char *std_file; |
62 |
char bad_atom[80],bad_contact[80],bad_param[80]; |
63 |
FILE *logfile; |
64 |
|
65 |
// check filetype |
66 |
ncount = 0; |
67 |
|
68 |
/* |
69 |
std_file = strdup(out_file); |
70 |
strcpy(Savebox.fname, out_file); |
71 |
*/ |
72 |
|
73 |
// build filenames |
74 |
strcpy(bad_atom,"_batom.sdf"); |
75 |
strcpy(bad_contact,"_bcontact.sdf"); |
76 |
strcpy(bad_param,"_bparam.sdf"); |
77 |
|
78 |
// open logfile |
79 |
logfile = pcmlogfile; |
80 |
|
81 |
files.nfiles = 0; |
82 |
// |
83 |
minim_control.field = MMFF94; |
84 |
//strcpy(line,paramfile); |
85 |
zero_data(); |
86 |
read_datafiles(paramfile); |
87 |
minim_control.method = 3; // just do first deriv minimization |
88 |
icount = 0; |
89 |
// |
90 |
ngood = nbparam = nbatom = nbcontact = 0; |
91 |
i = 0; |
92 |
while (1) |
93 |
// for (i=1; i <= files.nfiles; i++) |
94 |
{ |
95 |
++i; |
96 |
initialize(); |
97 |
nret = rd_sdf(infile); |
98 |
if (-1 == nret) return ; |
99 |
if (FALSE == nret) |
100 |
{ |
101 |
fprintf(logfile,"Strnum: %d Cannot handle aromatic bonds. Try Kekule version.\n",i); |
102 |
goto L_10; |
103 |
natom = 0; |
104 |
} |
105 |
|
106 |
if (natom == 1) |
107 |
{ |
108 |
fprintf(logfile,"Strnum: %d has only one atom\n",i); |
109 |
goto L_10; |
110 |
} |
111 |
if (nret == TRUE) |
112 |
{ |
113 |
// TJO |
114 |
if (flags & DO_ADDH) |
115 |
{ |
116 |
hdel(0); |
117 |
hadd(); |
118 |
} |
119 |
set_connectivity(); // count number of connections and total bond order for all atoms |
120 |
type(); |
121 |
|
122 |
#ifdef KEG_DEBUG |
123 |
minim_values.iprint = TRUE; |
124 |
for(i=1; i<= natom; i++) |
125 |
fprintf(pcmlogfile,"Atom: %d %d %ld\n",i,atom.type[i],atom.flags[i]); |
126 |
#endif |
127 |
nret = setup_calculation(); |
128 |
|
129 |
#ifdef KEG_DEBUG |
130 |
print_energy_totals(); |
131 |
for(i=1; i<= natom; i++) |
132 |
fprintf(pcmlogfile,"Atom: %d %d %ld\n",i,atom.type[i],atom.flags[i]); |
133 |
if (minim_values.iprint == TRUE) exit(0); |
134 |
minim_values.iprint = FALSE; |
135 |
#endif |
136 |
if (nret == TRUE) |
137 |
{ |
138 |
minimize(natom,atom.use,atom.x,atom.y,atom.z); |
139 |
// compute xlogp |
140 |
if (flags & DO_XLOGP) { |
141 |
if (VERBOSE) fprintf(pcmlogfile, " Evaluating XlogP\n"); |
142 |
logp = xlogp(natom,atom.atomnum,atom.iat,atom.bo,atom.flags,atom.type); |
143 |
logp_calc.logp = logp; |
144 |
} |
145 |
// compute dipole moment |
146 |
if (flags & DO_DIPOLE) { |
147 |
if (VERBOSE) fprintf(pcmlogfile, " Evaluating dipole moment\n"); |
148 |
charge_dipole(natom,atom.x,atom.y,atom.z,atom.atomwt,atom.charge); |
149 |
} |
150 |
// compute vib |
151 |
if (flags & DO_VIBRATION) { |
152 |
if (VERBOSE) fprintf(pcmlogfile, " Evaluating vibrational parameters\n"); |
153 |
vibrate(natom,atom.atomwt,atom.x,atom.y,atom.z,atom.charge); |
154 |
} |
155 |
write_sdf(flags); |
156 |
++ngood; |
157 |
++files.nfiles; |
158 |
#ifdef KEG_DEBUG |
159 |
print_energy_totals(); |
160 |
#endif |
161 |
} else |
162 |
{ |
163 |
strcpy(Savebox.fname,bad_param); |
164 |
files.append = TRUE; |
165 |
//write_sdf(flags); |
166 |
nbparam++; |
167 |
//strcpy(Savebox.fname,std_file); |
168 |
fprintf(logfile,"Strnum: %d CID: %s fails on missing parameters\n",files.nfiles,get_structure_title()); |
169 |
if (VERBOSE) |
170 |
fprintf(stdout, "Missing parameters - no calc - %s\n",get_structure_title()); |
171 |
} |
172 |
end_calculation(); |
173 |
if (VERBOSE) { |
174 |
fprintf(stdout, |
175 |
"Strnum: %d %s natom %d \n",files.nfiles,get_structure_title(), natom); |
176 |
} |
177 |
} else // read_sdf failed |
178 |
{ |
179 |
fprintf(logfile,"Strnum: %d CID: %s fails on unsupported atom\n",files.nfiles,get_structure_title()); |
180 |
strcpy(Savebox.fname,bad_atom); |
181 |
files.append = TRUE; |
182 |
if (batom == FALSE) |
183 |
{ |
184 |
batom = TRUE; |
185 |
files.append = FALSE; |
186 |
} |
187 |
//write_sdf(flags); |
188 |
nbatom++; |
189 |
//strcpy(Savebox.fname,std_file); |
190 |
} |
191 |
L_10: |
192 |
continue; |
193 |
} |
194 |
// fclose(infile); |
195 |
fprintf(logfile,"Ngood: %d Nbparam: %d Nbatom: %d Nbcontact: %d \n",ngood,nbparam,nbatom,nbcontact); |
196 |
if (VERBOSE) { |
197 |
fprintf(stdout, |
198 |
"Ngood: %d Nbparam: %d Nbatom: %d Nbcontact: %d \n", |
199 |
ngood,nbparam,nbatom,nbcontact); |
200 |
} |
201 |
|
202 |
// fclose(logfile); |
203 |
} |
204 |
// ================================================ |
205 |
int close_contact(int mode) |
206 |
{ |
207 |
int i, j; |
208 |
double dx,dy,dz; |
209 |
|
210 |
for (i=1; i < natom; i++) |
211 |
{ |
212 |
for (j=i+1; j <= natom; j++) |
213 |
{ |
214 |
dx = fabs(atom.x[i] - atom.x[j]); |
215 |
dy = fabs(atom.y[i] - atom.y[j]); |
216 |
dz = fabs(atom.z[i] - atom.z[j]); |
217 |
if ( (dx + dy + dz) < 0.1) |
218 |
{ |
219 |
return FALSE; |
220 |
} |
221 |
} |
222 |
} |
223 |
return TRUE; |
224 |
} |