From 33bd52ff41b7ec5d5c98fbb1f2523a92c57b4ac8 Mon Sep 17 00:00:00 2001 From: Berk Hess Date: Tue, 21 Jun 2016 12:00:24 +0200 Subject: [PATCH] Add gmx dump option for printing parameters For comparing topologies of the same system with different force fields it can be very useful to print bonded parameters instead of the parameter type index. Change-Id: I06a3b65b695ace4de1491f2d58abb6b7e752693d --- src/gromacs/tools/dump.cpp | 15 ++++++++++----- src/gromacs/topology/idef.cpp | 32 ++++++++++++++++++------------- src/gromacs/topology/idef.h | 9 ++++++--- src/gromacs/topology/topology.cpp | 17 +++++++++------- src/gromacs/topology/topology.h | 5 +++-- 5 files changed, 48 insertions(+), 30 deletions(-) diff --git a/src/gromacs/tools/dump.cpp b/src/gromacs/tools/dump.cpp index 07d1c0fc8d..1e9da4acfd 100644 --- a/src/gromacs/tools/dump.cpp +++ b/src/gromacs/tools/dump.cpp @@ -71,8 +71,11 @@ #include "gromacs/utility/smalloc.h" #include "gromacs/utility/txtdump.h" -static void list_tpx(const char *fn, gmx_bool bShowNumbers, const char *mdpfn, - gmx_bool bSysTop) +static void list_tpx(const char *fn, + gmx_bool bShowNumbers, + gmx_bool bShowParameters, + const char *mdpfn, + gmx_bool bSysTop) { FILE *gp; int indent, i, j, **gcount, atot; @@ -118,11 +121,11 @@ static void list_tpx(const char *fn, gmx_bool bShowNumbers, const char *mdpfn, if (!bSysTop) { - pr_mtop(stdout, indent, "topology", &(mtop), bShowNumbers); + pr_mtop(stdout, indent, "topology", &(mtop), bShowNumbers, bShowParameters); } else { - pr_top(stdout, indent, "topology", &(top), bShowNumbers); + pr_top(stdout, indent, "topology", &(top), bShowNumbers, bShowParameters); } pr_rvecs(stdout, indent, "box", tpx.bBox ? state.box : NULL, DIM); @@ -623,9 +626,11 @@ int gmx_dump(int argc, char *argv[]) gmx_output_env_t *oenv; /* Command line options */ static gmx_bool bShowNumbers = TRUE; + static gmx_bool bShowParams = FALSE; static gmx_bool bSysTop = FALSE; t_pargs pa[] = { { "-nr", FALSE, etBOOL, {&bShowNumbers}, "Show index numbers in output (leaving them out makes comparison easier, but creates a useless topology)" }, + { "-param", FALSE, etBOOL, {&bShowParams}, "Show parameters for each bonded interaction (for comparing dumps, it is useful to combine this with -nonr)" }, { "-sys", FALSE, etBOOL, {&bSysTop}, "List the atoms and bonded interactions for the whole system instead of for each molecule type" } }; @@ -638,7 +643,7 @@ int gmx_dump(int argc, char *argv[]) if (ftp2bSet(efTPR, NFILE, fnm)) { - list_tpx(ftp2fn(efTPR, NFILE, fnm), bShowNumbers, + list_tpx(ftp2fn(efTPR, NFILE, fnm), bShowNumbers, bShowParams, ftp2fn_null(efMDP, NFILE, fnm), bSysTop); } else if (ftp2bSet(efTRX, NFILE, fnm)) diff --git a/src/gromacs/topology/idef.cpp b/src/gromacs/topology/idef.cpp index 0ad6cbba67..4806c2a9ba 100644 --- a/src/gromacs/topology/idef.cpp +++ b/src/gromacs/topology/idef.cpp @@ -3,7 +3,7 @@ * * Copyright (c) 1991-2000, University of Groningen, The Netherlands. * Copyright (c) 2001-2004, The GROMACS development team. - * Copyright (c) 2013,2014,2015, by the GROMACS development team, led by + * Copyright (c) 2013,2014,2015,2016, by the GROMACS development team, led by * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, * and including many others, as listed in the AUTHORS file in the * top-level source directory and at http://www.gromacs.org. @@ -303,7 +303,9 @@ void pr_iparams(FILE *fp, t_functype ftype, const t_iparams *iparams) } void pr_ilist(FILE *fp, int indent, const char *title, - const t_functype *functype, const t_ilist *ilist, gmx_bool bShowNumbers) + const t_functype *functype, const t_ilist *ilist, + gmx_bool bShowNumbers, + gmx_bool bShowParameters, const t_iparams *iparams) { int i, j, k, type, ftype; t_iatom *iatoms; @@ -320,24 +322,26 @@ void pr_ilist(FILE *fp, int indent, const char *title, iatoms = ilist->iatoms; for (i = j = 0; i < ilist->nr; ) { -#ifndef DEBUG pr_indent(fp, indent+INDENT); type = *(iatoms++); ftype = functype[type]; - fprintf(fp, "%d type=%d (%s)", - bShowNumbers ? j : -1, bShowNumbers ? type : -1, - interaction_function[ftype].name); + if (bShowNumbers) + { + fprintf(fp, "%d type=%d ", j, type); + } j++; + printf("(%s)", interaction_function[ftype].name); for (k = 0; k < interaction_function[ftype].nratoms; k++) { - fprintf(fp, " %d", *(iatoms++)); + fprintf(fp, " %3d", *(iatoms++)); + } + if (bShowParameters) + { + fprintf(fp, " "); + pr_iparams(fp, ftype, &iparams[type]); } fprintf(fp, "\n"); i += 1+interaction_function[ftype].nratoms; -#else - fprintf(fp, "%5d%5d\n", i, iatoms[i]); - i++; -#endif } } } @@ -406,7 +410,8 @@ void pr_ffparams(FILE *fp, int indent, const char *title, pr_cmap(fp, indent, "cmap", &ffparams->cmap_grid, bShowNumbers); } -void pr_idef(FILE *fp, int indent, const char *title, const t_idef *idef, gmx_bool bShowNumbers) +void pr_idef(FILE *fp, int indent, const char *title, const t_idef *idef, + gmx_bool bShowNumbers, gmx_bool bShowParameters) { int i, j; @@ -430,7 +435,8 @@ void pr_idef(FILE *fp, int indent, const char *title, const t_idef *idef, gmx_bo for (j = 0; (j < F_NRE); j++) { pr_ilist(fp, indent, interaction_function[j].longname, - idef->functype, &idef->il[j], bShowNumbers); + idef->functype, &idef->il[j], bShowNumbers, + bShowParameters, idef->iparams); } } } diff --git a/src/gromacs/topology/idef.h b/src/gromacs/topology/idef.h index 4446a3f170..467d7d45f8 100644 --- a/src/gromacs/topology/idef.h +++ b/src/gromacs/topology/idef.h @@ -3,7 +3,7 @@ * * Copyright (c) 1991-2000, University of Groningen, The Netherlands. * Copyright (c) 2001-2004, The GROMACS development team. - * Copyright (c) 2013,2014,2015, by the GROMACS development team, led by + * Copyright (c) 2013,2014,2015,2016, by the GROMACS development team, led by * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, * and including many others, as listed in the AUTHORS file in the * top-level source directory and at http://www.gromacs.org. @@ -406,10 +406,13 @@ typedef struct t_idef void pr_iparams(FILE *fp, t_functype ftype, const t_iparams *iparams); void pr_ilist(FILE *fp, int indent, const char *title, - const t_functype *functype, const t_ilist *ilist, gmx_bool bShowNumbers); + const t_functype *functype, const t_ilist *ilist, + gmx_bool bShowNumbers, + gmx_bool bShowParameters, const t_iparams *iparams); void pr_ffparams(FILE *fp, int indent, const char *title, const gmx_ffparams_t *ffparams, gmx_bool bShowNumbers); -void pr_idef(FILE *fp, int indent, const char *title, const t_idef *idef, gmx_bool bShowNumbers); +void pr_idef(FILE *fp, int indent, const char *title, const t_idef *idef, + gmx_bool bShowNumbers, gmx_bool bShowParameters); #ifdef __cplusplus } diff --git a/src/gromacs/topology/topology.cpp b/src/gromacs/topology/topology.cpp index 27eed7a3c3..f113900e13 100644 --- a/src/gromacs/topology/topology.cpp +++ b/src/gromacs/topology/topology.cpp @@ -317,7 +317,7 @@ static void pr_groups(FILE *fp, int indent, static void pr_moltype(FILE *fp, int indent, const char *title, const gmx_moltype_t *molt, int n, const gmx_ffparams_t *ffparams, - gmx_bool bShowNumbers) + gmx_bool bShowNumbers, gmx_bool bShowParameters) { int j; @@ -330,7 +330,8 @@ static void pr_moltype(FILE *fp, int indent, const char *title, for (j = 0; (j < F_NRE); j++) { pr_ilist(fp, indent, interaction_function[j].longname, - ffparams->functype, &molt->ilist[j], bShowNumbers); + ffparams->functype, &molt->ilist[j], + bShowNumbers, bShowParameters, ffparams->iparams); } } @@ -357,7 +358,7 @@ static void pr_molblock(FILE *fp, int indent, const char *title, } void pr_mtop(FILE *fp, int indent, const char *title, const gmx_mtop_t *mtop, - gmx_bool bShowNumbers) + gmx_bool bShowNumbers, gmx_bool bShowParameters) { int mt, mb, j; @@ -380,7 +381,8 @@ void pr_mtop(FILE *fp, int indent, const char *title, const gmx_mtop_t *mtop, { pr_ilist(fp, indent, interaction_function[j].longname, mtop->ffparams.functype, - &mtop->intermolecular_ilist[j], bShowNumbers); + &mtop->intermolecular_ilist[j], + bShowNumbers, bShowParameters, mtop->ffparams.iparams); } } pr_ffparams(fp, indent, "ffparams", &(mtop->ffparams), bShowNumbers); @@ -388,13 +390,14 @@ void pr_mtop(FILE *fp, int indent, const char *title, const gmx_mtop_t *mtop, for (mt = 0; mt < mtop->nmoltype; mt++) { pr_moltype(fp, indent, "moltype", &mtop->moltype[mt], mt, - &mtop->ffparams, bShowNumbers); + &mtop->ffparams, bShowNumbers, bShowParameters); } pr_groups(fp, indent, &mtop->groups, bShowNumbers); } } -void pr_top(FILE *fp, int indent, const char *title, const t_topology *top, gmx_bool bShowNumbers) +void pr_top(FILE *fp, int indent, const char *title, const t_topology *top, + gmx_bool bShowNumbers, gmx_bool bShowParameters) { if (available(fp, top, indent, title)) { @@ -408,7 +411,7 @@ void pr_top(FILE *fp, int indent, const char *title, const t_topology *top, gmx_ pr_str(fp, indent, "bIntermolecularInteractions", gmx::boolToString(top->bIntermolecularInteractions)); pr_blocka(fp, indent, "excls", &top->excls, bShowNumbers); - pr_idef(fp, indent, "idef", &top->idef, bShowNumbers); + pr_idef(fp, indent, "idef", &top->idef, bShowNumbers, bShowParameters); } } diff --git a/src/gromacs/topology/topology.h b/src/gromacs/topology/topology.h index ed2b5ae7de..89c40bfc08 100644 --- a/src/gromacs/topology/topology.h +++ b/src/gromacs/topology/topology.h @@ -160,8 +160,9 @@ bool gmx_mtop_has_atomtypes(const gmx_mtop_t *mtop); bool gmx_mtop_has_pdbinfo(const gmx_mtop_t *mtop); void pr_mtop(FILE *fp, int indent, const char *title, const gmx_mtop_t *mtop, - gmx_bool bShowNumbers); -void pr_top(FILE *fp, int indent, const char *title, const t_topology *top, gmx_bool bShowNumbers); + gmx_bool bShowNumbers, gmx_bool bShowParameters); +void pr_top(FILE *fp, int indent, const char *title, const t_topology *top, + gmx_bool bShowNumbers, gmx_bool bShowParameters); void cmp_top(FILE *fp, const t_topology *t1, const t_topology *t2, real ftol, real abstol); void cmp_groups(FILE *fp, const gmx_groups_t *g0, const gmx_groups_t *g1, -- 2.22.0