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
, 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
, 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
,bool vec
[],int n
, 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
, 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
, 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
, 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
)
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)
479 static void pr_pullgrp(FILE *fp
,int indent
,int g
,t_pullgrp
*pg
)
481 pr_indent(fp
,indent
);
482 fprintf(fp
,"pull_group %d:\n",g
);
484 pr_ivec_block(fp
,indent
,"atom",pg
->ind
,pg
->nat
,TRUE
);
485 pr_rvec(fp
,indent
,"weight",pg
->weight
,pg
->nweight
,TRUE
);
486 PI("pbcatom",pg
->pbcatom
);
487 pr_rvec(fp
,indent
,"vec",pg
->vec
,DIM
,TRUE
);
488 pr_rvec(fp
,indent
,"init",pg
->init
,DIM
,TRUE
);
494 static void pr_pull(FILE *fp
,int indent
,t_pull
*pull
)
498 PS("pull_geometry",EPULLGEOM(pull
->eGeom
));
499 pr_ivec(fp
,indent
,"pull_dim",pull
->dim
,DIM
,TRUE
);
500 PR("pull_r1",pull
->cyl_r1
);
501 PR("pull_r0",pull
->cyl_r0
);
502 PR("pull_constr_tol",pull
->constr_tol
);
503 PI("pull_nstxout",pull
->nstxout
);
504 PI("pull_nstfout",pull
->nstfout
);
505 PI("pull_ngrp",pull
->ngrp
);
506 for(g
=0; g
<pull
->ngrp
+1; g
++)
507 pr_pullgrp(fp
,indent
,g
,&pull
->grp
[g
]);
510 void pr_inputrec(FILE *fp
,int indent
,const char *title
,t_inputrec
*ir
,
513 const char *infbuf
="inf";
516 if (available(fp
,ir
,indent
,title
)) {
518 indent
=pr_title(fp
,indent
,title
);
519 PS("integrator",EI(ir
->eI
));
520 PSTEP("nsteps",ir
->nsteps
);
521 PSTEP("init_step",ir
->init_step
);
522 PI("nstcalcenergy",ir
->nstcalcenergy
);
523 PS("ns_type",ENS(ir
->ns_type
));
524 PI("nstlist",ir
->nstlist
);
525 PI("ndelta",ir
->ndelta
);
526 PI("nstcomm",ir
->nstcomm
);
527 PS("comm_mode",ECOM(ir
->comm_mode
));
528 PI("nstlog",ir
->nstlog
);
529 PI("nstxout",ir
->nstxout
);
530 PI("nstvout",ir
->nstvout
);
531 PI("nstfout",ir
->nstfout
);
532 PI("nstenergy",ir
->nstenergy
);
533 PI("nstxtcout",ir
->nstxtcout
);
534 PR("init_t",ir
->init_t
);
535 PR("delta_t",ir
->delta_t
);
537 PR("xtcprec",ir
->xtcprec
);
541 PI("pme_order",ir
->pme_order
);
542 PR("ewald_rtol",ir
->ewald_rtol
);
543 PR("ewald_geometry",ir
->ewald_geometry
);
544 PR("epsilon_surface",ir
->epsilon_surface
);
545 PS("optimize_fft",BOOL(ir
->bOptFFT
));
546 PS("ePBC",EPBC(ir
->ePBC
));
547 PS("bPeriodicMols",BOOL(ir
->bPeriodicMols
));
548 PS("bContinuation",BOOL(ir
->bContinuation
));
549 PS("bShakeSOR",BOOL(ir
->bShakeSOR
));
550 PS("etc",ETCOUPLTYPE(ir
->etc
));
551 PS("epc",EPCOUPLTYPE(ir
->epc
));
552 PS("epctype",EPCOUPLTYPETYPE(ir
->epct
));
553 PR("tau_p",ir
->tau_p
);
554 pr_matrix(fp
,indent
,"ref_p",ir
->ref_p
,bMDPformat
);
555 pr_matrix(fp
,indent
,"compress",ir
->compress
,bMDPformat
);
556 PS("refcoord_scaling",EREFSCALINGTYPE(ir
->refcoord_scaling
));
558 fprintf(fp
,"posres_com = %g %g %g\n",ir
->posres_com
[XX
],
559 ir
->posres_com
[YY
],ir
->posres_com
[ZZ
]);
561 pr_rvec(fp
,indent
,"posres_com",ir
->posres_com
,DIM
,TRUE
);
563 fprintf(fp
,"posres_comB = %g %g %g\n",ir
->posres_comB
[XX
],
564 ir
->posres_comB
[YY
],ir
->posres_comB
[ZZ
]);
566 pr_rvec(fp
,indent
,"posres_comB",ir
->posres_comB
,DIM
,TRUE
);
567 PI("andersen_seed",ir
->andersen_seed
);
568 PR("rlist",ir
->rlist
);
569 PR("rlistlong",ir
->rlistlong
);
571 PS("coulombtype",EELTYPE(ir
->coulombtype
));
572 PR("rcoulomb_switch",ir
->rcoulomb_switch
);
573 PR("rcoulomb",ir
->rcoulomb
);
574 PS("vdwtype",EVDWTYPE(ir
->vdwtype
));
575 PR("rvdw_switch",ir
->rvdw_switch
);
577 if (ir
->epsilon_r
!= 0)
578 PR("epsilon_r",ir
->epsilon_r
);
580 PS("epsilon_r",infbuf
);
581 if (ir
->epsilon_rf
!= 0)
582 PR("epsilon_rf",ir
->epsilon_rf
);
584 PS("epsilon_rf",infbuf
);
585 PR("tabext",ir
->tabext
);
586 PS("implicit_solvent",EIMPLICITSOL(ir
->implicit_solvent
));
587 PS("gb_algorithm",EGBALGORITHM(ir
->gb_algorithm
));
588 PR("gb_epsilon_solvent",ir
->gb_epsilon_solvent
);
589 PI("nstgbradii",ir
->nstgbradii
);
590 PR("rgbradii",ir
->rgbradii
);
591 PR("gb_saltconc",ir
->gb_saltconc
);
592 PR("gb_obc_alpha",ir
->gb_obc_alpha
);
593 PR("gb_obc_beta",ir
->gb_obc_beta
);
594 PR("gb_obc_gamma",ir
->gb_obc_gamma
);
595 PR("gb_dielectric_offset",ir
->gb_dielectric_offset
);
596 PS("sa_algorithm",ESAALGORITHM(ir
->gb_algorithm
));
597 PR("sa_surface_tension",ir
->sa_surface_tension
);
599 PS("DispCorr",EDISPCORR(ir
->eDispCorr
));
600 PS("free_energy",EFEPTYPE(ir
->efep
));
601 PR("init_lambda",ir
->init_lambda
);
602 PR("delta_lambda",ir
->delta_lambda
);
605 PI("n_foreign_lambda",ir
->n_flambda
);
607 if (ir
->n_flambda
> 0)
609 pr_indent(fp
,indent
);
610 fprintf(fp
,"foreign_lambda%s",bMDPformat
? " = " : ":");
611 for(i
=0; i
<ir
->n_flambda
; i
++)
613 fprintf(fp
," %10g",ir
->flambda
[i
]);
617 PR("sc_alpha",ir
->sc_alpha
);
618 PI("sc_power",ir
->sc_power
);
619 PR("sc_sigma",ir
->sc_sigma
);
620 PI("nstdhdl", ir
->nstdhdl
);
622 PI("nwall",ir
->nwall
);
623 PS("wall_type",EWALLTYPE(ir
->wall_type
));
624 PI("wall_atomtype[0]",ir
->wall_atomtype
[0]);
625 PI("wall_atomtype[1]",ir
->wall_atomtype
[1]);
626 PR("wall_density[0]",ir
->wall_density
[0]);
627 PR("wall_density[1]",ir
->wall_density
[1]);
628 PR("wall_ewald_zfac",ir
->wall_ewald_zfac
);
630 PS("pull",EPULLTYPE(ir
->ePull
));
631 if (ir
->ePull
!= epullNO
)
632 pr_pull(fp
,indent
,ir
->pull
);
634 PS("disre",EDISRETYPE(ir
->eDisre
));
635 PS("disre_weighting",EDISREWEIGHTING(ir
->eDisreWeighting
));
636 PS("disre_mixed",BOOL(ir
->bDisreMixed
));
637 PR("dr_fc",ir
->dr_fc
);
638 PR("dr_tau",ir
->dr_tau
);
639 PR("nstdisreout",ir
->nstdisreout
);
640 PR("orires_fc",ir
->orires_fc
);
641 PR("orires_tau",ir
->orires_tau
);
642 PR("nstorireout",ir
->nstorireout
);
644 PR("dihre-fc",ir
->dihre_fc
);
646 PR("em_stepsize",ir
->em_stepsize
);
647 PR("em_tol",ir
->em_tol
);
648 PI("niter",ir
->niter
);
649 PR("fc_stepsize",ir
->fc_stepsize
);
650 PI("nstcgsteep",ir
->nstcgsteep
);
651 PI("nbfgscorr",ir
->nbfgscorr
);
653 PS("ConstAlg",ECONSTRTYPE(ir
->eConstrAlg
));
654 PR("shake_tol",ir
->shake_tol
);
655 PI("lincs_order",ir
->nProjOrder
);
656 PR("lincs_warnangle",ir
->LincsWarnAngle
);
657 PI("lincs_iter",ir
->nLincsIter
);
658 PR("bd_fric",ir
->bd_fric
);
659 PI("ld_seed",ir
->ld_seed
);
660 PR("cos_accel",ir
->cos_accel
);
661 pr_matrix(fp
,indent
,"deform",ir
->deform
,bMDPformat
);
663 PS("adress_type",EADRESSTYPE(ir
->adress_type
));
664 PS("adress_new_wf",BOOL(ir
->badress_new_wf
));
665 PR("adress_const_wf",ir
->adress_const_wf
);
666 PR("adress_ex_width",ir
->adress_ex_width
);
667 PR("adress_hy_width",ir
->adress_hy_width
);
668 PS("adress_interface_correction",EADRESSICTYPE(ir
->adress_icor
));
669 PS("adress_exvdw",EVDWTYPE(ir
->adress_ivdw
));
670 PS("adress_site",EADRESSSITETYPE(ir
->adress_site
));
671 pr_rvecs(fp
,indent
,"adress_reference_coords",&(ir
->adress_refs
),bMDPformat
);
673 PI("userint1",ir
->userint1
);
674 PI("userint2",ir
->userint2
);
675 PI("userint3",ir
->userint3
);
676 PI("userint4",ir
->userint4
);
677 PR("userreal1",ir
->userreal1
);
678 PR("userreal2",ir
->userreal2
);
679 PR("userreal3",ir
->userreal3
);
680 PR("userreal4",ir
->userreal4
);
681 pr_grp_opts(fp
,indent
,"grpopts",&(ir
->opts
),bMDPformat
);
682 pr_cosine(fp
,indent
,"efield-x",&(ir
->ex
[XX
]),bMDPformat
);
683 pr_cosine(fp
,indent
,"efield-xt",&(ir
->et
[XX
]),bMDPformat
);
684 pr_cosine(fp
,indent
,"efield-y",&(ir
->ex
[YY
]),bMDPformat
);
685 pr_cosine(fp
,indent
,"efield-yt",&(ir
->et
[YY
]),bMDPformat
);
686 pr_cosine(fp
,indent
,"efield-z",&(ir
->ex
[ZZ
]),bMDPformat
);
687 pr_cosine(fp
,indent
,"efield-zt",&(ir
->et
[ZZ
]),bMDPformat
);
688 PS("bQMMM",BOOL(ir
->bQMMM
));
689 PI("QMconstraints",ir
->QMconstraints
);
690 PI("QMMMscheme",ir
->QMMMscheme
);
691 PR("scalefactor",ir
->scalefactor
);
692 pr_qm_opts(fp
,indent
,"qm_opts",&(ir
->opts
));
699 static void pr_harm(FILE *fp
,t_iparams
*iparams
,const char *r
,const char *kr
)
701 fprintf(fp
,"%sA=%12.5e, %sA=%12.5e, %sB=%12.5e, %sB=%12.5e\n",
702 r
,iparams
->harmonic
.rA
,kr
,iparams
->harmonic
.krA
,
703 r
,iparams
->harmonic
.rB
,kr
,iparams
->harmonic
.krB
);
706 void pr_iparams(FILE *fp
,t_functype ftype
,t_iparams
*iparams
)
709 real VA
[4],VB
[4],*rbcA
,*rbcB
;
714 pr_harm(fp
,iparams
,"th","ct");
716 case F_CROSS_BOND_BONDS
:
717 fprintf(fp
,"r1e=%15.8e, r2e=%15.8e, krr=%15.8e\n",
718 iparams
->cross_bb
.r1e
,iparams
->cross_bb
.r2e
,
719 iparams
->cross_bb
.krr
);
721 case F_CROSS_BOND_ANGLES
:
722 fprintf(fp
,"r1e=%15.8e, r1e=%15.8e, r3e=%15.8e, krt=%15.8e\n",
723 iparams
->cross_ba
.r1e
,iparams
->cross_ba
.r2e
,
724 iparams
->cross_ba
.r3e
,iparams
->cross_ba
.krt
);
727 fprintf(fp
,"theta=%15.8e, ktheta=%15.8e, r13=%15.8e, kUB=%15.8e\n",
728 iparams
->u_b
.theta
,iparams
->u_b
.ktheta
,iparams
->u_b
.r13
,iparams
->u_b
.kUB
);
730 case F_QUARTIC_ANGLES
:
731 fprintf(fp
,"theta=%15.8e",iparams
->qangle
.theta
);
733 fprintf(fp
,", c%c=%15.8e",'0'+i
,iparams
->qangle
.c
[i
]);
737 fprintf(fp
,"a=%15.8e, b=%15.8e, c=%15.8e\n",
738 iparams
->bham
.a
,iparams
->bham
.b
,iparams
->bham
.c
);
743 pr_harm(fp
,iparams
,"b0","cb");
746 pr_harm(fp
,iparams
,"xi","cx");
749 fprintf(fp
,"b0=%15.8e, cb=%15.8e, beta=%15.8e\n",
750 iparams
->morse
.b0
,iparams
->morse
.cb
,iparams
->morse
.beta
);
753 fprintf(fp
,"b0=%15.8e, kb=%15.8e, kcub=%15.8e\n",
754 iparams
->cubic
.b0
,iparams
->cubic
.kb
,iparams
->cubic
.kcub
);
760 fprintf(fp
,"bm=%15.8e, kb=%15.8e\n",iparams
->fene
.bm
,iparams
->fene
.kb
);
766 fprintf(fp
,"tab=%d, kA=%15.8e, kB=%15.8e\n",
767 iparams
->tab
.table
,iparams
->tab
.kA
,iparams
->tab
.kB
);
770 fprintf(fp
,"alpha=%15.8e\n",iparams
->polarize
.alpha
);
773 fprintf(fp
,"a=%15.8e, alpha1=%15.8e, alpha2=%15.8e, rfac=%15.8e\n",
774 iparams
->thole
.a
,iparams
->thole
.alpha1
,iparams
->thole
.alpha2
,
775 iparams
->thole
.rfac
);
778 fprintf(fp
,"al_x=%15.8e, al_y=%15.8e, al_z=%15.8e, rOH=%9.6f, rHH=%9.6f, rOD=%9.6f\n",
779 iparams
->wpol
.al_x
,iparams
->wpol
.al_y
,iparams
->wpol
.al_z
,
780 iparams
->wpol
.rOH
,iparams
->wpol
.rHH
,iparams
->wpol
.rOD
);
783 fprintf(fp
,"c6=%15.8e, c12=%15.8e\n",iparams
->lj
.c6
,iparams
->lj
.c12
);
786 fprintf(fp
,"c6A=%15.8e, c12A=%15.8e, c6B=%15.8e, c12B=%15.8e\n",
787 iparams
->lj14
.c6A
,iparams
->lj14
.c12A
,
788 iparams
->lj14
.c6B
,iparams
->lj14
.c12B
);
791 fprintf(fp
,"fqq=%15.8e, qi=%15.8e, qj=%15.8e, c6=%15.8e, c12=%15.8e\n",
793 iparams
->ljc14
.qi
,iparams
->ljc14
.qj
,
794 iparams
->ljc14
.c6
,iparams
->ljc14
.c12
);
797 fprintf(fp
,"qi=%15.8e, qj=%15.8e, c6=%15.8e, c12=%15.8e\n",
798 iparams
->ljcnb
.qi
,iparams
->ljcnb
.qj
,
799 iparams
->ljcnb
.c6
,iparams
->ljcnb
.c12
);
804 fprintf(fp
,"phiA=%15.8e, cpA=%15.8e, phiB=%15.8e, cpB=%15.8e, mult=%d\n",
805 iparams
->pdihs
.phiA
,iparams
->pdihs
.cpA
,
806 iparams
->pdihs
.phiB
,iparams
->pdihs
.cpB
,
807 iparams
->pdihs
.mult
);
810 fprintf(fp
,"label=%4d, type=%1d, low=%15.8e, up1=%15.8e, up2=%15.8e, fac=%15.8e)\n",
811 iparams
->disres
.label
,iparams
->disres
.type
,
812 iparams
->disres
.low
,iparams
->disres
.up1
,
813 iparams
->disres
.up2
,iparams
->disres
.kfac
);
816 fprintf(fp
,"ex=%4d, label=%d, power=%4d, c=%15.8e, obs=%15.8e, kfac=%15.8e)\n",
817 iparams
->orires
.ex
,iparams
->orires
.label
,iparams
->orires
.power
,
818 iparams
->orires
.c
,iparams
->orires
.obs
,iparams
->orires
.kfac
);
821 fprintf(fp
,"label=%d, power=%4d phi=%15.8e, dphi=%15.8e, kfac=%15.8e)\n",
822 iparams
->dihres
.label
,iparams
->dihres
.power
,
823 iparams
->dihres
.phi
,iparams
->dihres
.dphi
,iparams
->dihres
.kfac
);
826 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",
827 iparams
->posres
.pos0A
[XX
],iparams
->posres
.pos0A
[YY
],
828 iparams
->posres
.pos0A
[ZZ
],iparams
->posres
.fcA
[XX
],
829 iparams
->posres
.fcA
[YY
],iparams
->posres
.fcA
[ZZ
],
830 iparams
->posres
.pos0B
[XX
],iparams
->posres
.pos0B
[YY
],
831 iparams
->posres
.pos0B
[ZZ
],iparams
->posres
.fcB
[XX
],
832 iparams
->posres
.fcB
[YY
],iparams
->posres
.fcB
[ZZ
]);
835 for (i
=0; i
<NR_RBDIHS
; i
++)
836 fprintf(fp
,"%srbcA[%d]=%15.8e",i
==0?"":", ",i
,iparams
->rbdihs
.rbcA
[i
]);
838 for (i
=0; i
<NR_RBDIHS
; i
++)
839 fprintf(fp
,"%srbcB[%d]=%15.8e",i
==0?"":", ",i
,iparams
->rbdihs
.rbcB
[i
]);
843 /* Use the OPLS -> Ryckaert-Bellemans formula backwards to get the
844 * OPLS potential constants back.
846 rbcA
= iparams
->rbdihs
.rbcA
;
847 rbcB
= iparams
->rbdihs
.rbcB
;
849 VA
[3] = -0.25*rbcA
[4];
850 VA
[2] = -0.5*rbcA
[3];
851 VA
[1] = 4.0*VA
[3]-rbcA
[2];
852 VA
[0] = 3.0*VA
[2]-2.0*rbcA
[1];
854 VB
[3] = -0.25*rbcB
[4];
855 VB
[2] = -0.5*rbcB
[3];
856 VB
[1] = 4.0*VB
[3]-rbcB
[2];
857 VB
[0] = 3.0*VB
[2]-2.0*rbcB
[1];
859 for (i
=0; i
<NR_FOURDIHS
; i
++)
860 fprintf(fp
,"%sFourA[%d]=%15.8e",i
==0?"":", ",i
,VA
[i
]);
862 for (i
=0; i
<NR_FOURDIHS
; i
++)
863 fprintf(fp
,"%sFourB[%d]=%15.8e",i
==0?"":", ",i
,VB
[i
]);
869 fprintf(fp
,"dA=%15.8e, dB=%15.8e\n",iparams
->constr
.dA
,iparams
->constr
.dB
);
872 fprintf(fp
,"doh=%15.8e, dhh=%15.8e\n",iparams
->settle
.doh
,
873 iparams
->settle
.dhh
);
876 fprintf(fp
,"a=%15.8e\n",iparams
->vsite
.a
);
881 fprintf(fp
,"a=%15.8e, b=%15.8e\n",iparams
->vsite
.a
,iparams
->vsite
.b
);
886 fprintf(fp
,"a=%15.8e, b=%15.8e, c=%15.8e\n",
887 iparams
->vsite
.a
,iparams
->vsite
.b
,iparams
->vsite
.c
);
890 fprintf(fp
,"n=%2d, a=%15.8e\n",iparams
->vsiten
.n
,iparams
->vsiten
.a
);
895 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
);
898 fprintf(fp
, "cmapA=%1d, cmapB=%1d\n",iparams
->cmap
.cmapA
, iparams
->cmap
.cmapB
);
901 gmx_fatal(FARGS
,"unknown function type %d (%s) in %s line %d",
902 ftype
,interaction_function
[ftype
].name
,__FILE__
,__LINE__
);
906 void pr_ilist(FILE *fp
,int indent
,const char *title
,
907 t_functype
*functype
,t_ilist
*ilist
, bool bShowNumbers
)
909 int i
,j
,k
,type
,ftype
;
912 if (available(fp
,ilist
,indent
,title
) && ilist
->nr
> 0)
914 indent
=pr_title(fp
,indent
,title
);
915 (void) pr_indent(fp
,indent
);
916 fprintf(fp
,"nr: %d\n",ilist
->nr
);
918 (void) pr_indent(fp
,indent
);
919 fprintf(fp
,"iatoms:\n");
920 iatoms
=ilist
->iatoms
;
921 for (i
=j
=0; i
<ilist
->nr
;) {
923 (void) pr_indent(fp
,indent
+INDENT
);
925 ftype
=functype
[type
];
926 (void) fprintf(fp
,"%d type=%d (%s)",
927 bShowNumbers
?j
:-1,bShowNumbers
?type
:-1,
928 interaction_function
[ftype
].name
);
930 for (k
=0; k
<interaction_function
[ftype
].nratoms
; k
++)
931 (void) fprintf(fp
," %u",*(iatoms
++));
932 (void) fprintf(fp
,"\n");
933 i
+=1+interaction_function
[ftype
].nratoms
;
935 fprintf(fp
,"%5d%5d\n",i
,iatoms
[i
]);
943 void pr_ffparams(FILE *fp
,int indent
,const char *title
,
944 gmx_ffparams_t
*ffparams
,
949 indent
=pr_title(fp
,indent
,title
);
950 (void) pr_indent(fp
,indent
);
951 (void) fprintf(fp
,"atnr=%d\n",ffparams
->atnr
);
952 (void) pr_indent(fp
,indent
);
953 (void) fprintf(fp
,"ntypes=%d\n",ffparams
->ntypes
);
954 for (i
=0; i
<ffparams
->ntypes
; i
++) {
955 (void) pr_indent(fp
,indent
+INDENT
);
956 (void) fprintf(fp
,"functype[%d]=%s, ",
958 interaction_function
[ffparams
->functype
[i
]].name
);
959 pr_iparams(fp
,ffparams
->functype
[i
],&ffparams
->iparams
[i
]);
961 (void) pr_double(fp
,indent
,"reppow",ffparams
->reppow
);
962 (void) pr_real(fp
,indent
,"fudgeQQ",ffparams
->fudgeQQ
);
965 void pr_idef(FILE *fp
,int indent
,const char *title
,t_idef
*idef
, bool bShowNumbers
)
969 if (available(fp
,idef
,indent
,title
)) {
970 indent
=pr_title(fp
,indent
,title
);
971 (void) pr_indent(fp
,indent
);
972 (void) fprintf(fp
,"atnr=%d\n",idef
->atnr
);
973 (void) pr_indent(fp
,indent
);
974 (void) fprintf(fp
,"ntypes=%d\n",idef
->ntypes
);
975 for (i
=0; i
<idef
->ntypes
; i
++) {
976 (void) pr_indent(fp
,indent
+INDENT
);
977 (void) fprintf(fp
,"functype[%d]=%s, ",
979 interaction_function
[idef
->functype
[i
]].name
);
980 pr_iparams(fp
,idef
->functype
[i
],&idef
->iparams
[i
]);
982 (void) pr_real(fp
,indent
,"fudgeQQ",idef
->fudgeQQ
);
984 for(j
=0; (j
<F_NRE
); j
++)
985 pr_ilist(fp
,indent
,interaction_function
[j
].longname
,
986 idef
->functype
,&idef
->il
[j
],bShowNumbers
);
990 static int pr_block_title(FILE *fp
,int indent
,const char *title
,t_block
*block
)
994 if (available(fp
,block
,indent
,title
))
996 indent
=pr_title(fp
,indent
,title
);
997 (void) pr_indent(fp
,indent
);
998 (void) fprintf(fp
,"nr=%d\n",block
->nr
);
1003 static int pr_blocka_title(FILE *fp
,int indent
,const char *title
,t_blocka
*block
)
1007 if (available(fp
,block
,indent
,title
))
1009 indent
=pr_title(fp
,indent
,title
);
1010 (void) pr_indent(fp
,indent
);
1011 (void) fprintf(fp
,"nr=%d\n",block
->nr
);
1012 (void) pr_indent(fp
,indent
);
1013 (void) fprintf(fp
,"nra=%d\n",block
->nra
);
1018 static void low_pr_block(FILE *fp
,int indent
,const char *title
,t_block
*block
, bool bShowNumbers
)
1022 if (available(fp
,block
,indent
,title
))
1024 indent
=pr_block_title(fp
,indent
,title
,block
);
1025 for (i
=0; i
<=block
->nr
; i
++)
1027 (void) pr_indent(fp
,indent
+INDENT
);
1028 (void) fprintf(fp
,"%s->index[%d]=%u\n",
1029 title
,bShowNumbers
?i
:-1,block
->index
[i
]);
1034 static void low_pr_blocka(FILE *fp
,int indent
,const char *title
,t_blocka
*block
, bool bShowNumbers
)
1038 if (available(fp
,block
,indent
,title
))
1040 indent
=pr_blocka_title(fp
,indent
,title
,block
);
1041 for (i
=0; i
<=block
->nr
; i
++)
1043 (void) pr_indent(fp
,indent
+INDENT
);
1044 (void) fprintf(fp
,"%s->index[%d]=%u\n",
1045 title
,bShowNumbers
?i
:-1,block
->index
[i
]);
1047 for (i
=0; i
<block
->nra
; i
++)
1049 (void) pr_indent(fp
,indent
+INDENT
);
1050 (void) fprintf(fp
,"%s->a[%d]=%u\n",
1051 title
,bShowNumbers
?i
:-1,block
->a
[i
]);
1056 void pr_block(FILE *fp
,int indent
,const char *title
,t_block
*block
,bool bShowNumbers
)
1058 int i
,j
,ok
,size
,start
,end
;
1060 if (available(fp
,block
,indent
,title
))
1062 indent
=pr_block_title(fp
,indent
,title
,block
);
1065 if ((ok
=(block
->index
[start
]==0))==0)
1066 (void) fprintf(fp
,"block->index[%d] should be 0\n",start
);
1068 for (i
=0; i
<block
->nr
; i
++)
1070 end
=block
->index
[i
+1];
1071 size
=pr_indent(fp
,indent
);
1073 size
+=fprintf(fp
,"%s[%d]={}\n",title
,i
);
1075 size
+=fprintf(fp
,"%s[%d]={%d..%d}\n",
1076 title
,bShowNumbers
?i
:-1,
1077 bShowNumbers
?start
:-1,bShowNumbers
?end
-1:-1);
1083 void pr_blocka(FILE *fp
,int indent
,const char *title
,t_blocka
*block
,bool bShowNumbers
)
1085 int i
,j
,ok
,size
,start
,end
;
1087 if (available(fp
,block
,indent
,title
))
1089 indent
=pr_blocka_title(fp
,indent
,title
,block
);
1092 if ((ok
=(block
->index
[start
]==0))==0)
1093 (void) fprintf(fp
,"block->index[%d] should be 0\n",start
);
1095 for (i
=0; i
<block
->nr
; i
++)
1097 end
=block
->index
[i
+1];
1098 size
=pr_indent(fp
,indent
);
1100 size
+=fprintf(fp
,"%s[%d]={",title
,i
);
1102 size
+=fprintf(fp
,"%s[%d][%d..%d]={",
1103 title
,bShowNumbers
?i
:-1,
1104 bShowNumbers
?start
:-1,bShowNumbers
?end
-1:-1);
1105 for (j
=start
; j
<end
; j
++)
1107 if (j
>start
) size
+=fprintf(fp
,", ");
1108 if ((size
)>(USE_WIDTH
))
1110 (void) fprintf(fp
,"\n");
1111 size
=pr_indent(fp
,indent
+INDENT
);
1113 size
+=fprintf(fp
,"%u",block
->a
[j
]);
1115 (void) fprintf(fp
,"}\n");
1118 if ((end
!=block
->nra
)||(!ok
))
1120 (void) pr_indent(fp
,indent
);
1121 (void) fprintf(fp
,"tables inconsistent, dumping complete tables:\n");
1122 low_pr_blocka(fp
,indent
,title
,block
,bShowNumbers
);
1127 static void pr_strings(FILE *fp
,int indent
,const char *title
,char ***nm
,int n
, bool bShowNumbers
)
1131 if (available(fp
,nm
,indent
,title
))
1133 indent
=pr_title_n(fp
,indent
,title
,n
);
1136 (void) pr_indent(fp
,indent
);
1137 (void) fprintf(fp
,"%s[%d]={name=\"%s\"}\n",
1138 title
,bShowNumbers
?i
:-1,*(nm
[i
]));
1143 static void pr_strings2(FILE *fp
,int indent
,const char *title
,
1144 char ***nm
,char ***nmB
,int n
, bool bShowNumbers
)
1148 if (available(fp
,nm
,indent
,title
))
1150 indent
=pr_title_n(fp
,indent
,title
,n
);
1153 (void) pr_indent(fp
,indent
);
1154 (void) fprintf(fp
,"%s[%d]={name=\"%s\",nameB=\"%s\"}\n",
1155 title
,bShowNumbers
?i
:-1,*(nm
[i
]),*(nmB
[i
]));
1160 static void pr_resinfo(FILE *fp
,int indent
,const char *title
,t_resinfo
*resinfo
,int n
, bool bShowNumbers
)
1164 if (available(fp
,resinfo
,indent
,title
))
1166 indent
=pr_title_n(fp
,indent
,title
,n
);
1169 (void) pr_indent(fp
,indent
);
1170 (void) fprintf(fp
,"%s[%d]={name=\"%s\", nr=%d, ic='%c'}\n",
1171 title
,bShowNumbers
?i
:-1,
1172 *(resinfo
[i
].name
),resinfo
[i
].nr
,
1173 (resinfo
[i
].ic
== '\0') ? ' ' : resinfo
[i
].ic
);
1178 static void pr_atom(FILE *fp
,int indent
,const char *title
,t_atom
*atom
,int n
)
1182 if (available(fp
,atom
,indent
,title
)) {
1183 indent
=pr_title_n(fp
,indent
,title
,n
);
1184 for (i
=0; i
<n
; i
++) {
1185 (void) pr_indent(fp
,indent
);
1186 fprintf(fp
,"%s[%6d]={type=%3d, typeB=%3d, ptype=%8s, m=%12.5e, "
1187 "q=%12.5e, mB=%12.5e, qB=%12.5e, resind=%5d, atomnumber=%3d}\n",
1188 title
,i
,atom
[i
].type
,atom
[i
].typeB
,ptype_str
[atom
[i
].ptype
],
1189 atom
[i
].m
,atom
[i
].q
,atom
[i
].mB
,atom
[i
].qB
,
1190 atom
[i
].resind
,atom
[i
].atomnumber
);
1195 static void pr_grps(FILE *fp
,int indent
,const char *title
,t_grps grps
[],int ngrp
,
1196 char **grpname
[], bool bShowNumbers
)
1200 for(i
=0; (i
<ngrp
); i
++) {
1201 fprintf(fp
,"%s[%d] nr=%d, name=[",title
,bShowNumbers
?i
:-1,grps
[i
].nr
);
1202 for(j
=0; (j
<grps
[i
].nr
); j
++)
1203 fprintf(fp
," %s",*(grpname
[grps
[i
].nm_ind
[j
]]));
1208 static void pr_groups(FILE *fp
,int indent
,const char *title
,
1209 gmx_groups_t
*groups
,
1215 pr_grps(fp
,indent
,"grp",groups
->grps
,egcNR
,groups
->grpname
,bShowNumbers
);
1216 pr_strings(fp
,indent
,"grpname",groups
->grpname
,groups
->ngrpname
,bShowNumbers
);
1218 (void) pr_indent(fp
,indent
);
1219 fprintf(fp
,"groups ");
1220 for(g
=0; g
<egcNR
; g
++)
1222 printf(" %5.5s",gtypes
[g
]);
1226 (void) pr_indent(fp
,indent
);
1227 fprintf(fp
,"allocated ");
1229 for(g
=0; g
<egcNR
; g
++)
1231 printf(" %5d",groups
->ngrpnr
[g
]);
1232 nat_max
= max(nat_max
,groups
->ngrpnr
[g
]);
1238 (void) pr_indent(fp
,indent
);
1239 fprintf(fp
,"groupnr[%5s] =","*");
1240 for(g
=0; g
<egcNR
; g
++)
1242 fprintf(fp
," %3d ",0);
1248 for(i
=0; i
<nat_max
; i
++)
1250 (void) pr_indent(fp
,indent
);
1251 fprintf(fp
,"groupnr[%5d] =",i
);
1252 for(g
=0; g
<egcNR
; g
++)
1255 groups
->grpnr
[g
] ? groups
->grpnr
[g
][i
] : 0);
1262 void pr_atoms(FILE *fp
,int indent
,const char *title
,t_atoms
*atoms
,
1265 if (available(fp
,atoms
,indent
,title
))
1267 indent
=pr_title(fp
,indent
,title
);
1268 pr_atom(fp
,indent
,"atom",atoms
->atom
,atoms
->nr
);
1269 pr_strings(fp
,indent
,"atom",atoms
->atomname
,atoms
->nr
,bShownumbers
);
1270 pr_strings2(fp
,indent
,"type",atoms
->atomtype
,atoms
->atomtypeB
,atoms
->nr
,bShownumbers
);
1271 pr_resinfo(fp
,indent
,"residue",atoms
->resinfo
,atoms
->nres
,bShownumbers
);
1275 void pr_cmap(FILE *fp
, int indent
, const char *title
, gmx_cmap_t
*cmap_grid
,
1281 dx
= 360.0 / cmap_grid
->grid_spacing
;
1282 nelem
= cmap_grid
->grid_spacing
*cmap_grid
->grid_spacing
;
1284 if(available(fp
,cmap_grid
,indent
,title
))
1286 fprintf(fp
,"%s\n",title
);
1288 for(i
=0;i
<cmap_grid
->ngrid
;i
++)
1291 fprintf(fp
,"%8s %8s %8s %8s\n","V","dVdx","dVdy","d2dV");
1293 fprintf(fp
,"grid[%3d]={\n",bShowNumbers
?i
:-1);
1295 for(j
=0;j
<nelem
;j
++)
1297 if( (j
%cmap_grid
->grid_spacing
)==0)
1299 fprintf(fp
,"%8.1f\n",idx
);
1303 fprintf(fp
,"%8.3f ",cmap_grid
->cmapdata
[i
].cmap
[j
*4]);
1304 fprintf(fp
,"%8.3f ",cmap_grid
->cmapdata
[i
].cmap
[j
*4+1]);
1305 fprintf(fp
,"%8.3f ",cmap_grid
->cmapdata
[i
].cmap
[j
*4+2]);
1306 fprintf(fp
,"%8.3f\n",cmap_grid
->cmapdata
[i
].cmap
[j
*4+3]);
1315 void pr_atomtypes(FILE *fp
,int indent
,const char *title
,t_atomtypes
*atomtypes
,
1319 if (available(fp
,atomtypes
,indent
,title
))
1321 indent
=pr_title(fp
,indent
,title
);
1322 for(i
=0;i
<atomtypes
->nr
;i
++) {
1323 pr_indent(fp
,indent
);
1325 "atomtype[%3d]={radius=%12.5e, volume=%12.5e, gb_radius=%12.5e, surftens=%12.5e, atomnumber=%4d, S_hct=%12.5e)}\n",
1326 bShowNumbers
?i
:-1,atomtypes
->radius
[i
],atomtypes
->vol
[i
],
1327 atomtypes
->gb_radius
[i
],
1328 atomtypes
->surftens
[i
],atomtypes
->atomnumber
[i
],atomtypes
->S_hct
[i
]);
1333 static void pr_moltype(FILE *fp
,int indent
,const char *title
,
1334 gmx_moltype_t
*molt
,int n
,
1335 gmx_ffparams_t
*ffparams
,
1340 indent
= pr_title_n(fp
,indent
,title
,n
);
1341 (void) pr_indent(fp
,indent
);
1342 (void) fprintf(fp
,"name=\"%s\"\n",*(molt
->name
));
1343 pr_atoms(fp
,indent
,"atoms",&(molt
->atoms
),bShowNumbers
);
1344 pr_block(fp
,indent
,"cgs",&molt
->cgs
, bShowNumbers
);
1345 pr_blocka(fp
,indent
,"excls",&molt
->excls
, bShowNumbers
);
1346 for(j
=0; (j
<F_NRE
); j
++) {
1347 pr_ilist(fp
,indent
,interaction_function
[j
].longname
,
1348 ffparams
->functype
,&molt
->ilist
[j
],bShowNumbers
);
1352 static void pr_molblock(FILE *fp
,int indent
,const char *title
,
1353 gmx_molblock_t
*molb
,int n
,
1354 gmx_moltype_t
*molt
,
1357 indent
= pr_title_n(fp
,indent
,title
,n
);
1358 (void) pr_indent(fp
,indent
);
1359 (void) fprintf(fp
,"%-20s = %d \"%s\"\n",
1360 "moltype",molb
->type
,*(molt
[molb
->type
].name
));
1361 pr_int(fp
,indent
,"#molecules",molb
->nmol
);
1362 pr_int(fp
,indent
,"#atoms_mol",molb
->natoms_mol
);
1363 pr_int(fp
,indent
,"#posres_xA",molb
->nposres_xA
);
1364 if (molb
->nposres_xA
> 0) {
1365 pr_rvecs(fp
,indent
,"posres_xA",molb
->posres_xA
,molb
->nposres_xA
);
1367 pr_int(fp
,indent
,"#posres_xB",molb
->nposres_xB
);
1368 if (molb
->nposres_xB
> 0) {
1369 pr_rvecs(fp
,indent
,"posres_xB",molb
->posres_xB
,molb
->nposres_xB
);
1373 void pr_mtop(FILE *fp
,int indent
,const char *title
,gmx_mtop_t
*mtop
,
1378 if (available(fp
,mtop
,indent
,title
)) {
1379 indent
=pr_title(fp
,indent
,title
);
1380 (void) pr_indent(fp
,indent
);
1381 (void) fprintf(fp
,"name=\"%s\"\n",*(mtop
->name
));
1382 pr_int(fp
,indent
,"#atoms",mtop
->natoms
);
1383 for(mb
=0; mb
<mtop
->nmolblock
; mb
++) {
1384 pr_molblock(fp
,indent
,"molblock",&mtop
->molblock
[mb
],mb
,
1385 mtop
->moltype
,bShowNumbers
);
1387 pr_ffparams(fp
,indent
,"ffparams",&(mtop
->ffparams
),bShowNumbers
);
1388 pr_atomtypes(fp
,indent
,"atomtypes",&(mtop
->atomtypes
),bShowNumbers
);
1389 for(mt
=0; mt
<mtop
->nmoltype
; mt
++) {
1390 pr_moltype(fp
,indent
,"moltype",&mtop
->moltype
[mt
],mt
,
1391 &mtop
->ffparams
,bShowNumbers
);
1393 pr_groups(fp
,indent
,"groups",&mtop
->groups
,bShowNumbers
);
1394 pr_cmap(fp
,indent
,"cmap",&mtop
->cmap_grid
,bShowNumbers
);
1398 void pr_top(FILE *fp
,int indent
,const char *title
,t_topology
*top
, bool bShowNumbers
)
1400 if (available(fp
,top
,indent
,title
)) {
1401 indent
=pr_title(fp
,indent
,title
);
1402 (void) pr_indent(fp
,indent
);
1403 (void) fprintf(fp
,"name=\"%s\"\n",*(top
->name
));
1404 pr_atoms(fp
,indent
,"atoms",&(top
->atoms
),bShowNumbers
);
1405 pr_atomtypes(fp
,indent
,"atomtypes",&(top
->atomtypes
),bShowNumbers
);
1406 pr_block(fp
,indent
,"cgs",&top
->cgs
, bShowNumbers
);
1407 pr_block(fp
,indent
,"mols",&top
->mols
, bShowNumbers
);
1408 pr_blocka(fp
,indent
,"excls",&top
->excls
, bShowNumbers
);
1409 pr_idef(fp
,indent
,"idef",&top
->idef
,bShowNumbers
);
1413 void pr_header(FILE *fp
,int indent
,const char *title
,t_tpxheader
*sh
)
1417 if (available(fp
,sh
,indent
,title
))
1419 indent
=pr_title(fp
,indent
,title
);
1420 pr_indent(fp
,indent
);
1421 fprintf(fp
,"bIr = %spresent\n",sh
->bIr
?"":"not ");
1422 pr_indent(fp
,indent
);
1423 fprintf(fp
,"bBox = %spresent\n",sh
->bBox
?"":"not ");
1424 pr_indent(fp
,indent
);
1425 fprintf(fp
,"bTop = %spresent\n",sh
->bTop
?"":"not ");
1426 pr_indent(fp
,indent
);
1427 fprintf(fp
,"bX = %spresent\n",sh
->bX
?"":"not ");
1428 pr_indent(fp
,indent
);
1429 fprintf(fp
,"bV = %spresent\n",sh
->bV
?"":"not ");
1430 pr_indent(fp
,indent
);
1431 fprintf(fp
,"bF = %spresent\n",sh
->bF
?"":"not ");
1433 pr_indent(fp
,indent
);
1434 fprintf(fp
,"natoms = %d\n",sh
->natoms
);
1435 pr_indent(fp
,indent
);
1436 fprintf(fp
,"lambda = %e\n",sh
->lambda
);
1440 void pr_commrec(FILE *fp
,int indent
,t_commrec
*cr
)
1442 pr_indent(fp
,indent
);
1443 fprintf(fp
,"commrec:\n");
1445 pr_indent(fp
,indent
);
1446 fprintf(fp
,"nodeid = %d\n",cr
->nodeid
);
1447 pr_indent(fp
,indent
);
1448 fprintf(fp
,"nnodes = %d\n",cr
->nnodes
);
1449 pr_indent(fp
,indent
);
1450 fprintf(fp
,"npmenodes = %d\n",cr
->npmenodes
);
1451 pr_indent(fp
,indent
);
1452 fprintf(fp
,"threadid = %d\n",cr
->threadid
);
1453 pr_indent(fp
,indent
);
1454 fprintf(fp
,"nthreads = %d\n",cr
->nthreads
);