?  search:
summary | log | graphiclog1 | graphiclog2 | commit | commitdiff | tree | refs | edit | fork
blame | history | raw | HEAD
Only issue FFT warning messages on changes
[gromacs/AngularHB.git] / src / contrib / testlr.c
blobf470e974c279eaf985b44eb302ad79aaa729a109
1 /*
2 *
3 * This source code is part of
4 *
5 * G R O M A C S
6 *
7 * GROningen MAchine for Chemical Simulations
8 *
9 * VERSION 3.2.0
10 * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
11 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
12 * Copyright (c) 2001-2004, The GROMACS development team,
13 * check out http://www.gromacs.org for more information.
14
15 * This program is free software; you can redistribute it and/or
16 * modify it under the terms of the GNU General Public License
17 * as published by the Free Software Foundation; either version 2
18 * of the License, or (at your option) any later version.
19 *
20 * If you want to redistribute modifications, please consider that
21 * scientific software is very special. Version control is crucial -
22 * bugs must be traceable. We will be happy to consider code for
23 * inclusion in the official distribution, but derived work must not
24 * be called official GROMACS. Details are found in the README & COPYING
25 * files - if they are missing, get the official version at www.gromacs.org.
26 *
27 * To help us fund GROMACS development, we humbly ask that you cite
28 * the papers on the package - you can find them in the top README file.
29 *
30 * For more info, check our website at http://www.gromacs.org
31 *
32 * And Hey:
33 * GROningen Mixture of Alchemy and Childrens' Stories
34 */
35 #include <math.h>
36 #include <string.h>
37 #include "typedefs.h"
38 #include "gromacs/math/vec.h"
39 #include "gromacs/math/units.h"
40 #include "macros.h"
41 #include "names.h"
42 #include "gromacs/utility/smalloc.h"
43 #include "gromacs/fileio/tpxio.h"
44 #include "gromacs/commandline/pargs.h"
45 #include "gromacs/fileio/writeps.h"
46 #include "copyrite.h"
47 #include "pppm.h"
48 #include "readinp.h"
49 #include "force.h"
50 #include "nrnb.h"
51 #include "coulomb.h"
52 #include "gromacs/pbcutil/mshift.h"
53 #include "poisson.h"
54 #include "mdatoms.h"
55
56 static real phi_sr(FILE *log,int nj,rvec x[],real charge[],real rc,real r1,
57 rvec box, real phi[],t_block *excl,rvec f_sr[],gmx_bool bOld)
58 {
59 int i,j,k,m,ni,i1,i2;
60 real pp,r2,R,R_1,R_2,rc2;
61 real qi,qj,vsr,eps,fscal;
62 rvec dx;
63
64 vsr = 0.0;
65 eps = ONE_4PI_EPS0;
66 rc2 = rc*rc;
67 ni=0;
68 for(i=0; (i<nj-1); i++) {
69 qi=charge[i];
70 for(j=i+1; (j<nj); j++) {
71 i1=excl->index[i];
72 i2=excl->index[i+1];
73 for(k=i1; (k<i2); k++)
74 if (excl->a[k] == j)
75 break;
76 if (k == i2) {
77 r2=calc_dx2dx(x[i],x[j],box,dx);
78 if (r2 < rc2) {
79 qj = charge[j];
80 R_1 = invsqrt(r2);
81 R_2 = R_1*R_1;
82 R = invsqrt(R_2);
83 if (bOld) {
84 fscal = old_f(R,rc,r1)*R_2;
85 pp = old_phi(R,rc,r1);
86 }
87 else {
88 fscal = new_f(R,rc)*R_2;
89 pp = new_phi(R,rc);
90 }
91 phi[i] += eps*qj*pp;
92 phi[j] += eps*qi*pp;
93 vsr += eps*qj*qi*pp;
94 for(m=0; (m<DIM); m++) {
95 f_sr[i][m] += dx[m]*fscal;
96 f_sr[j][m] -= dx[m]*fscal;
97 }
98 ni++;
99 }
100 }
101 }
102 }
103 fprintf(log,"There were %d short range interactions, vsr=%g\n",ni,vsr);
104
105 return vsr;
106 }
107
108 void calc_ener(FILE *fp,char *title,gmx_bool bHeader,int nmol,
109 int natoms,real phi[],real charge[],t_block *excl)
110 {
111 int i,i1,i2,j,k;
112 real qq,qi,vv,V,Vex,Vc,Vt;
113
114 qq = 0;
115 V = 0;
116 for(i=0; (i<natoms); i++) {
117 vv = charge[i]*phi[i];
118 V += vv;
119 qq += charge[i]*charge[i];
120 }
121 V = 0.5*V;
122 Vc = 0.5*C*ONE_4PI_EPS0*qq;
123
124 qq = 0;
125 for(i=0; (i<excl->nr); i++) {
126 i1 = excl->index[i];
127 i2 = excl->index[i+1];
128 qi = charge[i];
129 for(j=i1; (j<i2); j++) {
130 k = excl->a[j];
131 if (k != i)
132 qq+=qi*charge[k];
133 }
134 }
135 Vex = qq*0.5*C*ONE_4PI_EPS0;
136
137 Vt = V - Vc - Vex;
138
139 if (bHeader) {
140 fprintf(fp,"%12s %12s %12s %12s %12s %12s\n",
141 "","Vphi","Vself","Vexcl","Vtot","1Mol");
142
143 }
144 fprintf(fp,"%12s %12.5e %12.5e %12.5e %12.5e %12.5e\n",
145 title,V,Vc,Vex,Vt,Vt/nmol);
146 }
147
148 void write_pqr(char *fn,t_atoms *atoms,rvec x[],real phi[],real dx)
149 {
150 FILE *fp;
151 int i,rnr;
152
153 fp=gmx_fio_fopen(fn,"w");
154 for(i=0; (i<atoms->nr); i++) {
155 rnr=atoms->atom[i].resnr;
156 fprintf(fp,"%-6s%5d %-4.4s%3.3s %c%4d %8.3f%8.3f%8.3f%6.2f%6.2f\n",
157 "ATOM",(i+1),*atoms->atomname[i],*atoms->resname[rnr],' ',
158 (rnr+1) % 10000,
159 10*(dx+x[i][XX]),10*x[i][YY],10*(x[i][ZZ]),0.0,phi[i]);
160 }
161 gmx_fio_fclose(fp);
162 }
163
164 void write_grid_pqr(char *fn,int nx,int ny,int nz,real ***phi)
165 {
166 FILE *fp;
167 int i,j,k,rnr=0;
168 real fac=4.0;
169
170 fp=gmx_fio_fopen(fn,"w");
171 for(i=0; (i<nx); i++)
172 for(j=0; (j<ny); j++)
173 for(k=0; (k<nz); k++,rnr++)
174 fprintf(fp,"%-6s%5d %-4.4s%3.3s %c%4d %8.3f%8.3f%8.3f%6.2f%6.2f\n",
175 "ATOM",(i+1),"C","C",' ',
176 1+((rnr+1) % 10000),fac*i,fac*j,fac*k,0.0,phi[i][j][k]);
177 gmx_fio_fclose(fp);
178 }
179 void plot_phi(char *fn,rvec box,int natoms,rvec x[],real phi[])
180 {
181 t_psdata eps;
182 real phi_max,rr,gg,bb,fac,dx,x0,y0;
183 real offset;
184 int i;
185
186 phi_max=phi[0];
187 rr=gg=bb=0.0;
188 for(i=0; (i<natoms); i++)
189 phi_max=max(phi_max,fabs(phi[i]));
190
191 if (phi_max==0.0) {
192 fprintf(stderr,"All values zero, see .out file\n");
193 return;
194 }
195 offset=20.0;
196 fac=15.0;
197 #ifdef DEBUG
198 fprintf(stderr,"Scaling box by %g\n",fac);
199 #endif
200 eps=ps_open(fn,0,0,
201 (real)(fac*box[XX]+2*offset),(real)(fac*box[YY]+2*offset));
202 ps_translate(eps,offset,offset);
203 ps_color(eps,0,0,0);
204 ps_box(eps,1,1,(real)(fac*box[XX]-1),(real)(fac*box[YY]-1));
205 dx=0.15*fac;
206 for(i=0; (i<natoms); i++) {
207 rr=gg=bb=1.0;
208 if (phi[i] < 0)
209 gg=bb=(1.0+(phi[i]/phi_max));
210 else
211 rr=gg=(1.0-(phi[i]/phi_max));
212 rr=rgbset(rr);
213 gg=rgbset(gg);
214 bb=rgbset(bb);
215 ps_color(eps,rr,gg,bb);
216 x0=fac*x[i][XX];
217 y0=fac*x[i][YY];
218 ps_fillbox(eps,(real)(x0-dx),(real)(y0-dx),(real)(x0+dx),(real)(y0+dx));
219 }
220 ps_close(eps);
221 }
222
223 void plot_qtab(char *fn,int nx,int ny,int nz,real ***qtab)
224 {
225 rvec box;
226 rvec *xx;
227 real *phi;
228 int i,npt,ix,iy,iz;
229
230 box[XX]=nx;
231 box[YY]=ny;
232 box[ZZ]=nz;
233
234 npt=(box[XX]*box[YY]*box[ZZ]);
235 snew(xx,npt);
236 snew(phi,npt);
237 nx/=2;
238 ny/=2;
239 nz/=2;
240 i=0;
241 for(ix=-nx; (ix<nx); ix++)
242 for(iy=-ny; (iy<ny); iy++)
243 for(iz=-nz; (iz<nz); iz++,i++) {
244 xx[i][XX]=ix+nx+0.5;
245 xx[i][YY]=iy+ny+0.5;
246 xx[i][ZZ]=iz+nz+0.5; /* onzin */
247 phi[i]=qtab[ix+nx][iy+ny][iz+nz];
248 }
249
250 plot_phi(fn,box,npt,xx,phi);
251
252 sfree(xx);
253 sfree(phi);
254 }
255
256 void print_phi(char *fn,int natoms,rvec x[],real phi[])
257 {
258 FILE *fp;
259 int i;
260
261 fp=gmx_fio_fopen(fn,"w");
262 for(i=0; (i<natoms); i++)
263 fprintf(fp,"%10d %12.5e\n",i,phi[i]);
264 gmx_fio_fclose(fp);
265 }
266
267 void pr_f(char *fn,int natoms,rvec f[])
268 {
269 FILE *fp;
270 int i;
271
272 fp=gmx_fio_fopen(fn,"w");
273 for(i=0; (i<natoms); i++)
274 fprintf(fp," %12.5e\n %12.5e\n %12.5e\n",f[i][XX],f[i][YY],f[i][ZZ]);
275 gmx_fio_fclose(fp);
276 }
277
278 void test_pppm(FILE *log, gmx_bool bVerbose,
279 gmx_bool bGenerGhat, char *ghatfn,
280 t_atoms *atoms, t_inputrec *ir,
281 rvec x[], rvec f[],
282 real charge[], rvec box,
283 real phi[], real phi_s[],
284 int nmol, t_commrec *cr,
285 gmx_bool bOld, t_block *cgs)
286 {
287 char buf[256];
288 real ener;
289 int i;
290 t_nrnb nrnb;
291 t_nsborder nsb;
292
293 init_nrnb(&nrnb);
294 calc_nsb(cgs,1,&nsb,0);
295
296 /* First time only setup is done! */
297 init_pppm(log,cr,&nsb,bVerbose,bOld,box,ghatfn,ir);
298
299 ener = do_pppm(log,bVerbose,x,f,charge,box,phi,cr,&nsb,&nrnb);
300 fprintf(log,"Vpppm = %g\n",ener);
301
302 sprintf(buf,"PPPM-%d.pdb",ir->nkx);
303 write_pqr(buf,atoms,x,phi,0);
304
305 pr_f("pppm-force",atoms->nr,f);
306
307 calc_ener(log,buf,FALSE,nmol,atoms->nr,phi,charge,&atoms->excl);
308
309 for(i=0; (i<atoms->nr); i++)
310 phi[i]+=phi_s[i];
311 sprintf(buf,"PPPM-%d+SR",ir->nkx);
312 calc_ener(log,buf,FALSE,nmol,atoms->nr,phi,charge,&atoms->excl);
313 strcat(buf,".pdb");
314 write_pqr(buf,atoms,x,phi,0);
315 }
316
317 void test_poisson(FILE *log, gmx_bool bVerbose,
318 t_atoms *atoms, t_inputrec *ir,
319 rvec x[], rvec f[],
320 real charge[], rvec box,
321 real phi[], real phi_s[],
322 int nmol, t_commrec *cr,
323 gmx_bool bFour, rvec f_four[],
324 real phi_f[], gmx_bool bOld)
325 {
326 char buf[256];
327 real ener;
328 rvec beta;
329 int i,nit;
330 t_nrnb nrnb;
331
332 init_nrnb(&nrnb);
333
334 /* First time only setup is done! */
335 if (bFour) {
336 for(i=0; (i<atoms->nr); i++)
337 phi_f[i] -= phi_s[i];
338 ener = do_optimize_poisson(log,bVerbose,ir,atoms->nr,x,f,charge,
339 box,phi,cr,&nrnb,f_four,phi_f,beta,bOld);
340 for(i=0; (i<atoms->nr); i++)
341 phi_f[i] += phi_s[i];
342 nit = 0;
343 }
344 else {
345 ener = do_poisson(log,bVerbose,ir,atoms->nr,x,f,charge,box,phi,
346 cr,&nrnb,&nit,bOld);
347 }
348
349 fprintf(log,"Vpoisson = %g, nit = %d\n",ener,nit);
350
351 sprintf(buf,"POISSON-%d.pdb",ir->nkx);
352 write_pqr(buf,atoms,x,phi,0);
353
354 pr_f("poisson-force",atoms->nr,f);
355
356 calc_ener(log,buf,FALSE,nmol,atoms->nr,phi,charge,&atoms->excl);
357
358 for(i=0; (i<atoms->nr); i++)
359 phi[i]+=phi_s[i];
360 sprintf(buf,"POISSON-%d+SR",ir->nkx);
361 calc_ener(log,buf,FALSE,nmol,atoms->nr,phi,charge,&atoms->excl);
362 strcat(buf,".pdb");
363 write_pqr(buf,atoms,x,phi,0);
364 }
365
366 void test_four(FILE *log,int NFILE,t_filenm fnm[],t_atoms *atoms,
367 t_inputrec *ir,rvec x[],rvec f[],rvec box,real charge[],
368 real phi_f[],real phi_s[],int nmol,t_commrec *cr,
369 gmx_bool bOld,gmx_bool bOldEwald)
370 {
371 int i;
372 real energy;
373 ewald_tab_t et;
374
375 init_ewald_tab(&et, NULL, ir, log);
376
377 if (bOldEwald)
378 energy = do_ewald(log,ir,atoms->nr,x,f,charge,box,phi_f,cr,bOld,et);
379 else
380 energy = do_ewald_new(log,ir,atoms->nr,x,f,charge,box,phi_f,cr,bOld,et);
381
382 /*symmetrize_phi(log,atoms->nr,phi_f,bVerbose);*/
383
384 /* Plot the long range fourier solution in a matrix */
385 plot_phi(opt2fn("-elf",NFILE,fnm),box,atoms->nr,x,phi_f);
386 print_phi(opt2fn("-of",NFILE,fnm),atoms->nr,x,phi_f);
387 calc_ener(log,"Fourier",FALSE,nmol,atoms->nr,phi_f,charge,&atoms->excl);
388 write_pqr(opt2fn("-pf",NFILE,fnm),atoms,x,phi_f,0.0/*1.5*box[XX]*/);
389 pr_f("four-force",atoms->nr,f);
390
391 /* Calc and plot the total potential */
392 for(i=0; (i<atoms->nr); i++) {
393 phi_f[i]+=phi_s[i];
394 /*clear_rvec(f[i]);*/
395 }
396 calc_ener(log,"Fourier+SR",FALSE,nmol,atoms->nr,phi_f,charge,&atoms->excl);
397 }
398
399 static void print_opts(FILE *fp,t_inputrec *ir,gmx_bool bFour)
400 {
401 fprintf(fp,"Ewald solution: %s\n",gmx_bool_names[bFour]);
402 fprintf(fp,"r1: %10.3e\n",ir->rcoulomb_switch);
403 fprintf(fp,"rc: %10.3e\n",ir->rcoulomb);
404 if (bFour)
405 fprintf(fp,"KVectors: %10d %10d %10d\n",ir->nkx,ir->nky,ir->nkz);
406 fprintf(fp,"\n");
407 }
408
409 int main(int argc,char *argv[])
410 {
411 static char *desc[] = {
412 "[TT]testlr[tt] tests the PPPM and Ewald method for the",
413 "long range electrostatics problem."
414 };
415 static t_filenm fnm[] = {
416 { efTPX, NULL, NULL, ffREAD },
417 { efHAT, "-g", "ghat", ffOPTRD },
418 { efOUT, "-o", "rho", ffOPTWR },
419 { efOUT, "-op", "lr-pb", ffOPTWR },
420 { efOUT, "-of", "lr-four", ffOPTWR },
421 { efOUT, "-opt", "tot-pb", ffOPTWR },
422 { efOUT, "-oft", "tot-four", ffOPTWR },
423 { efOUT, "-fin", "lr-four", ffOPTWR },
424 { efEPS, "-es", "sr", ffOPTWR },
425 { efEPS, "-elf", "lr-four", ffOPTWR },
426 { efEPS, "-etf", "tot-four", ffOPTWR },
427 { efEPS, "-qr", "qk-real", ffOPTWR },
428 { efEPS, "-qi", "qk-im", ffOPTWR },
429 { efEPS, "-elp", "lr-pb", ffOPTWR },
430 { efEPS, "-etp", "tot-pb", ffOPTWR },
431 { efEPS, "-rho", "rho", ffOPTWR },
432 { efEPS, "-qq", "charge", ffOPTWR },
433 { efXVG, "-gt", "gk-tab", ffOPTWR },
434 { efXVG, "-fcorr","fcorr", ffWRITE },
435 { efXVG, "-pcorr","pcorr", ffWRITE },
436 { efXVG, "-ftotcorr","ftotcorr", ffWRITE },
437 { efXVG, "-ptotcorr","ptotcorr", ffWRITE },
438 { efLOG, "-l", "fptest", ffWRITE },
439 { efXVG, "-gr", "spread", ffOPTWR },
440 { efPDB, "-pf", "pqr-four", ffOPTWR },
441 { efPDB, "-phitot", "pppm-phitot", ffOPTWR }
442 };
443 #define NFILE asize(fnm)
444 FILE *log;
445 t_topology top;
446 t_tpxheader stath;
447 t_inputrec ir;
448 t_block *excl;
449 t_forcerec *fr;
450 t_commrec *cr;
451 t_mdatoms *mdatoms;
452 t_graph *graph;
453 int i,step,nre,natoms,nmol;
454 rvec *x,*f_sr,*f_excl,*f_four,*f_pppm,*f_pois,box_size,hbox;
455 matrix box;
456 real t,lambda,vsr,*charge,*phi_f,*phi_pois,*phi_s,*phi_p3m,*rho;
457 output_env_t oenv;
458
459 static gmx_bool bFour=FALSE,bVerbose=FALSE,bGGhat=FALSE,bPPPM=TRUE,
460 bPoisson=FALSE,bOld=FALSE,bOldEwald=TRUE;
461 static int nprocs = 1;
462 static t_pargs pa[] = {
463 { "-np", FALSE, etINT, &nprocs, "Do it in parallel" },
464 { "-ewald", FALSE, etBOOL, &bFour, "Do an Ewald solution"},
465 { "-pppm", FALSE, etBOOL, &bPPPM, "Do a PPPM solution" },
466 { "-poisson",FALSE, etBOOL, &bPoisson,"Do a Poisson solution" },
467 { "-v", FALSE, etBOOL, &bVerbose,"Verbose on"},
468 { "-ghat", FALSE, etBOOL, &bGGhat, "Generate Ghat function"},
469 { "-old", FALSE, etBOOL, &bOld, "Use old function types"},
470 { "-oldewald",FALSE,etBOOL, &bOldEwald,"Use old Ewald code"}
471 };
472
473 CopyRight(stderr,argv[0]);
474 parse_common_args_r(&argc,argv,PCA_CAN_TIME | PCA_CAN_VIEW,
475 NFILE,fnm,asize(pa),pa,asize(desc),desc,0,NULL,&oenv);
476
477 if (nprocs > 1) {
478 cr = init_par(&argc,argv);
479 open_log(ftp2fn(efLOG,NFILE,fnm),cr);
480 log = stdlog;
481 }
482 else {
483 cr = init_par(&argc,argv);
484 log = ftp2FILE(efLOG,NFILE,fnm,"w");
485 stdlog = log; }
486
487
488 /* Read topology and coordinates */
489 read_tpxheader(ftp2fn(efTPX,NFILE,fnm),&stath,FALSE);
490 snew(x,stath.natoms);
491 snew(f_sr,stath.natoms);
492 snew(f_excl,stath.natoms);
493 snew(f_four,stath.natoms);
494 snew(f_pppm,stath.natoms);
495 snew(f_pois,stath.natoms);
496 read_tpx(ftp2fn(efTPX,NFILE,fnm),&step,&t,&lambda,&ir,
497 box,&natoms,x,NULL,NULL,&top);
498 excl=&(top.atoms.excl);
499 nmol=top.blocks[ebMOLS].nr;
500
501 /* Allocate space for potential, charges and rho (charge density) */
502 snew(charge,stath.natoms);
503 snew(phi_f,stath.natoms);
504 snew(phi_p3m,stath.natoms);
505 snew(phi_pois,stath.natoms);
506 snew(phi_s,stath.natoms);
507 snew(rho,stath.natoms);
508
509 /* Set the charges */
510 for(i=0; (i<natoms); i++)
511 charge[i]=top.atoms.atom[i].q;
512
513 /* Make a simple box vector instead of tensor */
514 for(i=0; (i<DIM); i++)
515 box_size[i]=box[i][i];
516
517 /* Set some constants */
518 fr = mk_forcerec();
519 mdatoms = atoms2md(&(top.atoms),FALSE,FALSE);
520
521 set_LRconsts(log,ir.rcoulomb_switch,ir.rcoulomb,box_size,fr);
522 init_forcerec(log,fr,&ir,&(top.blocks[ebMOLS]),cr,
523 &(top.blocks[ebCGS]),&(top.idef),mdatoms,box,FALSE);
524 calc_shifts(box,box_size,fr->shift_vec,FALSE);
525
526 /* Periodicity stuff */
527 graph = mk_graph(&(top.idef),top.atoms.nr,FALSE,FALSE);
528 shift_self(graph,fr->shift_vec,x);
529
530 calc_LRcorrections(log,0,natoms,ir.rcoulomb_switch,
531 ir.rcoulomb,charge,excl,x,f_excl,bOld);
532 pr_f("f_excl.dat",natoms,f_excl);
533
534 /* Compute the short range potential */
535 put_atoms_in_box(natoms,box,x);
536 vsr=phi_sr(log,natoms,x,charge,ir.rcoulomb,
537 ir.rcoulomb_switch,box_size,phi_s,excl,f_sr,bOld);
538 pr_f("f_sr.dat",natoms,f_sr);
539
540 /* Plot the short range potential in a matrix */
541 calc_ener(log,"Short Range",TRUE,nmol,natoms,phi_s,charge,excl);
542
543
544 if (bFour)
545 test_four(log,NFILE,fnm,&(top.atoms),&ir,x,f_four,box_size,charge,phi_f,
546 phi_s,nmol,cr,bOld,bOldEwald);
547
548 if (bPPPM)
549 test_pppm(log,bVerbose,bGGhat,opt2fn("-g",NFILE,fnm),
550 &(top.atoms),&ir,x,f_pppm,charge,box_size,phi_p3m,phi_s,nmol,
551 cr,bOld,&(top.blocks[ebCGS]));
552
553 if (bPoisson)
554 test_poisson(log,bVerbose,
555 &(top.atoms),&ir,x,f_pois,charge,box_size,phi_pois,
556 phi_s,nmol,cr,bFour,f_four,phi_f,bOld);
557
558 if (bPPPM && bFour)
559 analyse_diff(log,"PPPM",oenv,
560 top.atoms.nr,f_four,f_pppm,phi_f,phi_p3m,phi_s,
561 opt2fn("-fcorr",NFILE,fnm),
562 opt2fn("-pcorr",NFILE,fnm),
563 opt2fn("-ftotcorr",NFILE,fnm),
564 opt2fn("-ptotcorr",NFILE,fnm));
565
566 if (bPoisson && bFour)
567 analyse_diff(log,"Poisson",oenv,
568 top.atoms.nr,f_four,f_pois,phi_f,phi_pois,phi_s,
569 opt2fn("-fcorr",NFILE,fnm),
570 opt2fn("-pcorr",NFILE,fnm),
571 opt2fn("-ftotcorr",NFILE,fnm),
572 opt2fn("-ptotcorr",NFILE,fnm));
573
574 gmx_fio_fclose(log);
575
576 gmx_thanx(stderr);
577
578 return 0;
579 }
580