ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/freemol/trunk/src/mengine/src/gmmx_run.c
Revision: 85
Committed: Mon Jan 12 15:05:44 2009 UTC (11 years, 9 months ago) by gilbertke
File size: 5427 byte(s)
Log Message:
more cleanup and code removal
Line User Rev File contents
1 tjod 48 /* NOTICE: this source code file has been modified for use with FreeMOL */
2 tjod 3 #define EXTERN extern
3     #include "pcwin.h"
4     #include "pcmod.h"
5     #include "utility.h"
6     #include "energies.h"
7 gilbertke 63 #include "dipmom.h"
8     #include "vibrate.h"
9 tjod 3
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 gilbertke 63 static int close_contact(int);
34 tjod 3 void initialize(void);
35     void type(void);
36 tjod 16 void read_gmmxinp(char*,FILE*,int);
37 tjod 3 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 tjod 16 void read_gmmxinp(char *paramfile,
48 tjod 15 FILE *infile,
49 tjod 3 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 tjod 15 //char *std_file;
59 tjod 3 char bad_atom[80],bad_contact[80],bad_param[80];
60 tjod 15 FILE *logfile;
61 tjod 3
62     // check filetype
63     ncount = 0;
64    
65 tjod 15 /*
66 tjod 3 std_file = strdup(out_file);
67     strcpy(Savebox.fname, out_file);
68 tjod 15 */
69 tjod 3
70     // build filenames
71 tjod 15 strcpy(bad_atom,"_batom.sdf");
72     strcpy(bad_contact,"_bcontact.sdf");
73     strcpy(bad_param,"_bparam.sdf");
74 tjod 3
75     // open logfile
76 wdelano 58 logfile = pcmlogfile;
77 tjod 3
78     files.nfiles = 0;
79     //
80     minim_control.field = MMFF94;
81 tjod 25 //strcpy(line,paramfile);
82 tjod 3 zero_data();
83 tjod 25 read_datafiles(paramfile);
84 tjod 3 minim_control.method = 1; // just do first deriv minimization
85     icount = 0;
86     //
87     ngood = nbparam = nbatom = nbcontact = 0;
88 tjod 31 i = 0;
89 tjod 15 while (1)
90     // for (i=1; i <= files.nfiles; i++)
91 tjod 3 {
92 tjod 31 ++i;
93 tjod 3 initialize();
94     nret = rd_sdf(infile);
95 tjod 15 if (-1 == nret) return ;
96 gilbertke 64 if (FALSE == nret)
97     {
98 tjod 31 fprintf(logfile,"Strnum: %d Cannot handle aromatic bonds. Try Kekule version.\n",i);
99     goto L_10;
100     natom = 0;
101 gilbertke 64 }
102 tjod 3
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 tjod 20 // TJO
111 gilbertke 64 if (flags & DO_ADDH)
112     {
113 tjod 20 hdel(0);
114     hadd();
115     }
116 gilbertke 65 type();
117 tjod 3
118 gilbertke 75 // minim_values.iprint = TRUE;
119 tjod 3 nret = setup_calculation();
120 gilbertke 75 // for(i=1; i<= natom; i++)
121     //printf("Atom: %d %d %x\n",i,atom[i].type,atom[i].flags);
122     //if (minim_values.iprint == TRUE) exit(0);
123 tjod 3 if (nret == TRUE)
124     {
125     minimize();
126     // compute xlogp
127     if (energies.total < 1000.0)
128     {
129     if (flags & DO_XLOGP) {
130 wdelano 58 if (VERBOSE) fprintf(pcmlogfile, " Evaluating XlogP\n");
131 tjod 3 xlogp(&logp);
132     logp_calc.logp = logp;
133     }
134     // compute dipole moment
135     dipolemom.xdipole = 0.0;
136     dipolemom.ydipole = 0.0;
137     dipolemom.zdipole = 0.0;
138     if (flags & DO_DIPOLE) {
139 wdelano 58 if (VERBOSE) fprintf(pcmlogfile, " Evaluating dipole moment\n");
140 tjod 3 charge_dipole();
141     dipolemom.total = sqrt(dipolemom.xdipole*dipolemom.xdipole +
142     dipolemom.ydipole*dipolemom.ydipole +
143     dipolemom.zdipole*dipolemom.zdipole);
144     }
145     // compute vib
146     if (flags & DO_VIBRATION) {
147 wdelano 58 if (VERBOSE) fprintf(pcmlogfile, " Evaluating vibrational parameters\n");
148 tjod 3 vibrate();
149     }
150     write_sdf(flags);
151 tjod 15 ++ngood;
152     ++files.nfiles;
153 tjod 3 }
154     } else
155     {
156     strcpy(Savebox.fname,bad_param);
157     files.append = TRUE;
158 tjod 33 //write_sdf(flags);
159 tjod 3 nbparam++;
160 tjod 15 //strcpy(Savebox.fname,std_file);
161     fprintf(logfile,"Strnum: %d CID: %s fails on missing parameters\n",files.nfiles,Struct_Title);
162 tjod 3 if (VERBOSE)
163     fprintf(stdout, "Missing parameters - no calc - %s\n",Struct_Title);
164     }
165     end_calculation();
166     if (VERBOSE) {
167     fprintf(stdout,
168 gilbertke 63 "Strnum: %d %s natom %d \n",files.nfiles,Struct_Title, natom);
169 tjod 3 }
170     } else // read_sdf failed
171     {
172 tjod 15 fprintf(logfile,"Strnum: %d CID: %s fails on unsupported atom\n",files.nfiles,Struct_Title);
173 tjod 3 strcpy(Savebox.fname,bad_atom);
174     files.append = TRUE;
175     if (batom == FALSE)
176     {
177     batom = TRUE;
178     files.append = FALSE;
179     }
180 tjod 33 //write_sdf(flags);
181 tjod 3 nbatom++;
182 tjod 15 //strcpy(Savebox.fname,std_file);
183 tjod 3 }
184     L_10:
185     continue;
186     }
187 wdelano 58 // fclose(infile);
188 tjod 3 fprintf(logfile,"Ngood: %d Nbparam: %d Nbatom: %d Nbcontact: %d \n",ngood,nbparam,nbatom,nbcontact);
189     if (VERBOSE) {
190     fprintf(stdout,
191     "Ngood: %d Nbparam: %d Nbatom: %d Nbcontact: %d \n",
192     ngood,nbparam,nbatom,nbcontact);
193     }
194    
195 wdelano 58 // fclose(logfile);
196 tjod 3 }
197     // ================================================
198     int close_contact(int mode)
199     {
200     int i, j;
201     double dx,dy,dz;
202    
203     for (i=1; i < natom; i++)
204     {
205     for (j=i+1; j <= natom; j++)
206     {
207     dx = fabs(atom[i].x - atom[j].x);
208     dy = fabs(atom[i].y - atom[j].y);
209     dz = fabs(atom[i].z - atom[j].z);
210     if ( (dx + dy + dz) < 0.1)
211     {
212     return FALSE;
213     }
214     }
215     }
216     return TRUE;
217     }