1 /* -*- mode: c; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4; c-file-style: "stroustrup"; -*-
4 * This source code is part of
8 * GROningen MAchine for Chemical Simulations
11 * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
12 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
13 * Copyright (c) 2001-2004, The GROMACS development team,
14 * check out http://www.gromacs.org for more information.
16 * This program is free software; you can redistribute it and/or
17 * modify it under the terms of the GNU General Public License
18 * as published by the Free Software Foundation; either version 2
19 * of the License, or (at your option) any later version.
21 * If you want to redistribute modifications, please consider that
22 * scientific software is very special. Version control is crucial -
23 * bugs must be traceable. We will be happy to consider code for
24 * inclusion in the official distribution, but derived work must not
25 * be called official GROMACS. Details are found in the README & COPYING
26 * files - if they are missing, get the official version at www.gromacs.org.
28 * To help us fund GROMACS development, we humbly ask that you cite
29 * the papers on the package - you can find them in the top README file.
31 * For more info, check our website at http://www.gromacs.org
34 * GROningen Mixture of Alchemy and Childrens' Stories
40 /* This file is completely threadsafe - please keep it that way! */
42 #include <thread_mpi.h>
55 int pr_indent(FILE *fp,int n)
59 for (i=0; i<n; i++) (void) fprintf(fp," ");
63 int available(FILE *fp,void *p,int indent,const char *title)
68 (void) fprintf(fp,"%s: not available\n",title);
73 int pr_title(FILE *fp,int indent,const char *title)
75 (void) pr_indent(fp,indent);
76 (void) fprintf(fp,"%s:\n",title);
77 return (indent+INDENT);
80 int pr_title_n(FILE *fp,int indent,const char *title,int n)
82 (void) pr_indent(fp,indent);
83 (void) fprintf(fp,"%s (%d):\n",title,n);
84 return (indent+INDENT);
87 int pr_title_nxn(FILE *fp,int indent,const char *title,int n1,int n2)
89 (void) pr_indent(fp,indent);
90 (void) fprintf(fp,"%s (%dx%d):\n",title,n1,n2);
91 return (indent+INDENT);
94 void pr_ivec(FILE *fp,int indent,const char *title,int vec[],int n, gmx_bool bShowNumbers)
98 if (available(fp,vec,indent,title))
100 indent=pr_title_n(fp,indent,title,n);
103 (void) pr_indent(fp,indent);
104 (void) fprintf(fp,"%s[%d]=%d\n",title,bShowNumbers?i:-1,vec[i]);
109 void pr_ivec_block(FILE *fp,int indent,const char *title,int vec[],int n, gmx_bool bShowNumbers)
113 if (available(fp,vec,indent,title))
115 indent=pr_title_n(fp,indent,title,n);
120 while (j < n && vec[j] == vec[j-1]+1)
124 /* Print consecutive groups of 3 or more as blocks */
129 (void) pr_indent(fp,indent);
130 (void) fprintf(fp,"%s[%d]=%d\n",
131 title,bShowNumbers?i:-1,vec[i]);
137 (void) pr_indent(fp,indent);
138 (void) fprintf(fp,"%s[%d,...,%d] = {%d,...,%d}\n",
149 void pr_bvec(FILE *fp,int indent,const char *title,gmx_bool vec[],int n, gmx_bool bShowNumbers)
153 if (available(fp,vec,indent,title))
155 indent=pr_title_n(fp,indent,title,n);
158 (void) pr_indent(fp,indent);
159 (void) fprintf(fp,"%s[%d]=%s\n",title,bShowNumbers?i:-1,
165 void pr_ivecs(FILE *fp,int indent,const char *title,ivec vec[],int n, gmx_bool bShowNumbers)
169 if (available(fp,vec,indent,title))
171 indent=pr_title_nxn(fp,indent,title,n,DIM);
174 (void) pr_indent(fp,indent);
175 (void) fprintf(fp,"%s[%d]={",title,bShowNumbers?i:-1);
176 for (j=0; j<DIM; j++)
178 if (j!=0) (void) fprintf(fp,", ");
179 fprintf(fp,"%d",vec[i][j]);
181 (void) fprintf(fp,"}\n");
186 void pr_rvec(FILE *fp,int indent,const char *title,real vec[],int n, gmx_bool bShowNumbers)
190 if (available(fp,vec,indent,title))
192 indent=pr_title_n(fp,indent,title,n);
195 pr_indent(fp,indent);
196 fprintf(fp,"%s[%d]=%12.5e\n",title,bShowNumbers?i:-1,vec[i]);
201 void pr_dvec(FILE *fp,int indent,const char *title,double vec[],int n, gmx_bool bShowNumbers)
205 if (available(fp,vec,indent,title))
207 indent=pr_title_n(fp,indent,title,n);
210 pr_indent(fp,indent);
211 fprintf(fp,"%s[%d]=%12.5e\n",title,bShowNumbers?i:-1,vec[i]);
218 void pr_mat(FILE *fp,int indent,char *title,matrix m)
222 if (available(fp,m,indent,title)) {
223 indent=pr_title_n(fp,indent,title,n);
225 pr_indent(fp,indent);
226 fprintf(fp,"%s[%d]=%12.5e %12.5e %12.5e\n",
227 title,bShowNumbers?i:-1,m[i][XX],m[i][YY],m[i][ZZ]);
233 void pr_rvecs_len(FILE *fp,int indent,const char *title,rvec vec[],int n)
237 if (available(fp,vec,indent,title)) {
238 indent=pr_title_nxn(fp,indent,title,n,DIM);
239 for (i=0; i<n; i++) {
240 (void) pr_indent(fp,indent);
241 (void) fprintf(fp,"%s[%5d]={",title,i);
242 for (j=0; j<DIM; j++) {
244 (void) fprintf(fp,", ");
245 (void) fprintf(fp,"%12.5e",vec[i][j]);
247 (void) fprintf(fp,"} len=%12.5e\n",norm(vec[i]));
252 void pr_rvecs(FILE *fp,int indent,const char *title,rvec vec[],int n)
254 const char *fshort = "%12.5e";
255 const char *flong = "%15.8e";
259 if (getenv("LONGFORMAT") != NULL)
264 if (available(fp,vec,indent,title)) {
265 indent=pr_title_nxn(fp,indent,title,n,DIM);
266 for (i=0; i<n; i++) {
267 (void) pr_indent(fp,indent);
268 (void) fprintf(fp,"%s[%5d]={",title,i);
269 for (j=0; j<DIM; j++) {
271 (void) fprintf(fp,", ");
272 (void) fprintf(fp,format,vec[i][j]);
274 (void) fprintf(fp,"}\n");
280 void pr_reals(FILE *fp,int indent,const char *title,real *vec,int n)
284 if (available(fp,vec,indent,title)) {
285 (void) pr_indent(fp,indent);
286 (void) fprintf(fp,"%s:\t",title);
288 fprintf(fp," %10g",vec[i]);
289 (void) fprintf(fp,"\n");
293 void pr_doubles(FILE *fp,int indent,const char *title,double *vec,int n)
297 if (available(fp,vec,indent,title)) {
298 (void) pr_indent(fp,indent);
299 (void) fprintf(fp,"%s:\t",title);
301 fprintf(fp," %10g",vec[i]);
302 (void) fprintf(fp,"\n");
306 static void pr_int(FILE *fp,int indent,const char *title,int i)
308 pr_indent(fp,indent);
309 fprintf(fp,"%-20s = %d\n",title,i);
312 static void pr_gmx_large_int(FILE *fp,int indent,const char *title,gmx_large_int_t i)
314 char buf[STEPSTRSIZE];
316 pr_indent(fp,indent);
317 fprintf(fp,"%-20s = %s\n",title,gmx_step_str(i,buf));
320 static void pr_real(FILE *fp,int indent,const char *title,real r)
322 pr_indent(fp,indent);
323 fprintf(fp,"%-20s = %g\n",title,r);
326 static void pr_double(FILE *fp,int indent,const char *title,double d)
328 pr_indent(fp,indent);
329 fprintf(fp,"%-20s = %g\n",title,d);
332 static void pr_str(FILE *fp,int indent,const char *title,const char *s)
334 pr_indent(fp,indent);
335 fprintf(fp,"%-20s = %s\n",title,s);
338 void pr_qm_opts(FILE *fp,int indent,const char *title,t_grpopts *opts)
342 fprintf(fp,"%s:\n",title);
344 pr_int(fp,indent,"ngQM",opts->ngQM);
345 if (opts->ngQM > 0) {
346 pr_ivec(fp,indent,"QMmethod",opts->QMmethod,opts->ngQM,FALSE);
347 pr_ivec(fp,indent,"QMbasis",opts->QMbasis,opts->ngQM,FALSE);
348 pr_ivec(fp,indent,"QMcharge",opts->QMcharge,opts->ngQM,FALSE);
349 pr_ivec(fp,indent,"QMmult",opts->QMmult,opts->ngQM,FALSE);
350 pr_bvec(fp,indent,"bSH",opts->bSH,opts->ngQM,FALSE);
351 pr_ivec(fp,indent,"CASorbitals",opts->CASorbitals,opts->ngQM,FALSE);
352 pr_ivec(fp,indent,"CASelectrons",opts->CASelectrons,opts->ngQM,FALSE);
353 pr_rvec(fp,indent,"SAon",opts->SAon,opts->ngQM,FALSE);
354 pr_rvec(fp,indent,"SAon",opts->SAon,opts->ngQM,FALSE);
355 pr_ivec(fp,indent,"SAsteps",opts->SAsteps,opts->ngQM,FALSE);
356 pr_bvec(fp,indent,"bOPT",opts->bOPT,opts->ngQM,FALSE);
357 pr_bvec(fp,indent,"bTS",opts->bTS,opts->ngQM,FALSE);
361 static void pr_grp_opts(FILE *out,int indent,const char *title,t_grpopts *opts,
367 fprintf(out,"%s:\n",title);
369 pr_indent(out,indent);
370 fprintf(out,"nrdf%s",bMDPformat ? " = " : ":");
371 for(i=0; (i<opts->ngtc); i++)
372 fprintf(out," %10g",opts->nrdf[i]);
375 pr_indent(out,indent);
376 fprintf(out,"ref-t%s",bMDPformat ? " = " : ":");
377 for(i=0; (i<opts->ngtc); i++)
378 fprintf(out," %10g",opts->ref_t[i]);
381 pr_indent(out,indent);
382 fprintf(out,"tau-t%s",bMDPformat ? " = " : ":");
383 for(i=0; (i<opts->ngtc); i++)
384 fprintf(out," %10g",opts->tau_t[i]);
387 /* Pretty-print the simulated annealing info */
388 fprintf(out,"anneal%s",bMDPformat ? " = " : ":");
389 for(i=0; (i<opts->ngtc); i++)
390 fprintf(out," %10s",EANNEAL(opts->annealing[i]));
393 fprintf(out,"ann-npoints%s",bMDPformat ? " = " : ":");
394 for(i=0; (i<opts->ngtc); i++)
395 fprintf(out," %10d",opts->anneal_npoints[i]);
398 for(i=0; (i<opts->ngtc); i++) {
399 if(opts->anneal_npoints[i]>0) {
400 fprintf(out,"ann. times [%d]:\t",i);
401 for(j=0; (j<opts->anneal_npoints[i]); j++)
402 fprintf(out," %10.1f",opts->anneal_time[i][j]);
404 fprintf(out,"ann. temps [%d]:\t",i);
405 for(j=0; (j<opts->anneal_npoints[i]); j++)
406 fprintf(out," %10.1f",opts->anneal_temp[i][j]);
411 pr_indent(out,indent);
412 fprintf(out,"acc:\t");
413 for(i=0; (i<opts->ngacc); i++)
414 for(m=0; (m<DIM); m++)
415 fprintf(out," %10g",opts->acc[i][m]);
418 pr_indent(out,indent);
419 fprintf(out,"nfreeze:");
420 for(i=0; (i<opts->ngfrz); i++)
421 for(m=0; (m<DIM); m++)
422 fprintf(out," %10s",opts->nFreeze[i][m] ? "Y" : "N");
426 for(i=0; (i<opts->ngener); i++) {
427 pr_indent(out,indent);
428 fprintf(out,"energygrp-flags[%3d]:",i);
429 for(m=0; (m<opts->ngener); m++)
430 fprintf(out," %d",opts->egp_flags[opts->ngener*i+m]);
437 static void pr_matrix(FILE *fp,int indent,const char *title,rvec *m,
441 fprintf(fp,"%-10s = %g %g %g %g %g %g\n",title,
442 m[XX][XX],m[YY][YY],m[ZZ][ZZ],m[XX][YY],m[XX][ZZ],m[YY][ZZ]);
444 pr_rvecs(fp,indent,title,m,DIM);
447 static void pr_cosine(FILE *fp,int indent,const char *title,t_cosines *cos,
453 fprintf(fp,"%s = %d\n",title,cos->n);
456 indent=pr_title(fp,indent,title);
457 (void) pr_indent(fp,indent);
458 fprintf(fp,"n = %d\n",cos->n);
460 (void) pr_indent(fp,indent+2);
462 for(j=0; (j<cos->n); j++)
463 fprintf(fp," %e",cos->a[j]);
465 (void) pr_indent(fp,indent+2);
467 for(j=0; (j<cos->n); j++)
468 fprintf(fp," %e",cos->phi[j]);
474 #define PS(t,s) pr_str(fp,indent,t,s)
475 #define PI(t,s) pr_int(fp,indent,t,s)
476 #define PSTEP(t,s) pr_gmx_large_int(fp,indent,t,s)
477 #define PR(t,s) pr_real(fp,indent,t,s)
478 #define PD(t,s) pr_double(fp,indent,t,s)
480 static void pr_pullgrp(FILE *fp,int indent,int g,t_pullgrp *pg)
482 pr_indent(fp,indent);
483 fprintf(fp,"pull-group %d:\n",g);
485 pr_ivec_block(fp,indent,"atom",pg->ind,pg->nat,TRUE);
486 pr_rvec(fp,indent,"weight",pg->weight,pg->nweight,TRUE);
487 PI("pbcatom",pg->pbcatom);
488 pr_rvec(fp,indent,"vec",pg->vec,DIM,TRUE);
489 pr_rvec(fp,indent,"init",pg->init,DIM,TRUE);
495 static void pr_pull(FILE *fp,int indent,t_pull *pull)
499 PS("pull-geometry",EPULLGEOM(pull->eGeom));
500 pr_ivec(fp,indent,"pull-dim",pull->dim,DIM,TRUE);
501 PR("pull-r1",pull->cyl_r1);
502 PR("pull-r0",pull->cyl_r0);
503 PR("pull-constr-tol",pull->constr_tol);
504 PI("pull-nstxout",pull->nstxout);
505 PI("pull-nstfout",pull->nstfout);
506 PI("pull-ngrp",pull->ngrp);
507 for(g=0; g<pull->ngrp+1; g++)
508 pr_pullgrp(fp,indent,g,&pull->grp[g]);
511 static void pr_rotgrp(FILE *fp,int indent,int g,t_rotgrp *rotg)
513 pr_indent(fp,indent);
514 fprintf(fp,"rotation_group %d:\n",g);
516 PS("type",EROTGEOM(rotg->eType));
517 PS("massw",BOOL(rotg->bMassW));
518 pr_ivec_block(fp,indent,"atom",rotg->ind,rotg->nat,TRUE);
519 pr_rvecs(fp,indent,"x_ref",rotg->x_ref,rotg->nat);
520 pr_rvec(fp,indent,"vec",rotg->vec,DIM,TRUE);
521 pr_rvec(fp,indent,"pivot",rotg->pivot,DIM,TRUE);
522 PR("rate",rotg->rate);
524 PR("slab_dist",rotg->slab_dist);
525 PR("min_gaussian",rotg->min_gaussian);
526 PR("epsilon",rotg->eps);
527 PS("fit_method",EROTFIT(rotg->eFittype));
528 PI("potfitangle_nstep",rotg->PotAngle_nstep);
529 PR("potfitangle_step",rotg->PotAngle_step);
532 static void pr_rot(FILE *fp,int indent,t_rot *rot)
536 PI("rot_nstrout",rot->nstrout);
537 PI("rot_nstsout",rot->nstsout);
538 PI("rot_ngrp",rot->ngrp);
539 for(g=0; g<rot->ngrp; g++)
540 pr_rotgrp(fp,indent,g,&rot->grp[g]);
543 void pr_inputrec(FILE *fp,int indent,const char *title,t_inputrec *ir,
546 const char *infbuf="inf";
549 if (available(fp,ir,indent,title)) {
551 indent=pr_title(fp,indent,title);
552 PS("integrator",EI(ir->eI));
553 PSTEP("nsteps",ir->nsteps);
554 PSTEP("init-step",ir->init_step);
555 PS("ns-type",ENS(ir->ns_type));
556 PI("nstlist",ir->nstlist);
557 PI("ndelta",ir->ndelta);
558 PI("nstcomm",ir->nstcomm);
559 PS("comm-mode",ECOM(ir->comm_mode));
560 PI("nstlog",ir->nstlog);
561 PI("nstxout",ir->nstxout);
562 PI("nstvout",ir->nstvout);
563 PI("nstfout",ir->nstfout);
564 PI("nstcalcenergy",ir->nstcalcenergy);
565 PI("nstenergy",ir->nstenergy);
566 PI("nstxtcout",ir->nstxtcout);
567 PR("init-t",ir->init_t);
568 PR("delta-t",ir->delta_t);
570 PR("xtcprec",ir->xtcprec);
574 PI("pme-order",ir->pme_order);
575 PR("ewald-rtol",ir->ewald_rtol);
576 PR("ewald-geometry",ir->ewald_geometry);
577 PR("epsilon-surface",ir->epsilon_surface);
578 PS("optimize-fft",BOOL(ir->bOptFFT));
579 PS("ePBC",EPBC(ir->ePBC));
580 PS("bPeriodicMols",BOOL(ir->bPeriodicMols));
581 PS("bContinuation",BOOL(ir->bContinuation));
582 PS("bShakeSOR",BOOL(ir->bShakeSOR));
583 PS("etc",ETCOUPLTYPE(ir->etc));
584 PI("nsttcouple",ir->nsttcouple);
585 PS("epc",EPCOUPLTYPE(ir->epc));
586 PS("epctype",EPCOUPLTYPETYPE(ir->epct));
587 PI("nstpcouple",ir->nstpcouple);
588 PR("tau-p",ir->tau_p);
589 pr_matrix(fp,indent,"ref-p",ir->ref_p,bMDPformat);
590 pr_matrix(fp,indent,"compress",ir->compress,bMDPformat);
591 PS("refcoord-scaling",EREFSCALINGTYPE(ir->refcoord_scaling));
593 fprintf(fp,"posres-com = %g %g %g\n",ir->posres_com[XX],
594 ir->posres_com[YY],ir->posres_com[ZZ]);
596 pr_rvec(fp,indent,"posres-com",ir->posres_com,DIM,TRUE);
598 fprintf(fp,"posres-comB = %g %g %g\n",ir->posres_comB[XX],
599 ir->posres_comB[YY],ir->posres_comB[ZZ]);
601 pr_rvec(fp,indent,"posres-comB",ir->posres_comB,DIM,TRUE);
602 PI("andersen-seed",ir->andersen_seed);
603 PR("rlist",ir->rlist);
604 PR("rlistlong",ir->rlistlong);
606 PS("coulombtype",EELTYPE(ir->coulombtype));
607 PR("rcoulomb-switch",ir->rcoulomb_switch);
608 PR("rcoulomb",ir->rcoulomb);
609 PS("vdwtype",EVDWTYPE(ir->vdwtype));
610 PR("rvdw-switch",ir->rvdw_switch);
612 if (ir->epsilon_r != 0)
613 PR("epsilon-r",ir->epsilon_r);
615 PS("epsilon-r",infbuf);
616 if (ir->epsilon_rf != 0)
617 PR("epsilon-rf",ir->epsilon_rf);
619 PS("epsilon-rf",infbuf);
620 PR("tabext",ir->tabext);
621 PS("implicit-solvent",EIMPLICITSOL(ir->implicit_solvent));
622 PS("gb-algorithm",EGBALGORITHM(ir->gb_algorithm));
623 PR("gb-epsilon-solvent",ir->gb_epsilon_solvent);
624 PI("nstgbradii",ir->nstgbradii);
625 PR("rgbradii",ir->rgbradii);
626 PR("gb-saltconc",ir->gb_saltconc);
627 PR("gb-obc-alpha",ir->gb_obc_alpha);
628 PR("gb-obc-beta",ir->gb_obc_beta);
629 PR("gb-obc-gamma",ir->gb_obc_gamma);
630 PR("gb-dielectric-offset",ir->gb_dielectric_offset);
631 PS("sa-algorithm",ESAALGORITHM(ir->gb_algorithm));
632 PR("sa-surface-tension",ir->sa_surface_tension);
634 PS("DispCorr",EDISPCORR(ir->eDispCorr));
635 PS("free-energy",EFEPTYPE(ir->efep));
636 PR("init-lambda",ir->init_lambda);
637 PR("delta-lambda",ir->delta_lambda);
640 PI("n-foreign-lambda",ir->n_flambda);
642 if (ir->n_flambda > 0)
644 pr_indent(fp,indent);
645 fprintf(fp,"foreign-lambda%s",bMDPformat ? " = " : ":");
646 for(i=0; i<ir->n_flambda; i++)
648 fprintf(fp," %10g",ir->flambda[i]);
652 PR("sc-alpha",ir->sc_alpha);
653 PI("sc-power",ir->sc_power);
654 PR("sc-sigma",ir->sc_sigma);
655 PR("sc-sigma-min",ir->sc_sigma_min);
656 PI("nstdhdl", ir->nstdhdl);
657 PS("separate-dhdl-file", SEPDHDLFILETYPE(ir->separate_dhdl_file));
658 PS("dhdl-derivatives", DHDLDERIVATIVESTYPE(ir->dhdl_derivatives));
659 PI("dh-hist-size", ir->dh_hist_size);
660 PD("dh-hist-spacing", ir->dh_hist_spacing);
662 PI("nwall",ir->nwall);
663 PS("wall-type",EWALLTYPE(ir->wall_type));
664 PI("wall-atomtype[0]",ir->wall_atomtype[0]);
665 PI("wall-atomtype[1]",ir->wall_atomtype[1]);
666 PR("wall-density[0]",ir->wall_density[0]);
667 PR("wall-density[1]",ir->wall_density[1]);
668 PR("wall-ewald-zfac",ir->wall_ewald_zfac);
670 PS("pull",EPULLTYPE(ir->ePull));
671 if (ir->ePull != epullNO)
672 pr_pull(fp,indent,ir->pull);
674 PS("rotation",BOOL(ir->bRot));
676 pr_rot(fp,indent,ir->rot);
678 PS("disre",EDISRETYPE(ir->eDisre));
679 PS("disre-weighting",EDISREWEIGHTING(ir->eDisreWeighting));
680 PS("disre-mixed",BOOL(ir->bDisreMixed));
681 PR("dr-fc",ir->dr_fc);
682 PR("dr-tau",ir->dr_tau);
683 PR("nstdisreout",ir->nstdisreout);
684 PR("orires-fc",ir->orires_fc);
685 PR("orires-tau",ir->orires_tau);
686 PR("nstorireout",ir->nstorireout);
688 PR("dihre-fc",ir->dihre_fc);
690 PR("em-stepsize",ir->em_stepsize);
691 PR("em-tol",ir->em_tol);
692 PI("niter",ir->niter);
693 PR("fc-stepsize",ir->fc_stepsize);
694 PI("nstcgsteep",ir->nstcgsteep);
695 PI("nbfgscorr",ir->nbfgscorr);
697 PS("ConstAlg",ECONSTRTYPE(ir->eConstrAlg));
698 PR("shake-tol",ir->shake_tol);
699 PI("lincs-order",ir->nProjOrder);
700 PR("lincs-warnangle",ir->LincsWarnAngle);
701 PI("lincs-iter",ir->nLincsIter);
702 PR("bd-fric",ir->bd_fric);
703 PI("ld-seed",ir->ld_seed);
704 PR("cos-accel",ir->cos_accel);
705 pr_matrix(fp,indent,"deform",ir->deform,bMDPformat);
707 PS("adress",BOOL(ir->bAdress));
709 PS("adress_type",EADRESSTYPE(ir->adress->type));
710 PR("adress_const_wf",ir->adress->const_wf);
711 PR("adress_ex_width",ir->adress->ex_width);
712 PR("adress_hy_width",ir->adress->hy_width);
713 PS("adress_interface_correction",EADRESSICTYPE(ir->adress->icor));
714 PS("adress_site",EADRESSSITETYPE(ir->adress->site));
715 PR("adress_ex_force_cap",ir->adress->ex_forcecap);
716 PS("adress_do_hybridpairs", BOOL(ir->adress->do_hybridpairs));
718 pr_rvec(fp,indent,"adress_reference_coords",ir->adress->refs,DIM,TRUE);
720 PI("userint1",ir->userint1);
721 PI("userint2",ir->userint2);
722 PI("userint3",ir->userint3);
723 PI("userint4",ir->userint4);
724 PR("userreal1",ir->userreal1);
725 PR("userreal2",ir->userreal2);
726 PR("userreal3",ir->userreal3);
727 PR("userreal4",ir->userreal4);
728 pr_grp_opts(fp,indent,"grpopts",&(ir->opts),bMDPformat);
729 pr_cosine(fp,indent,"efield-x",&(ir->ex[XX]),bMDPformat);
730 pr_cosine(fp,indent,"efield-xt",&(ir->et[XX]),bMDPformat);
731 pr_cosine(fp,indent,"efield-y",&(ir->ex[YY]),bMDPformat);
732 pr_cosine(fp,indent,"efield-yt",&(ir->et[YY]),bMDPformat);
733 pr_cosine(fp,indent,"efield-z",&(ir->ex[ZZ]),bMDPformat);
734 pr_cosine(fp,indent,"efield-zt",&(ir->et[ZZ]),bMDPformat);
735 PS("bQMMM",BOOL(ir->bQMMM));
736 PI("QMconstraints",ir->QMconstraints);
737 PI("QMMMscheme",ir->QMMMscheme);
738 PR("scalefactor",ir->scalefactor);
739 pr_qm_opts(fp,indent,"qm-opts",&(ir->opts));
746 static void pr_harm(FILE *fp,t_iparams *iparams,const char *r,const char *kr)
748 fprintf(fp,"%sA=%12.5e, %sA=%12.5e, %sB=%12.5e, %sB=%12.5e\n",
749 r,iparams->harmonic.rA,kr,iparams->harmonic.krA,
750 r,iparams->harmonic.rB,kr,iparams->harmonic.krB);
753 void pr_iparams(FILE *fp,t_functype ftype,t_iparams *iparams)
756 real VA[4],VB[4],*rbcA,*rbcB;
761 pr_harm(fp,iparams,"th","ct");
763 case F_CROSS_BOND_BONDS:
764 fprintf(fp,"r1e=%15.8e, r2e=%15.8e, krr=%15.8e\n",
765 iparams->cross_bb.r1e,iparams->cross_bb.r2e,
766 iparams->cross_bb.krr);
768 case F_CROSS_BOND_ANGLES:
769 fprintf(fp,"r1e=%15.8e, r1e=%15.8e, r3e=%15.8e, krt=%15.8e\n",
770 iparams->cross_ba.r1e,iparams->cross_ba.r2e,
771 iparams->cross_ba.r3e,iparams->cross_ba.krt);
774 fprintf(fp,"theta=%15.8e, ktheta=%15.8e, r13=%15.8e, kUB=%15.8e\n",
775 iparams->u_b.theta,iparams->u_b.ktheta,iparams->u_b.r13,iparams->u_b.kUB);
777 case F_QUARTIC_ANGLES:
778 fprintf(fp,"theta=%15.8e",iparams->qangle.theta);
780 fprintf(fp,", c%c=%15.8e",'0'+i,iparams->qangle.c[i]);
784 fprintf(fp,"a=%15.8e, b=%15.8e, c=%15.8e\n",
785 iparams->bham.a,iparams->bham.b,iparams->bham.c);
790 pr_harm(fp,iparams,"b0","cb");
793 pr_harm(fp,iparams,"xi","cx");
796 fprintf(fp,"b0=%15.8e, cb=%15.8e, beta=%15.8e\n",
797 iparams->morse.b0,iparams->morse.cb,iparams->morse.beta);
800 fprintf(fp,"b0=%15.8e, kb=%15.8e, kcub=%15.8e\n",
801 iparams->cubic.b0,iparams->cubic.kb,iparams->cubic.kcub);
807 fprintf(fp,"bm=%15.8e, kb=%15.8e\n",iparams->fene.bm,iparams->fene.kb);
810 fprintf(fp,"lowA=%15.8e, up1A=%15.8e, up2A=%15.8e, kA=%15.8e, lowB=%15.8e, up1B=%15.8e, up2B=%15.8e, kB=%15.8e,\n",
811 iparams->restraint.lowA,iparams->restraint.up1A,
812 iparams->restraint.up2A,iparams->restraint.kA,
813 iparams->restraint.lowB,iparams->restraint.up1B,
814 iparams->restraint.up2B,iparams->restraint.kB);
820 fprintf(fp,"tab=%d, kA=%15.8e, kB=%15.8e\n",
821 iparams->tab.table,iparams->tab.kA,iparams->tab.kB);
824 fprintf(fp,"alpha=%15.8e\n",iparams->polarize.alpha);
827 fprintf(fp,"a=%15.8e, alpha1=%15.8e, alpha2=%15.8e, rfac=%15.8e\n",
828 iparams->thole.a,iparams->thole.alpha1,iparams->thole.alpha2,
829 iparams->thole.rfac);
832 fprintf(fp,"al_x=%15.8e, al_y=%15.8e, al_z=%15.8e, rOH=%9.6f, rHH=%9.6f, rOD=%9.6f\n",
833 iparams->wpol.al_x,iparams->wpol.al_y,iparams->wpol.al_z,
834 iparams->wpol.rOH,iparams->wpol.rHH,iparams->wpol.rOD);
837 fprintf(fp,"c6=%15.8e, c12=%15.8e\n",iparams->lj.c6,iparams->lj.c12);
840 fprintf(fp,"c6A=%15.8e, c12A=%15.8e, c6B=%15.8e, c12B=%15.8e\n",
841 iparams->lj14.c6A,iparams->lj14.c12A,
842 iparams->lj14.c6B,iparams->lj14.c12B);
845 fprintf(fp,"fqq=%15.8e, qi=%15.8e, qj=%15.8e, c6=%15.8e, c12=%15.8e\n",
847 iparams->ljc14.qi,iparams->ljc14.qj,
848 iparams->ljc14.c6,iparams->ljc14.c12);
851 fprintf(fp,"qi=%15.8e, qj=%15.8e, c6=%15.8e, c12=%15.8e\n",
852 iparams->ljcnb.qi,iparams->ljcnb.qj,
853 iparams->ljcnb.c6,iparams->ljcnb.c12);
859 fprintf(fp,"phiA=%15.8e, cpA=%15.8e, phiB=%15.8e, cpB=%15.8e, mult=%d\n",
860 iparams->pdihs.phiA,iparams->pdihs.cpA,
861 iparams->pdihs.phiB,iparams->pdihs.cpB,
862 iparams->pdihs.mult);
865 fprintf(fp,"label=%4d, type=%1d, low=%15.8e, up1=%15.8e, up2=%15.8e, fac=%15.8e)\n",
866 iparams->disres.label,iparams->disres.type,
867 iparams->disres.low,iparams->disres.up1,
868 iparams->disres.up2,iparams->disres.kfac);
871 fprintf(fp,"ex=%4d, label=%d, power=%4d, c=%15.8e, obs=%15.8e, kfac=%15.8e)\n",
872 iparams->orires.ex,iparams->orires.label,iparams->orires.power,
873 iparams->orires.c,iparams->orires.obs,iparams->orires.kfac);
876 fprintf(fp,"label=%d, power=%4d phi=%15.8e, dphi=%15.8e, kfac=%15.8e)\n",
877 iparams->dihres.label,iparams->dihres.power,
878 iparams->dihres.phi,iparams->dihres.dphi,iparams->dihres.kfac);
881 fprintf(fp,"pos0A=(%15.8e,%15.8e,%15.8e), fcA=(%15.8e,%15.8e,%15.8e), pos0B=(%15.8e,%15.8e,%15.8e), fcB=(%15.8e,%15.8e,%15.8e)\n",
882 iparams->posres.pos0A[XX],iparams->posres.pos0A[YY],
883 iparams->posres.pos0A[ZZ],iparams->posres.fcA[XX],
884 iparams->posres.fcA[YY],iparams->posres.fcA[ZZ],
885 iparams->posres.pos0B[XX],iparams->posres.pos0B[YY],
886 iparams->posres.pos0B[ZZ],iparams->posres.fcB[XX],
887 iparams->posres.fcB[YY],iparams->posres.fcB[ZZ]);
890 for (i=0; i<NR_RBDIHS; i++)
891 fprintf(fp,"%srbcA[%d]=%15.8e",i==0?"":", ",i,iparams->rbdihs.rbcA[i]);
893 for (i=0; i<NR_RBDIHS; i++)
894 fprintf(fp,"%srbcB[%d]=%15.8e",i==0?"":", ",i,iparams->rbdihs.rbcB[i]);
898 /* Use the OPLS -> Ryckaert-Bellemans formula backwards to get the
899 * OPLS potential constants back.
901 rbcA = iparams->rbdihs.rbcA;
902 rbcB = iparams->rbdihs.rbcB;
904 VA[3] = -0.25*rbcA[4];
905 VA[2] = -0.5*rbcA[3];
906 VA[1] = 4.0*VA[3]-rbcA[2];
907 VA[0] = 3.0*VA[2]-2.0*rbcA[1];
909 VB[3] = -0.25*rbcB[4];
910 VB[2] = -0.5*rbcB[3];
911 VB[1] = 4.0*VB[3]-rbcB[2];
912 VB[0] = 3.0*VB[2]-2.0*rbcB[1];
914 for (i=0; i<NR_FOURDIHS; i++)
915 fprintf(fp,"%sFourA[%d]=%15.8e",i==0?"":", ",i,VA[i]);
917 for (i=0; i<NR_FOURDIHS; i++)
918 fprintf(fp,"%sFourB[%d]=%15.8e",i==0?"":", ",i,VB[i]);
924 fprintf(fp,"dA=%15.8e, dB=%15.8e\n",iparams->constr.dA,iparams->constr.dB);
927 fprintf(fp,"doh=%15.8e, dhh=%15.8e\n",iparams->settle.doh,
928 iparams->settle.dhh);
931 fprintf(fp,"a=%15.8e\n",iparams->vsite.a);
936 fprintf(fp,"a=%15.8e, b=%15.8e\n",iparams->vsite.a,iparams->vsite.b);
941 fprintf(fp,"a=%15.8e, b=%15.8e, c=%15.8e\n",
942 iparams->vsite.a,iparams->vsite.b,iparams->vsite.c);
945 fprintf(fp,"n=%2d, a=%15.8e\n",iparams->vsiten.n,iparams->vsiten.a);
950 fprintf(fp, "sar=%15.8e, st=%15.8e, pi=%15.8e, gbr=%15.8e, bmlt=%15.8e\n",iparams->gb.sar,iparams->gb.st,iparams->gb.pi,iparams->gb.gbr,iparams->gb.bmlt);
953 fprintf(fp, "cmapA=%1d, cmapB=%1d\n",iparams->cmap.cmapA, iparams->cmap.cmapB);
956 gmx_fatal(FARGS,"unknown function type %d (%s) in %s line %d",
957 ftype,interaction_function[ftype].name,__FILE__,__LINE__);
961 void pr_ilist(FILE *fp,int indent,const char *title,
962 t_functype *functype,t_ilist *ilist, gmx_bool bShowNumbers)
964 int i,j,k,type,ftype;
967 if (available(fp,ilist,indent,title) && ilist->nr > 0)
969 indent=pr_title(fp,indent,title);
970 (void) pr_indent(fp,indent);
971 fprintf(fp,"nr: %d\n",ilist->nr);
973 (void) pr_indent(fp,indent);
974 fprintf(fp,"iatoms:\n");
975 iatoms=ilist->iatoms;
976 for (i=j=0; i<ilist->nr;) {
978 (void) pr_indent(fp,indent+INDENT);
980 ftype=functype[type];
981 (void) fprintf(fp,"%d type=%d (%s)",
982 bShowNumbers?j:-1,bShowNumbers?type:-1,
983 interaction_function[ftype].name);
985 for (k=0; k<interaction_function[ftype].nratoms; k++)
986 (void) fprintf(fp," %u",*(iatoms++));
987 (void) fprintf(fp,"\n");
988 i+=1+interaction_function[ftype].nratoms;
990 fprintf(fp,"%5d%5d\n",i,iatoms[i]);
998 static void pr_cmap(FILE *fp, int indent, const char *title,
999 gmx_cmap_t *cmap_grid, gmx_bool bShowNumbers)
1004 dx = 360.0 / cmap_grid->grid_spacing;
1005 nelem = cmap_grid->grid_spacing*cmap_grid->grid_spacing;
1007 if(available(fp,cmap_grid,indent,title))
1009 fprintf(fp,"%s\n",title);
1011 for(i=0;i<cmap_grid->ngrid;i++)
1014 fprintf(fp,"%8s %8s %8s %8s\n","V","dVdx","dVdy","d2dV");
1016 fprintf(fp,"grid[%3d]={\n",bShowNumbers?i:-1);
1018 for(j=0;j<nelem;j++)
1020 if( (j%cmap_grid->grid_spacing)==0)
1022 fprintf(fp,"%8.1f\n",idx);
1026 fprintf(fp,"%8.3f ",cmap_grid->cmapdata[i].cmap[j*4]);
1027 fprintf(fp,"%8.3f ",cmap_grid->cmapdata[i].cmap[j*4+1]);
1028 fprintf(fp,"%8.3f ",cmap_grid->cmapdata[i].cmap[j*4+2]);
1029 fprintf(fp,"%8.3f\n",cmap_grid->cmapdata[i].cmap[j*4+3]);
1037 void pr_ffparams(FILE *fp,int indent,const char *title,
1038 gmx_ffparams_t *ffparams,
1039 gmx_bool bShowNumbers)
1043 indent=pr_title(fp,indent,title);
1044 (void) pr_indent(fp,indent);
1045 (void) fprintf(fp,"atnr=%d\n",ffparams->atnr);
1046 (void) pr_indent(fp,indent);
1047 (void) fprintf(fp,"ntypes=%d\n",ffparams->ntypes);
1048 for (i=0; i<ffparams->ntypes; i++) {
1049 (void) pr_indent(fp,indent+INDENT);
1050 (void) fprintf(fp,"functype[%d]=%s, ",
1052 interaction_function[ffparams->functype[i]].name);
1053 pr_iparams(fp,ffparams->functype[i],&ffparams->iparams[i]);
1055 (void) pr_double(fp,indent,"reppow",ffparams->reppow);
1056 (void) pr_real(fp,indent,"fudgeQQ",ffparams->fudgeQQ);
1057 pr_cmap(fp,indent,"cmap",&ffparams->cmap_grid,bShowNumbers);
1060 void pr_idef(FILE *fp,int indent,const char *title,t_idef *idef, gmx_bool bShowNumbers)
1064 if (available(fp,idef,indent,title)) {
1065 indent=pr_title(fp,indent,title);
1066 (void) pr_indent(fp,indent);
1067 (void) fprintf(fp,"atnr=%d\n",idef->atnr);
1068 (void) pr_indent(fp,indent);
1069 (void) fprintf(fp,"ntypes=%d\n",idef->ntypes);
1070 for (i=0; i<idef->ntypes; i++) {
1071 (void) pr_indent(fp,indent+INDENT);
1072 (void) fprintf(fp,"functype[%d]=%s, ",
1074 interaction_function[idef->functype[i]].name);
1075 pr_iparams(fp,idef->functype[i],&idef->iparams[i]);
1077 (void) pr_real(fp,indent,"fudgeQQ",idef->fudgeQQ);
1079 for(j=0; (j<F_NRE); j++)
1080 pr_ilist(fp,indent,interaction_function[j].longname,
1081 idef->functype,&idef->il[j],bShowNumbers);
1085 static int pr_block_title(FILE *fp,int indent,const char *title,t_block *block)
1089 if (available(fp,block,indent,title))
1091 indent=pr_title(fp,indent,title);
1092 (void) pr_indent(fp,indent);
1093 (void) fprintf(fp,"nr=%d\n",block->nr);
1098 static int pr_blocka_title(FILE *fp,int indent,const char *title,t_blocka *block)
1102 if (available(fp,block,indent,title))
1104 indent=pr_title(fp,indent,title);
1105 (void) pr_indent(fp,indent);
1106 (void) fprintf(fp,"nr=%d\n",block->nr);
1107 (void) pr_indent(fp,indent);
1108 (void) fprintf(fp,"nra=%d\n",block->nra);
1113 static void low_pr_block(FILE *fp,int indent,const char *title,t_block *block, gmx_bool bShowNumbers)
1117 if (available(fp,block,indent,title))
1119 indent=pr_block_title(fp,indent,title,block);
1120 for (i=0; i<=block->nr; i++)
1122 (void) pr_indent(fp,indent+INDENT);
1123 (void) fprintf(fp,"%s->index[%d]=%u\n",
1124 title,bShowNumbers?i:-1,block->index[i]);
1129 static void low_pr_blocka(FILE *fp,int indent,const char *title,t_blocka *block, gmx_bool bShowNumbers)
1133 if (available(fp,block,indent,title))
1135 indent=pr_blocka_title(fp,indent,title,block);
1136 for (i=0; i<=block->nr; i++)
1138 (void) pr_indent(fp,indent+INDENT);
1139 (void) fprintf(fp,"%s->index[%d]=%u\n",
1140 title,bShowNumbers?i:-1,block->index[i]);
1142 for (i=0; i<block->nra; i++)
1144 (void) pr_indent(fp,indent+INDENT);
1145 (void) fprintf(fp,"%s->a[%d]=%u\n",
1146 title,bShowNumbers?i:-1,block->a[i]);
1151 void pr_block(FILE *fp,int indent,const char *title,t_block *block,gmx_bool bShowNumbers)
1153 int i,j,ok,size,start,end;
1155 if (available(fp,block,indent,title))
1157 indent=pr_block_title(fp,indent,title,block);
1160 if ((ok=(block->index[start]==0))==0)
1161 (void) fprintf(fp,"block->index[%d] should be 0\n",start);
1163 for (i=0; i<block->nr; i++)
1165 end=block->index[i+1];
1166 size=pr_indent(fp,indent);
1168 size+=fprintf(fp,"%s[%d]={}\n",title,i);
1170 size+=fprintf(fp,"%s[%d]={%d..%d}\n",
1171 title,bShowNumbers?i:-1,
1172 bShowNumbers?start:-1,bShowNumbers?end-1:-1);
1178 void pr_blocka(FILE *fp,int indent,const char *title,t_blocka *block,gmx_bool bShowNumbers)
1180 int i,j,ok,size,start,end;
1182 if (available(fp,block,indent,title))
1184 indent=pr_blocka_title(fp,indent,title,block);
1187 if ((ok=(block->index[start]==0))==0)
1188 (void) fprintf(fp,"block->index[%d] should be 0\n",start);
1190 for (i=0; i<block->nr; i++)
1192 end=block->index[i+1];
1193 size=pr_indent(fp,indent);
1195 size+=fprintf(fp,"%s[%d]={",title,i);
1197 size+=fprintf(fp,"%s[%d][%d..%d]={",
1198 title,bShowNumbers?i:-1,
1199 bShowNumbers?start:-1,bShowNumbers?end-1:-1);
1200 for (j=start; j<end; j++)
1202 if (j>start) size+=fprintf(fp,", ");
1203 if ((size)>(USE_WIDTH))
1205 (void) fprintf(fp,"\n");
1206 size=pr_indent(fp,indent+INDENT);
1208 size+=fprintf(fp,"%u",block->a[j]);
1210 (void) fprintf(fp,"}\n");
1213 if ((end!=block->nra)||(!ok))
1215 (void) pr_indent(fp,indent);
1216 (void) fprintf(fp,"tables inconsistent, dumping complete tables:\n");
1217 low_pr_blocka(fp,indent,title,block,bShowNumbers);
1222 static void pr_strings(FILE *fp,int indent,const char *title,char ***nm,int n, gmx_bool bShowNumbers)
1226 if (available(fp,nm,indent,title))
1228 indent=pr_title_n(fp,indent,title,n);
1231 (void) pr_indent(fp,indent);
1232 (void) fprintf(fp,"%s[%d]={name=\"%s\"}\n",
1233 title,bShowNumbers?i:-1,*(nm[i]));
1238 static void pr_strings2(FILE *fp,int indent,const char *title,
1239 char ***nm,char ***nmB,int n, gmx_bool bShowNumbers)
1243 if (available(fp,nm,indent,title))
1245 indent=pr_title_n(fp,indent,title,n);
1248 (void) pr_indent(fp,indent);
1249 (void) fprintf(fp,"%s[%d]={name=\"%s\",nameB=\"%s\"}\n",
1250 title,bShowNumbers?i:-1,*(nm[i]),*(nmB[i]));
1255 static void pr_resinfo(FILE *fp,int indent,const char *title,t_resinfo *resinfo,int n, gmx_bool bShowNumbers)
1259 if (available(fp,resinfo,indent,title))
1261 indent=pr_title_n(fp,indent,title,n);
1264 (void) pr_indent(fp,indent);
1265 (void) fprintf(fp,"%s[%d]={name=\"%s\", nr=%d, ic='%c'}\n",
1266 title,bShowNumbers?i:-1,
1267 *(resinfo[i].name),resinfo[i].nr,
1268 (resinfo[i].ic == '\0') ? ' ' : resinfo[i].ic);
1273 static void pr_atom(FILE *fp,int indent,const char *title,t_atom *atom,int n)
1277 if (available(fp,atom,indent,title)) {
1278 indent=pr_title_n(fp,indent,title,n);
1279 for (i=0; i<n; i++) {
1280 (void) pr_indent(fp,indent);
1281 fprintf(fp,"%s[%6d]={type=%3d, typeB=%3d, ptype=%8s, m=%12.5e, "
1282 "q=%12.5e, mB=%12.5e, qB=%12.5e, resind=%5d, atomnumber=%3d}\n",
1283 title,i,atom[i].type,atom[i].typeB,ptype_str[atom[i].ptype],
1284 atom[i].m,atom[i].q,atom[i].mB,atom[i].qB,
1285 atom[i].resind,atom[i].atomnumber);
1290 static void pr_grps(FILE *fp,int indent,const char *title,t_grps grps[],
1291 char **grpname[], gmx_bool bShowNumbers)
1295 for(i=0; (i<egcNR); i++)
1297 fprintf(fp,"%s[%-12s] nr=%d, name=[",title,gtypes[i],grps[i].nr);
1298 for(j=0; (j<grps[i].nr); j++)
1300 fprintf(fp," %s",*(grpname[grps[i].nm_ind[j]]));
1306 static void pr_groups(FILE *fp,int indent,const char *title,
1307 gmx_groups_t *groups,
1308 gmx_bool bShowNumbers)
1313 pr_grps(fp,indent,"grp",groups->grps,groups->grpname,bShowNumbers);
1314 pr_strings(fp,indent,"grpname",groups->grpname,groups->ngrpname,bShowNumbers);
1316 (void) pr_indent(fp,indent);
1317 fprintf(fp,"groups ");
1318 for(g=0; g<egcNR; g++)
1320 printf(" %5.5s",gtypes[g]);
1324 (void) pr_indent(fp,indent);
1325 fprintf(fp,"allocated ");
1327 for(g=0; g<egcNR; g++)
1329 printf(" %5d",groups->ngrpnr[g]);
1330 nat_max = max(nat_max,groups->ngrpnr[g]);
1336 (void) pr_indent(fp,indent);
1337 fprintf(fp,"groupnr[%5s] =","*");
1338 for(g=0; g<egcNR; g++)
1340 fprintf(fp," %3d ",0);
1346 for(i=0; i<nat_max; i++)
1348 (void) pr_indent(fp,indent);
1349 fprintf(fp,"groupnr[%5d] =",i);
1350 for(g=0; g<egcNR; g++)
1353 groups->grpnr[g] ? groups->grpnr[g][i] : 0);
1360 void pr_atoms(FILE *fp,int indent,const char *title,t_atoms *atoms,
1361 gmx_bool bShownumbers)
1363 if (available(fp,atoms,indent,title))
1365 indent=pr_title(fp,indent,title);
1366 pr_atom(fp,indent,"atom",atoms->atom,atoms->nr);
1367 pr_strings(fp,indent,"atom",atoms->atomname,atoms->nr,bShownumbers);
1368 pr_strings2(fp,indent,"type",atoms->atomtype,atoms->atomtypeB,atoms->nr,bShownumbers);
1369 pr_resinfo(fp,indent,"residue",atoms->resinfo,atoms->nres,bShownumbers);
1374 void pr_atomtypes(FILE *fp,int indent,const char *title,t_atomtypes *atomtypes,
1375 gmx_bool bShowNumbers)
1378 if (available(fp,atomtypes,indent,title))
1380 indent=pr_title(fp,indent,title);
1381 for(i=0;i<atomtypes->nr;i++) {
1382 pr_indent(fp,indent);
1384 "atomtype[%3d]={radius=%12.5e, volume=%12.5e, gb_radius=%12.5e, surftens=%12.5e, atomnumber=%4d, S_hct=%12.5e)}\n",
1385 bShowNumbers?i:-1,atomtypes->radius[i],atomtypes->vol[i],
1386 atomtypes->gb_radius[i],
1387 atomtypes->surftens[i],atomtypes->atomnumber[i],atomtypes->S_hct[i]);
1392 static void pr_moltype(FILE *fp,int indent,const char *title,
1393 gmx_moltype_t *molt,int n,
1394 gmx_ffparams_t *ffparams,
1395 gmx_bool bShowNumbers)
1399 indent = pr_title_n(fp,indent,title,n);
1400 (void) pr_indent(fp,indent);
1401 (void) fprintf(fp,"name=\"%s\"\n",*(molt->name));
1402 pr_atoms(fp,indent,"atoms",&(molt->atoms),bShowNumbers);
1403 pr_block(fp,indent,"cgs",&molt->cgs, bShowNumbers);
1404 pr_blocka(fp,indent,"excls",&molt->excls, bShowNumbers);
1405 for(j=0; (j<F_NRE); j++) {
1406 pr_ilist(fp,indent,interaction_function[j].longname,
1407 ffparams->functype,&molt->ilist[j],bShowNumbers);
1411 static void pr_molblock(FILE *fp,int indent,const char *title,
1412 gmx_molblock_t *molb,int n,
1413 gmx_moltype_t *molt,
1414 gmx_bool bShowNumbers)
1416 indent = pr_title_n(fp,indent,title,n);
1417 (void) pr_indent(fp,indent);
1418 (void) fprintf(fp,"%-20s = %d \"%s\"\n",
1419 "moltype",molb->type,*(molt[molb->type].name));
1420 pr_int(fp,indent,"#molecules",molb->nmol);
1421 pr_int(fp,indent,"#atoms_mol",molb->natoms_mol);
1422 pr_int(fp,indent,"#posres_xA",molb->nposres_xA);
1423 if (molb->nposres_xA > 0) {
1424 pr_rvecs(fp,indent,"posres_xA",molb->posres_xA,molb->nposres_xA);
1426 pr_int(fp,indent,"#posres_xB",molb->nposres_xB);
1427 if (molb->nposres_xB > 0) {
1428 pr_rvecs(fp,indent,"posres_xB",molb->posres_xB,molb->nposres_xB);
1432 void pr_mtop(FILE *fp,int indent,const char *title,gmx_mtop_t *mtop,
1433 gmx_bool bShowNumbers)
1437 if (available(fp,mtop,indent,title)) {
1438 indent=pr_title(fp,indent,title);
1439 (void) pr_indent(fp,indent);
1440 (void) fprintf(fp,"name=\"%s\"\n",*(mtop->name));
1441 pr_int(fp,indent,"#atoms",mtop->natoms);
1442 for(mb=0; mb<mtop->nmolblock; mb++) {
1443 pr_molblock(fp,indent,"molblock",&mtop->molblock[mb],mb,
1444 mtop->moltype,bShowNumbers);
1446 pr_ffparams(fp,indent,"ffparams",&(mtop->ffparams),bShowNumbers);
1447 pr_atomtypes(fp,indent,"atomtypes",&(mtop->atomtypes),bShowNumbers);
1448 for(mt=0; mt<mtop->nmoltype; mt++) {
1449 pr_moltype(fp,indent,"moltype",&mtop->moltype[mt],mt,
1450 &mtop->ffparams,bShowNumbers);
1452 pr_groups(fp,indent,"groups",&mtop->groups,bShowNumbers);
1456 void pr_top(FILE *fp,int indent,const char *title,t_topology *top, gmx_bool bShowNumbers)
1458 if (available(fp,top,indent,title)) {
1459 indent=pr_title(fp,indent,title);
1460 (void) pr_indent(fp,indent);
1461 (void) fprintf(fp,"name=\"%s\"\n",*(top->name));
1462 pr_atoms(fp,indent,"atoms",&(top->atoms),bShowNumbers);
1463 pr_atomtypes(fp,indent,"atomtypes",&(top->atomtypes),bShowNumbers);
1464 pr_block(fp,indent,"cgs",&top->cgs, bShowNumbers);
1465 pr_block(fp,indent,"mols",&top->mols, bShowNumbers);
1466 pr_blocka(fp,indent,"excls",&top->excls, bShowNumbers);
1467 pr_idef(fp,indent,"idef",&top->idef,bShowNumbers);
1471 void pr_header(FILE *fp,int indent,const char *title,t_tpxheader *sh)
1475 if (available(fp,sh,indent,title))
1477 indent=pr_title(fp,indent,title);
1478 pr_indent(fp,indent);
1479 fprintf(fp,"bIr = %spresent\n",sh->bIr?"":"not ");
1480 pr_indent(fp,indent);
1481 fprintf(fp,"bBox = %spresent\n",sh->bBox?"":"not ");
1482 pr_indent(fp,indent);
1483 fprintf(fp,"bTop = %spresent\n",sh->bTop?"":"not ");
1484 pr_indent(fp,indent);
1485 fprintf(fp,"bX = %spresent\n",sh->bX?"":"not ");
1486 pr_indent(fp,indent);
1487 fprintf(fp,"bV = %spresent\n",sh->bV?"":"not ");
1488 pr_indent(fp,indent);
1489 fprintf(fp,"bF = %spresent\n",sh->bF?"":"not ");
1491 pr_indent(fp,indent);
1492 fprintf(fp,"natoms = %d\n",sh->natoms);
1493 pr_indent(fp,indent);
1494 fprintf(fp,"lambda = %e\n",sh->lambda);
1498 void pr_commrec(FILE *fp,int indent,t_commrec *cr)
1500 pr_indent(fp,indent);
1501 fprintf(fp,"commrec:\n");
1503 pr_indent(fp,indent);
1504 fprintf(fp,"nodeid = %d\n",cr->nodeid);
1505 pr_indent(fp,indent);
1506 fprintf(fp,"nnodes = %d\n",cr->nnodes);
1507 pr_indent(fp,indent);
1508 fprintf(fp,"npmenodes = %d\n",cr->npmenodes);
1510 pr_indent(fp,indent);
1511 fprintf(fp,"threadid = %d\n",cr->threadid);
1512 pr_indent(fp,indent);
1513 fprintf(fp,"nthreads = %d\n",cr->nthreads);