2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
5 * Copyright (c) 2001-2004, The GROMACS development team.
6 * Copyright (c) 2013,2014,2015,2016,2018, by the GROMACS development team, led by
7 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
8 * and including many others, as listed in the AUTHORS file in the
9 * top-level source directory and at http://www.gromacs.org.
11 * GROMACS is free software; you can redistribute it and/or
12 * modify it under the terms of the GNU Lesser General Public License
13 * as published by the Free Software Foundation; either version 2.1
14 * of the License, or (at your option) any later version.
16 * GROMACS is distributed in the hope that it will be useful,
17 * but WITHOUT ANY WARRANTY; without even the implied warranty of
18 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
19 * Lesser General Public License for more details.
21 * You should have received a copy of the GNU Lesser General Public
22 * License along with GROMACS; if not, see
23 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
24 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
26 * If you want to redistribute modifications to GROMACS, please
27 * consider that scientific software is very special. Version
28 * control is crucial - bugs must be traceable. We will be happy to
29 * consider code for inclusion in the official distribution, but
30 * derived work must not be called official GROMACS. Details are found
31 * in the README & COPYING files - if they are missing, get the
32 * official version at http://www.gromacs.org.
34 * To help us fund GROMACS development, we humbly ask that you cite
35 * the research papers on the package. Check out http://www.gromacs.org.
43 #include "gromacs/topology/ifunc.h"
44 #include "gromacs/utility/fatalerror.h"
45 #include "gromacs/utility/txtdump.h"
47 static void pr_harm(FILE *fp
, const t_iparams
*iparams
, const char *r
, const char *kr
)
49 fprintf(fp
, "%sA=%12.5e, %sA=%12.5e, %sB=%12.5e, %sB=%12.5e\n",
50 r
, iparams
->harmonic
.rA
, kr
, iparams
->harmonic
.krA
,
51 r
, iparams
->harmonic
.rB
, kr
, iparams
->harmonic
.krB
);
54 void pr_iparams(FILE *fp
, t_functype ftype
, const t_iparams
*iparams
)
60 pr_harm(fp
, iparams
, "th", "ct");
62 case F_CROSS_BOND_BONDS
:
63 fprintf(fp
, "r1e=%15.8e, r2e=%15.8e, krr=%15.8e\n",
64 iparams
->cross_bb
.r1e
, iparams
->cross_bb
.r2e
,
65 iparams
->cross_bb
.krr
);
67 case F_CROSS_BOND_ANGLES
:
68 fprintf(fp
, "r1e=%15.8e, r1e=%15.8e, r3e=%15.8e, krt=%15.8e\n",
69 iparams
->cross_ba
.r1e
, iparams
->cross_ba
.r2e
,
70 iparams
->cross_ba
.r3e
, iparams
->cross_ba
.krt
);
73 fprintf(fp
, "klinA=%15.8e, aA=%15.8e, klinB=%15.8e, aB=%15.8e\n",
74 iparams
->linangle
.klinA
, iparams
->linangle
.aA
,
75 iparams
->linangle
.klinB
, iparams
->linangle
.aB
);
78 fprintf(fp
, "thetaA=%15.8e, kthetaA=%15.8e, r13A=%15.8e, kUBA=%15.8e, thetaB=%15.8e, kthetaB=%15.8e, r13B=%15.8e, kUBB=%15.8e\n", iparams
->u_b
.thetaA
, iparams
->u_b
.kthetaA
, iparams
->u_b
.r13A
, iparams
->u_b
.kUBA
, iparams
->u_b
.thetaB
, iparams
->u_b
.kthetaB
, iparams
->u_b
.r13B
, iparams
->u_b
.kUBB
);
80 case F_QUARTIC_ANGLES
:
81 fprintf(fp
, "theta=%15.8e", iparams
->qangle
.theta
);
82 for (int i
= 0; i
< 5; i
++)
84 fprintf(fp
, ", c%c=%15.8e", '0'+i
, iparams
->qangle
.c
[i
]);
89 fprintf(fp
, "a=%15.8e, b=%15.8e, c=%15.8e\n",
90 iparams
->bham
.a
, iparams
->bham
.b
, iparams
->bham
.c
);
95 pr_harm(fp
, iparams
, "b0", "cb");
98 pr_harm(fp
, iparams
, "xi", "cx");
101 fprintf(fp
, "b0A=%15.8e, cbA=%15.8e, betaA=%15.8e, b0B=%15.8e, cbB=%15.8e, betaB=%15.8e\n",
102 iparams
->morse
.b0A
, iparams
->morse
.cbA
, iparams
->morse
.betaA
,
103 iparams
->morse
.b0B
, iparams
->morse
.cbB
, iparams
->morse
.betaB
);
106 fprintf(fp
, "b0=%15.8e, kb=%15.8e, kcub=%15.8e\n",
107 iparams
->cubic
.b0
, iparams
->cubic
.kb
, iparams
->cubic
.kcub
);
113 fprintf(fp
, "bm=%15.8e, kb=%15.8e\n", iparams
->fene
.bm
, iparams
->fene
.kb
);
116 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",
117 iparams
->restraint
.lowA
, iparams
->restraint
.up1A
,
118 iparams
->restraint
.up2A
, iparams
->restraint
.kA
,
119 iparams
->restraint
.lowB
, iparams
->restraint
.up1B
,
120 iparams
->restraint
.up2B
, iparams
->restraint
.kB
);
126 fprintf(fp
, "tab=%d, kA=%15.8e, kB=%15.8e\n",
127 iparams
->tab
.table
, iparams
->tab
.kA
, iparams
->tab
.kB
);
130 fprintf(fp
, "alpha=%15.8e\n", iparams
->polarize
.alpha
);
133 fprintf(fp
, "alpha=%15.8e drcut=%15.8e khyp=%15.8e\n",
134 iparams
->anharm_polarize
.alpha
,
135 iparams
->anharm_polarize
.drcut
,
136 iparams
->anharm_polarize
.khyp
);
139 fprintf(fp
, "a=%15.8e, alpha1=%15.8e, alpha2=%15.8e, rfac=%15.8e\n",
140 iparams
->thole
.a
, iparams
->thole
.alpha1
, iparams
->thole
.alpha2
,
141 iparams
->thole
.rfac
);
144 fprintf(fp
, "al_x=%15.8e, al_y=%15.8e, al_z=%15.8e, rOH=%9.6f, rHH=%9.6f, rOD=%9.6f\n",
145 iparams
->wpol
.al_x
, iparams
->wpol
.al_y
, iparams
->wpol
.al_z
,
146 iparams
->wpol
.rOH
, iparams
->wpol
.rHH
, iparams
->wpol
.rOD
);
149 fprintf(fp
, "c6=%15.8e, c12=%15.8e\n", iparams
->lj
.c6
, iparams
->lj
.c12
);
152 fprintf(fp
, "c6A=%15.8e, c12A=%15.8e, c6B=%15.8e, c12B=%15.8e\n",
153 iparams
->lj14
.c6A
, iparams
->lj14
.c12A
,
154 iparams
->lj14
.c6B
, iparams
->lj14
.c12B
);
157 fprintf(fp
, "fqq=%15.8e, qi=%15.8e, qj=%15.8e, c6=%15.8e, c12=%15.8e\n",
159 iparams
->ljc14
.qi
, iparams
->ljc14
.qj
,
160 iparams
->ljc14
.c6
, iparams
->ljc14
.c12
);
163 fprintf(fp
, "qi=%15.8e, qj=%15.8e, c6=%15.8e, c12=%15.8e\n",
164 iparams
->ljcnb
.qi
, iparams
->ljcnb
.qj
,
165 iparams
->ljcnb
.c6
, iparams
->ljcnb
.c12
);
171 fprintf(fp
, "phiA=%15.8e, cpA=%15.8e, phiB=%15.8e, cpB=%15.8e, mult=%d\n",
172 iparams
->pdihs
.phiA
, iparams
->pdihs
.cpA
,
173 iparams
->pdihs
.phiB
, iparams
->pdihs
.cpB
,
174 iparams
->pdihs
.mult
);
177 fprintf(fp
, "label=%4d, type=%1d, low=%15.8e, up1=%15.8e, up2=%15.8e, fac=%15.8e)\n",
178 iparams
->disres
.label
, iparams
->disres
.type
,
179 iparams
->disres
.low
, iparams
->disres
.up1
,
180 iparams
->disres
.up2
, iparams
->disres
.kfac
);
183 fprintf(fp
, "ex=%4d, label=%d, power=%4d, c=%15.8e, obs=%15.8e, kfac=%15.8e)\n",
184 iparams
->orires
.ex
, iparams
->orires
.label
, iparams
->orires
.power
,
185 iparams
->orires
.c
, iparams
->orires
.obs
, iparams
->orires
.kfac
);
188 fprintf(fp
, "phiA=%15.8e, dphiA=%15.8e, kfacA=%15.8e, phiB=%15.8e, dphiB=%15.8e, kfacB=%15.8e\n",
189 iparams
->dihres
.phiA
, iparams
->dihres
.dphiA
, iparams
->dihres
.kfacA
,
190 iparams
->dihres
.phiB
, iparams
->dihres
.dphiB
, iparams
->dihres
.kfacB
);
193 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",
194 iparams
->posres
.pos0A
[XX
], iparams
->posres
.pos0A
[YY
],
195 iparams
->posres
.pos0A
[ZZ
], iparams
->posres
.fcA
[XX
],
196 iparams
->posres
.fcA
[YY
], iparams
->posres
.fcA
[ZZ
],
197 iparams
->posres
.pos0B
[XX
], iparams
->posres
.pos0B
[YY
],
198 iparams
->posres
.pos0B
[ZZ
], iparams
->posres
.fcB
[XX
],
199 iparams
->posres
.fcB
[YY
], iparams
->posres
.fcB
[ZZ
]);
202 fprintf(fp
, "pos0=(%15.8e,%15.8e,%15.8e), geometry=%d, r=%15.8e, k=%15.8e\n",
203 iparams
->fbposres
.pos0
[XX
], iparams
->fbposres
.pos0
[YY
],
204 iparams
->fbposres
.pos0
[ZZ
], iparams
->fbposres
.geom
,
205 iparams
->fbposres
.r
, iparams
->fbposres
.k
);
208 for (int i
= 0; i
< NR_RBDIHS
; i
++)
210 fprintf(fp
, "%srbcA[%d]=%15.8e", i
== 0 ? "" : ", ", i
, iparams
->rbdihs
.rbcA
[i
]);
213 for (int i
= 0; i
< NR_RBDIHS
; i
++)
215 fprintf(fp
, "%srbcB[%d]=%15.8e", i
== 0 ? "" : ", ", i
, iparams
->rbdihs
.rbcB
[i
]);
221 /* Use the OPLS -> Ryckaert-Bellemans formula backwards to get
222 * the OPLS potential constants back.
224 const real
*rbcA
= iparams
->rbdihs
.rbcA
;
225 const real
*rbcB
= iparams
->rbdihs
.rbcB
;
228 VA
[3] = -0.25*rbcA
[4];
229 VA
[2] = -0.5*rbcA
[3];
230 VA
[1] = 4.0*VA
[3]-rbcA
[2];
231 VA
[0] = 3.0*VA
[2]-2.0*rbcA
[1];
233 VB
[3] = -0.25*rbcB
[4];
234 VB
[2] = -0.5*rbcB
[3];
235 VB
[1] = 4.0*VB
[3]-rbcB
[2];
236 VB
[0] = 3.0*VB
[2]-2.0*rbcB
[1];
238 for (int i
= 0; i
< NR_FOURDIHS
; i
++)
240 fprintf(fp
, "%sFourA[%d]=%15.8e", i
== 0 ? "" : ", ", i
, VA
[i
]);
243 for (int i
= 0; i
< NR_FOURDIHS
; i
++)
245 fprintf(fp
, "%sFourB[%d]=%15.8e", i
== 0 ? "" : ", ", i
, VB
[i
]);
253 fprintf(fp
, "dA=%15.8e, dB=%15.8e\n", iparams
->constr
.dA
, iparams
->constr
.dB
);
256 fprintf(fp
, "doh=%15.8e, dhh=%15.8e\n", iparams
->settle
.doh
,
257 iparams
->settle
.dhh
);
260 fprintf(fp
, "a=%15.8e\n", iparams
->vsite
.a
);
265 fprintf(fp
, "a=%15.8e, b=%15.8e\n", iparams
->vsite
.a
, iparams
->vsite
.b
);
270 fprintf(fp
, "a=%15.8e, b=%15.8e, c=%15.8e\n",
271 iparams
->vsite
.a
, iparams
->vsite
.b
, iparams
->vsite
.c
);
274 fprintf(fp
, "n=%2d, a=%15.8e\n", iparams
->vsiten
.n
, iparams
->vsiten
.a
);
276 case F_GB12_NOLONGERUSED
:
277 case F_GB13_NOLONGERUSED
:
278 case F_GB14_NOLONGERUSED
:
279 // These could only be generated by grompp, not written in
280 // a .top file. Now that implicit solvent is not
281 // supported, they can't be generated, and the values are
282 // ignored if read from an old .tpr file. So there is
286 fprintf(fp
, "cmapA=%1d, cmapB=%1d\n", iparams
->cmap
.cmapA
, iparams
->cmap
.cmapB
);
289 pr_harm(fp
, iparams
, "ktheta", "costheta0");
292 fprintf(fp
, "phiA=%15.8e, cpA=%15.8e",
293 iparams
->pdihs
.phiA
, iparams
->pdihs
.cpA
);
296 fprintf(fp
, "kphi=%15.8e", iparams
->cbtdihs
.cbtcA
[0]);
297 for (int i
= 1; i
< NR_CBTDIHS
; i
++)
299 fprintf(fp
, ", cbtcA[%d]=%15.8e", i
-1, iparams
->cbtdihs
.cbtcA
[i
]);
304 gmx_fatal(FARGS
, "unknown function type %d (%s) in %s line %d",
305 ftype
, interaction_function
[ftype
].name
, __FILE__
, __LINE__
);
309 void pr_ilist(FILE *fp
, int indent
, const char *title
,
310 const t_functype
*functype
, const t_ilist
*ilist
,
311 gmx_bool bShowNumbers
,
312 gmx_bool bShowParameters
, const t_iparams
*iparams
)
314 int i
, j
, k
, type
, ftype
;
317 if (available(fp
, ilist
, indent
, title
) && ilist
->nr
> 0)
319 indent
= pr_title(fp
, indent
, title
);
320 pr_indent(fp
, indent
);
321 fprintf(fp
, "nr: %d\n", ilist
->nr
);
324 pr_indent(fp
, indent
);
325 fprintf(fp
, "iatoms:\n");
326 iatoms
= ilist
->iatoms
;
327 for (i
= j
= 0; i
< ilist
->nr
; )
329 pr_indent(fp
, indent
+INDENT
);
331 ftype
= functype
[type
];
334 fprintf(fp
, "%d type=%d ", j
, type
);
337 printf("(%s)", interaction_function
[ftype
].name
);
338 for (k
= 0; k
< interaction_function
[ftype
].nratoms
; k
++)
340 fprintf(fp
, " %3d", *(iatoms
++));
345 pr_iparams(fp
, ftype
, &iparams
[type
]);
348 i
+= 1+interaction_function
[ftype
].nratoms
;
354 static void pr_cmap(FILE *fp
, int indent
, const char *title
,
355 const gmx_cmap_t
*cmap_grid
, gmx_bool bShowNumbers
)
360 dx
= 360.0 / cmap_grid
->grid_spacing
;
361 nelem
= cmap_grid
->grid_spacing
*cmap_grid
->grid_spacing
;
363 if (available(fp
, cmap_grid
, indent
, title
))
365 fprintf(fp
, "%s\n", title
);
367 for (i
= 0; i
< cmap_grid
->ngrid
; i
++)
370 fprintf(fp
, "%8s %8s %8s %8s\n", "V", "dVdx", "dVdy", "d2dV");
372 fprintf(fp
, "grid[%3d]={\n", bShowNumbers
? i
: -1);
374 for (j
= 0; j
< nelem
; j
++)
376 if ( (j
%cmap_grid
->grid_spacing
) == 0)
378 fprintf(fp
, "%8.1f\n", idx
);
382 fprintf(fp
, "%8.3f ", cmap_grid
->cmapdata
[i
].cmap
[j
*4]);
383 fprintf(fp
, "%8.3f ", cmap_grid
->cmapdata
[i
].cmap
[j
*4+1]);
384 fprintf(fp
, "%8.3f ", cmap_grid
->cmapdata
[i
].cmap
[j
*4+2]);
385 fprintf(fp
, "%8.3f\n", cmap_grid
->cmapdata
[i
].cmap
[j
*4+3]);
393 void pr_ffparams(FILE *fp
, int indent
, const char *title
,
394 const gmx_ffparams_t
*ffparams
,
395 gmx_bool bShowNumbers
)
399 indent
= pr_title(fp
, indent
, title
);
400 pr_indent(fp
, indent
);
401 fprintf(fp
, "atnr=%d\n", ffparams
->atnr
);
402 pr_indent(fp
, indent
);
403 fprintf(fp
, "ntypes=%d\n", ffparams
->ntypes
);
404 for (i
= 0; i
< ffparams
->ntypes
; i
++)
406 pr_indent(fp
, indent
+INDENT
);
407 fprintf(fp
, "functype[%d]=%s, ",
408 bShowNumbers
? i
: -1,
409 interaction_function
[ffparams
->functype
[i
]].name
);
410 pr_iparams(fp
, ffparams
->functype
[i
], &ffparams
->iparams
[i
]);
412 pr_double(fp
, indent
, "reppow", ffparams
->reppow
);
413 pr_real(fp
, indent
, "fudgeQQ", ffparams
->fudgeQQ
);
414 pr_cmap(fp
, indent
, "cmap", &ffparams
->cmap_grid
, bShowNumbers
);
417 void pr_idef(FILE *fp
, int indent
, const char *title
, const t_idef
*idef
,
418 gmx_bool bShowNumbers
, gmx_bool bShowParameters
)
422 if (available(fp
, idef
, indent
, title
))
424 indent
= pr_title(fp
, indent
, title
);
425 pr_indent(fp
, indent
);
426 fprintf(fp
, "atnr=%d\n", idef
->atnr
);
427 pr_indent(fp
, indent
);
428 fprintf(fp
, "ntypes=%d\n", idef
->ntypes
);
429 for (i
= 0; i
< idef
->ntypes
; i
++)
431 pr_indent(fp
, indent
+INDENT
);
432 fprintf(fp
, "functype[%d]=%s, ",
433 bShowNumbers
? i
: -1,
434 interaction_function
[idef
->functype
[i
]].name
);
435 pr_iparams(fp
, idef
->functype
[i
], &idef
->iparams
[i
]);
437 pr_real(fp
, indent
, "fudgeQQ", idef
->fudgeQQ
);
439 for (j
= 0; (j
< F_NRE
); j
++)
441 pr_ilist(fp
, indent
, interaction_function
[j
].longname
,
442 idef
->functype
, &idef
->il
[j
], bShowNumbers
,
443 bShowParameters
, idef
->iparams
);