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, 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.
37 /* This file is completely threadsafe - keep it that way! */
44 #include "gromacs/fileio/enxio.h"
45 #include "gromacs/fileio/tpxio.h"
46 #include "gromacs/fileio/trxio.h"
47 #include "gromacs/legacyheaders/macros.h"
48 #include "gromacs/legacyheaders/names.h"
49 #include "gromacs/legacyheaders/txtdump.h"
50 #include "gromacs/topology/mtop_util.h"
51 #include "gromacs/utility/cstringutil.h"
52 #include "gromacs/utility/fatalerror.h"
53 #include "gromacs/utility/futil.h"
54 #include "gromacs/utility/smalloc.h"
56 static void cmp_int(FILE *fp
, const char *s
, int index
, int i1
, int i2
)
62 fprintf(fp
, "%s[%d] (%d - %d)\n", s
, index
, i1
, i2
);
66 fprintf(fp
, "%s (%d - %d)\n", s
, i1
, i2
);
71 static void cmp_int64(FILE *fp
, const char *s
, gmx_int64_t i1
, gmx_int64_t i2
)
75 fprintf(fp
, "%s (", s
);
76 fprintf(fp
, "%"GMX_PRId64
, i1
);
78 fprintf(fp
, "%"GMX_PRId64
, i2
);
83 static void cmp_us(FILE *fp
, const char *s
, int index
, unsigned short i1
, unsigned short i2
)
89 fprintf(fp
, "%s[%d] (%hu - %hu)\n", s
, index
, i1
, i2
);
93 fprintf(fp
, "%s (%hu - %hu)\n", s
, i1
, i2
);
98 static void cmp_uc(FILE *fp
, const char *s
, int index
, unsigned char i1
, unsigned char i2
)
104 fprintf(fp
, "%s[%d] (%d - %d)\n", s
, index
, i1
, i2
);
108 fprintf(fp
, "%s (%d - %d)\n", s
, i1
, i2
);
113 static gmx_bool
cmp_bool(FILE *fp
, const char *s
, int index
, gmx_bool b1
, gmx_bool b2
)
135 fprintf(fp
, "%s[%d] (%s - %s)\n", s
, index
,
136 bool_names
[b1
], bool_names
[b2
]);
140 fprintf(fp
, "%s (%s - %s)\n", s
,
141 bool_names
[b1
], bool_names
[b2
]);
147 static void cmp_str(FILE *fp
, const char *s
, int index
,
148 const char *s1
, const char *s2
)
150 if (strcmp(s1
, s2
) != 0)
154 fprintf(fp
, "%s[%d] (%s - %s)\n", s
, index
, s1
, s2
);
158 fprintf(fp
, "%s (%s - %s)\n", s
, s1
, s2
);
163 static gmx_bool
equal_real(real i1
, real i2
, real ftol
, real abstol
)
165 return ( ( 2*fabs(i1
- i2
) <= (fabs(i1
) + fabs(i2
))*ftol
) || fabs(i1
-i2
) <= abstol
);
168 static gmx_bool
equal_float(float i1
, float i2
, float ftol
, float abstol
)
170 return ( ( 2*fabs(i1
- i2
) <= (fabs(i1
) + fabs(i2
))*ftol
) || fabs(i1
-i2
) <= abstol
);
173 static gmx_bool
equal_double(double i1
, double i2
, real ftol
, real abstol
)
175 return ( ( 2*fabs(i1
- i2
) <= (fabs(i1
) + fabs(i2
))*ftol
) || fabs(i1
-i2
) <= abstol
);
179 cmp_real(FILE *fp
, const char *s
, int index
, real i1
, real i2
, real ftol
, real abstol
)
181 if (!equal_real(i1
, i2
, ftol
, abstol
))
185 fprintf(fp
, "%s[%2d] (%e - %e)\n", s
, index
, i1
, i2
);
189 fprintf(fp
, "%s (%e - %e)\n", s
, i1
, i2
);
195 cmp_float(FILE *fp
, const char *s
, int index
, float i1
, float i2
, float ftol
, float abstol
)
197 if (!equal_float(i1
, i2
, ftol
, abstol
))
201 fprintf(fp
, "%s[%2d] (%e - %e)\n", s
, index
, i1
, i2
);
205 fprintf(fp
, "%s (%e - %e)\n", s
, i1
, i2
);
213 cmp_double(FILE *fp
, const char *s
, int index
, double i1
, double i2
, double ftol
, double abstol
)
215 if (!equal_double(i1
, i2
, ftol
, abstol
))
219 fprintf(fp
, "%s[%2d] (%16.9e - %16.9e)\n", s
, index
, i1
, i2
);
223 fprintf(fp
, "%s (%16.9e - %16.9e)\n", s
, i1
, i2
);
228 static void cmp_rvec(FILE *fp
, const char *s
, int index
, rvec i1
, rvec i2
, real ftol
, real abstol
)
230 if (!equal_real(i1
[XX
], i2
[XX
], ftol
, abstol
) ||
231 !equal_real(i1
[YY
], i2
[YY
], ftol
, abstol
) ||
232 !equal_real(i1
[ZZ
], i2
[ZZ
], ftol
, abstol
))
236 fprintf(fp
, "%s[%5d] (%12.5e %12.5e %12.5e) - (%12.5e %12.5e %12.5e)\n",
237 s
, index
, i1
[XX
], i1
[YY
], i1
[ZZ
], i2
[XX
], i2
[YY
], i2
[ZZ
]);
241 fprintf(fp
, "%s (%12.5e %12.5e %12.5e) - (%12.5e %12.5e %12.5e)\n",
242 s
, i1
[XX
], i1
[YY
], i1
[ZZ
], i2
[XX
], i2
[YY
], i2
[ZZ
]);
247 static void cmp_ivec(FILE *fp
, const char *s
, int index
, ivec i1
, ivec i2
)
249 if ((i1
[XX
] != i2
[XX
]) || (i1
[YY
] != i2
[YY
]) || (i1
[ZZ
] != i2
[ZZ
]))
253 fprintf(fp
, "%s[%5d] (%8d,%8d,%8d - %8d,%8d,%8d)\n", s
, index
,
254 i1
[XX
], i1
[YY
], i1
[ZZ
], i2
[XX
], i2
[YY
], i2
[ZZ
]);
258 fprintf(fp
, "%s (%8d,%8d,%8d - %8d,%8d,%8d)\n", s
,
259 i1
[XX
], i1
[YY
], i1
[ZZ
], i2
[XX
], i2
[YY
], i2
[ZZ
]);
264 static void cmp_ilist(FILE *fp
, int ftype
, t_ilist
*il1
, t_ilist
*il2
)
269 fprintf(fp
, "comparing ilist %s\n", interaction_function
[ftype
].name
);
270 sprintf(buf
, "%s->nr", interaction_function
[ftype
].name
);
271 cmp_int(fp
, buf
, -1, il1
->nr
, il2
->nr
);
272 sprintf(buf
, "%s->iatoms", interaction_function
[ftype
].name
);
273 if (((il1
->nr
> 0) && (!il1
->iatoms
)) ||
274 ((il2
->nr
> 0) && (!il2
->iatoms
)) ||
275 ((il1
->nr
!= il2
->nr
)))
277 fprintf(fp
, "Comparing radically different topologies - %s is different\n",
282 for (i
= 0; (i
< il1
->nr
); i
++)
284 cmp_int(fp
, buf
, i
, il1
->iatoms
[i
], il2
->iatoms
[i
]);
289 void cmp_iparm(FILE *fp
, const char *s
, t_functype ft
,
290 t_iparams ip1
, t_iparams ip2
, real ftol
, real abstol
)
296 for (i
= 0; i
< MAXFORCEPARAM
&& !bDiff
; i
++)
298 bDiff
= !equal_real(ip1
.generic
.buf
[i
], ip2
.generic
.buf
[i
], ftol
, abstol
);
302 fprintf(fp
, "%s1: ", s
);
303 pr_iparams(fp
, ft
, &ip1
);
304 fprintf(fp
, "%s2: ", s
);
305 pr_iparams(fp
, ft
, &ip2
);
309 void cmp_iparm_AB(FILE *fp
, const char *s
, t_functype ft
, t_iparams ip1
, real ftol
, real abstol
)
311 int nrfpA
, nrfpB
, p0
, i
;
314 /* Normally the first parameter is perturbable */
316 nrfpA
= interaction_function
[ft
].nrfpA
;
317 nrfpB
= interaction_function
[ft
].nrfpB
;
322 else if (interaction_function
[ft
].flags
& IF_TABULATED
)
324 /* For tabulated interactions only the second parameter is perturbable */
329 for (i
= 0; i
< nrfpB
&& !bDiff
; i
++)
331 bDiff
= !equal_real(ip1
.generic
.buf
[p0
+i
], ip1
.generic
.buf
[nrfpA
+i
], ftol
, abstol
);
335 fprintf(fp
, "%s: ", s
);
336 pr_iparams(fp
, ft
, &ip1
);
340 static void cmp_cmap(FILE *fp
, const gmx_cmap_t
*cmap1
, const gmx_cmap_t
*cmap2
, real ftol
, real abstol
)
342 cmp_int(fp
, "cmap ngrid", -1, cmap1
->ngrid
, cmap2
->ngrid
);
343 cmp_int(fp
, "cmap grid_spacing", -1, cmap1
->grid_spacing
, cmap2
->grid_spacing
);
344 if (cmap1
->ngrid
== cmap2
->ngrid
&&
345 cmap1
->grid_spacing
== cmap2
->grid_spacing
)
349 for (g
= 0; g
< cmap1
->ngrid
; g
++)
353 fprintf(fp
, "comparing cmap %d\n", g
);
355 for (i
= 0; i
< 4*cmap1
->grid_spacing
*cmap1
->grid_spacing
; i
++)
357 cmp_real(fp
, "", i
, cmap1
->cmapdata
[g
].cmap
[i
], cmap2
->cmapdata
[g
].cmap
[i
], ftol
, abstol
);
363 static void cmp_idef(FILE *fp
, t_idef
*id1
, t_idef
*id2
, real ftol
, real abstol
)
366 char buf1
[64], buf2
[64];
368 fprintf(fp
, "comparing idef\n");
371 cmp_int(fp
, "idef->ntypes", -1, id1
->ntypes
, id2
->ntypes
);
372 cmp_int(fp
, "idef->atnr", -1, id1
->atnr
, id2
->atnr
);
373 for (i
= 0; (i
< min(id1
->ntypes
, id2
->ntypes
)); i
++)
375 sprintf(buf1
, "idef->functype[%d]", i
);
376 sprintf(buf2
, "idef->iparam[%d]", i
);
377 cmp_int(fp
, buf1
, i
, (int)id1
->functype
[i
], (int)id2
->functype
[i
]);
378 cmp_iparm(fp
, buf2
, id1
->functype
[i
],
379 id1
->iparams
[i
], id2
->iparams
[i
], ftol
, abstol
);
381 cmp_real(fp
, "fudgeQQ", -1, id1
->fudgeQQ
, id2
->fudgeQQ
, ftol
, abstol
);
382 cmp_cmap(fp
, &id1
->cmap_grid
, &id2
->cmap_grid
, ftol
, abstol
);
383 for (i
= 0; (i
< F_NRE
); i
++)
385 cmp_ilist(fp
, i
, &(id1
->il
[i
]), &(id2
->il
[i
]));
390 for (i
= 0; (i
< id1
->ntypes
); i
++)
392 cmp_iparm_AB(fp
, "idef->iparam", id1
->functype
[i
], id1
->iparams
[i
], ftol
, abstol
);
397 static void cmp_block(FILE *fp
, t_block
*b1
, t_block
*b2
, const char *s
)
402 fprintf(fp
, "comparing block %s\n", s
);
403 sprintf(buf
, "%s.nr", s
);
404 cmp_int(fp
, buf
, -1, b1
->nr
, b2
->nr
);
407 static void cmp_blocka(FILE *fp
, t_blocka
*b1
, t_blocka
*b2
, const char *s
)
412 fprintf(fp
, "comparing blocka %s\n", s
);
413 sprintf(buf
, "%s.nr", s
);
414 cmp_int(fp
, buf
, -1, b1
->nr
, b2
->nr
);
415 sprintf(buf
, "%s.nra", s
);
416 cmp_int(fp
, buf
, -1, b1
->nra
, b2
->nra
);
419 static void cmp_atom(FILE *fp
, int index
, t_atom
*a1
, t_atom
*a2
, real ftol
, real abstol
)
426 cmp_us(fp
, "atom.type", index
, a1
->type
, a2
->type
);
427 cmp_us(fp
, "atom.ptype", index
, a1
->ptype
, a2
->ptype
);
428 cmp_int(fp
, "atom.resind", index
, a1
->resind
, a2
->resind
);
429 cmp_int(fp
, "atom.atomnumber", index
, a1
->atomnumber
, a2
->atomnumber
);
430 cmp_real(fp
, "atom.m", index
, a1
->m
, a2
->m
, ftol
, abstol
);
431 cmp_real(fp
, "atom.q", index
, a1
->q
, a2
->q
, ftol
, abstol
);
432 cmp_us(fp
, "atom.typeB", index
, a1
->typeB
, a2
->typeB
);
433 cmp_real(fp
, "atom.mB", index
, a1
->mB
, a2
->mB
, ftol
, abstol
);
434 cmp_real(fp
, "atom.qB", index
, a1
->qB
, a2
->qB
, ftol
, abstol
);
438 cmp_us(fp
, "atom.type", index
, a1
->type
, a1
->typeB
);
439 cmp_real(fp
, "atom.m", index
, a1
->m
, a1
->mB
, ftol
, abstol
);
440 cmp_real(fp
, "atom.q", index
, a1
->q
, a1
->qB
, ftol
, abstol
);
444 static void cmp_atoms(FILE *fp
, t_atoms
*a1
, t_atoms
*a2
, real ftol
, real abstol
)
448 fprintf(fp
, "comparing atoms\n");
452 cmp_int(fp
, "atoms->nr", -1, a1
->nr
, a2
->nr
);
453 for (i
= 0; (i
< a1
->nr
); i
++)
455 cmp_atom(fp
, i
, &(a1
->atom
[i
]), &(a2
->atom
[i
]), ftol
, abstol
);
460 for (i
= 0; (i
< a1
->nr
); i
++)
462 cmp_atom(fp
, i
, &(a1
->atom
[i
]), NULL
, ftol
, abstol
);
467 static void cmp_top(FILE *fp
, t_topology
*t1
, t_topology
*t2
, real ftol
, real abstol
)
471 fprintf(fp
, "comparing top\n");
474 cmp_idef(fp
, &(t1
->idef
), &(t2
->idef
), ftol
, abstol
);
475 cmp_atoms(fp
, &(t1
->atoms
), &(t2
->atoms
), ftol
, abstol
);
476 cmp_block(fp
, &t1
->cgs
, &t2
->cgs
, "cgs");
477 cmp_block(fp
, &t1
->mols
, &t2
->mols
, "mols");
478 cmp_bool(fp
, "bIntermolecularInteractions", -1, t1
->bIntermolecularInteractions
, t2
->bIntermolecularInteractions
);
479 cmp_blocka(fp
, &t1
->excls
, &t2
->excls
, "excls");
483 cmp_idef(fp
, &(t1
->idef
), NULL
, ftol
, abstol
);
484 cmp_atoms(fp
, &(t1
->atoms
), NULL
, ftol
, abstol
);
488 static void cmp_groups(FILE *fp
, gmx_groups_t
*g0
, gmx_groups_t
*g1
,
489 int natoms0
, int natoms1
)
494 fprintf(fp
, "comparing groups\n");
496 for (i
= 0; i
< egcNR
; i
++)
498 sprintf(buf
, "grps[%d].nr", i
);
499 cmp_int(fp
, buf
, -1, g0
->grps
[i
].nr
, g1
->grps
[i
].nr
);
500 if (g0
->grps
[i
].nr
== g1
->grps
[i
].nr
)
502 for (j
= 0; j
< g0
->grps
[i
].nr
; j
++)
504 sprintf(buf
, "grps[%d].name[%d]", i
, j
);
506 *g0
->grpname
[g0
->grps
[i
].nm_ind
[j
]],
507 *g1
->grpname
[g1
->grps
[i
].nm_ind
[j
]]);
510 cmp_int(fp
, "ngrpnr", i
, g0
->ngrpnr
[i
], g1
->ngrpnr
[i
]);
511 if (g0
->ngrpnr
[i
] == g1
->ngrpnr
[i
] && natoms0
== natoms1
&&
512 (g0
->grpnr
[i
] != NULL
|| g1
->grpnr
[i
] != NULL
))
514 for (j
= 0; j
< natoms0
; j
++)
516 cmp_int(fp
, gtypes
[i
], j
, ggrpnr(g0
, i
, j
), ggrpnr(g1
, i
, j
));
520 /* We have compared the names in the groups lists,
521 * so we can skip the grpname list comparison.
525 static void cmp_rvecs_rmstol(FILE *fp
, const char *title
, int n
, rvec x1
[], rvec x2
[],
526 real ftol
, real abstol
)
531 /* For a vector you are usally not interested in a relative difference
532 * on a component that is very small compared to the other components.
533 * Therefore we do the relative comparision relative to the RMS component.
536 for (i
= 0; (i
< n
); i
++)
538 for (m
= 0; m
< DIM
; m
++)
540 rms
+= x1
[i
][m
]*x1
[i
][m
] + x2
[i
][m
]*x2
[i
][m
];
543 rms
= sqrt(rms
/(2*n
*DIM
));
545 /* Convert the relative tolerance into an absolute tolerance */
546 if (ftol
*rms
< abstol
)
551 /* And now do the actual comparision */
552 for (i
= 0; (i
< n
); i
++)
554 cmp_rvec(fp
, title
, i
, x1
[i
], x2
[i
], 0.0, abstol
);
558 static void cmp_rvecs(FILE *fp
, const char *title
, int n
, rvec x1
[], rvec x2
[],
559 gmx_bool bRMSD
, real ftol
, real abstol
)
567 for (i
= 0; (i
< n
); i
++)
569 for (m
= 0; m
< DIM
; m
++)
571 d
= x1
[i
][m
] - x2
[i
][m
];
575 fprintf(fp
, "%s RMSD %g\n", title
, sqrt(ssd
/n
));
579 cmp_rvecs_rmstol(fp
, title
, n
, x1
, x2
, ftol
, abstol
);
583 static void cmp_grpopts(FILE *fp
, t_grpopts
*opt1
, t_grpopts
*opt2
, real ftol
, real abstol
)
586 char buf1
[256], buf2
[256];
588 cmp_int(fp
, "inputrec->grpopts.ngtc", -1, opt1
->ngtc
, opt2
->ngtc
);
589 cmp_int(fp
, "inputrec->grpopts.ngacc", -1, opt1
->ngacc
, opt2
->ngacc
);
590 cmp_int(fp
, "inputrec->grpopts.ngfrz", -1, opt1
->ngfrz
, opt2
->ngfrz
);
591 cmp_int(fp
, "inputrec->grpopts.ngener", -1, opt1
->ngener
, opt2
->ngener
);
592 for (i
= 0; (i
< min(opt1
->ngtc
, opt2
->ngtc
)); i
++)
594 cmp_real(fp
, "inputrec->grpopts.nrdf", i
, opt1
->nrdf
[i
], opt2
->nrdf
[i
], ftol
, abstol
);
595 cmp_real(fp
, "inputrec->grpopts.ref_t", i
, opt1
->ref_t
[i
], opt2
->ref_t
[i
], ftol
, abstol
);
596 cmp_real(fp
, "inputrec->grpopts.tau_t", i
, opt1
->tau_t
[i
], opt2
->tau_t
[i
], ftol
, abstol
);
597 cmp_int(fp
, "inputrec->grpopts.annealing", i
, opt1
->annealing
[i
], opt2
->annealing
[i
]);
598 cmp_int(fp
, "inputrec->grpopts.anneal_npoints", i
,
599 opt1
->anneal_npoints
[i
], opt2
->anneal_npoints
[i
]);
600 if (opt1
->anneal_npoints
[i
] == opt2
->anneal_npoints
[i
])
602 sprintf(buf1
, "inputrec->grpopts.anneal_time[%d]", i
);
603 sprintf(buf2
, "inputrec->grpopts.anneal_temp[%d]", i
);
604 for (j
= 0; j
< opt1
->anneal_npoints
[i
]; j
++)
606 cmp_real(fp
, buf1
, j
, opt1
->anneal_time
[i
][j
], opt2
->anneal_time
[i
][j
], ftol
, abstol
);
607 cmp_real(fp
, buf2
, j
, opt1
->anneal_temp
[i
][j
], opt2
->anneal_temp
[i
][j
], ftol
, abstol
);
611 if (opt1
->ngener
== opt2
->ngener
)
613 for (i
= 0; i
< opt1
->ngener
; i
++)
615 for (j
= i
; j
< opt1
->ngener
; j
++)
617 sprintf(buf1
, "inputrec->grpopts.egp_flags[%d]", i
);
619 opt1
->egp_flags
[opt1
->ngener
*i
+j
],
620 opt2
->egp_flags
[opt1
->ngener
*i
+j
]);
624 for (i
= 0; (i
< min(opt1
->ngacc
, opt2
->ngacc
)); i
++)
626 cmp_rvec(fp
, "inputrec->grpopts.acc", i
, opt1
->acc
[i
], opt2
->acc
[i
], ftol
, abstol
);
628 for (i
= 0; (i
< min(opt1
->ngfrz
, opt2
->ngfrz
)); i
++)
630 cmp_ivec(fp
, "inputrec->grpopts.nFreeze", i
, opt1
->nFreeze
[i
], opt2
->nFreeze
[i
]);
634 static void cmp_cosines(FILE *fp
, const char *s
, t_cosines c1
[DIM
], t_cosines c2
[DIM
], real ftol
, real abstol
)
639 for (m
= 0; (m
< DIM
); m
++)
641 sprintf(buf
, "inputrec->%s[%d]", s
, m
);
642 cmp_int(fp
, buf
, 0, c1
->n
, c2
->n
);
643 for (i
= 0; (i
< min(c1
->n
, c2
->n
)); i
++)
645 cmp_real(fp
, buf
, i
, c1
->a
[i
], c2
->a
[i
], ftol
, abstol
);
646 cmp_real(fp
, buf
, i
, c1
->phi
[i
], c2
->phi
[i
], ftol
, abstol
);
650 static void cmp_adress(FILE *fp
, t_adress
*ad1
, t_adress
*ad2
,
651 real ftol
, real abstol
)
653 cmp_int(fp
, "ir->adress->type", -1, ad1
->type
, ad2
->type
);
654 cmp_real(fp
, "ir->adress->const_wf", -1, ad1
->const_wf
, ad2
->const_wf
, ftol
, abstol
);
655 cmp_real(fp
, "ir->adress->ex_width", -1, ad1
->ex_width
, ad2
->ex_width
, ftol
, abstol
);
656 cmp_real(fp
, "ir->adress->hy_width", -1, ad1
->hy_width
, ad2
->hy_width
, ftol
, abstol
);
657 cmp_int(fp
, "ir->adress->icor", -1, ad1
->icor
, ad2
->icor
);
658 cmp_int(fp
, "ir->adress->site", -1, ad1
->site
, ad2
->site
);
659 cmp_rvec(fp
, "ir->adress->refs", -1, ad1
->refs
, ad2
->refs
, ftol
, abstol
);
660 cmp_real(fp
, "ir->adress->ex_forcecap", -1, ad1
->ex_forcecap
, ad2
->ex_forcecap
, ftol
, abstol
);
663 static void cmp_pull(FILE *fp
)
665 fprintf(fp
, "WARNING: Both files use COM pulling, but comparing of the pull struct is not implemented (yet). The pull parameters could be the same or different.\n");
668 static void cmp_simtempvals(FILE *fp
, t_simtemp
*simtemp1
, t_simtemp
*simtemp2
, int n_lambda
, real ftol
, real abstol
)
671 cmp_int(fp
, "inputrec->simtempvals->eSimTempScale", -1, simtemp1
->eSimTempScale
, simtemp2
->eSimTempScale
);
672 cmp_real(fp
, "inputrec->simtempvals->simtemp_high", -1, simtemp1
->simtemp_high
, simtemp2
->simtemp_high
, ftol
, abstol
);
673 cmp_real(fp
, "inputrec->simtempvals->simtemp_low", -1, simtemp1
->simtemp_low
, simtemp2
->simtemp_low
, ftol
, abstol
);
674 for (i
= 0; i
< n_lambda
; i
++)
676 cmp_real(fp
, "inputrec->simtempvals->temperatures", -1, simtemp1
->temperatures
[i
], simtemp2
->temperatures
[i
], ftol
, abstol
);
680 static void cmp_expandedvals(FILE *fp
, t_expanded
*expand1
, t_expanded
*expand2
, int n_lambda
, real ftol
, real abstol
)
684 cmp_bool(fp
, "inputrec->fepvals->bInit_weights", -1, expand1
->bInit_weights
, expand2
->bInit_weights
);
685 cmp_bool(fp
, "inputrec->fepvals->bWLoneovert", -1, expand1
->bWLoneovert
, expand2
->bWLoneovert
);
687 for (i
= 0; i
< n_lambda
; i
++)
689 cmp_real(fp
, "inputrec->expandedvals->init_lambda_weights", -1,
690 expand1
->init_lambda_weights
[i
], expand2
->init_lambda_weights
[i
], ftol
, abstol
);
693 cmp_int(fp
, "inputrec->expandedvals->lambda-stats", -1, expand1
->elamstats
, expand2
->elamstats
);
694 cmp_int(fp
, "inputrec->expandedvals->lambda-mc-move", -1, expand1
->elmcmove
, expand2
->elmcmove
);
695 cmp_int(fp
, "inputrec->expandedvals->lmc-repeats", -1, expand1
->lmc_repeats
, expand2
->lmc_repeats
);
696 cmp_int(fp
, "inputrec->expandedvals->lmc-gibbsdelta", -1, expand1
->gibbsdeltalam
, expand2
->gibbsdeltalam
);
697 cmp_int(fp
, "inputrec->expandedvals->lmc-forced-nstart", -1, expand1
->lmc_forced_nstart
, expand2
->lmc_forced_nstart
);
698 cmp_int(fp
, "inputrec->expandedvals->lambda-weights-equil", -1, expand1
->elmceq
, expand2
->elmceq
);
699 cmp_int(fp
, "inputrec->expandedvals->,weight-equil-number-all-lambda", -1, expand1
->equil_n_at_lam
, expand2
->equil_n_at_lam
);
700 cmp_int(fp
, "inputrec->expandedvals->weight-equil-number-samples", -1, expand1
->equil_samples
, expand2
->equil_samples
);
701 cmp_int(fp
, "inputrec->expandedvals->weight-equil-number-steps", -1, expand1
->equil_steps
, expand2
->equil_steps
);
702 cmp_real(fp
, "inputrec->expandedvals->weight-equil-wl-delta", -1, expand1
->equil_wl_delta
, expand2
->equil_wl_delta
, ftol
, abstol
);
703 cmp_real(fp
, "inputrec->expandedvals->weight-equil-count-ratio", -1, expand1
->equil_ratio
, expand2
->equil_ratio
, ftol
, abstol
);
704 cmp_bool(fp
, "inputrec->expandedvals->symmetrized-transition-matrix", -1, expand1
->bSymmetrizedTMatrix
, expand2
->bSymmetrizedTMatrix
);
705 cmp_int(fp
, "inputrec->expandedvals->nstTij", -1, expand1
->nstTij
, expand2
->nstTij
);
706 cmp_int(fp
, "inputrec->expandedvals->mininum-var-min", -1, expand1
->minvarmin
, expand2
->minvarmin
); /*default is reasonable */
707 cmp_int(fp
, "inputrec->expandedvals->weight-c-range", -1, expand1
->c_range
, expand2
->c_range
); /* default is just C=0 */
708 cmp_real(fp
, "inputrec->expandedvals->wl-scale", -1, expand1
->wl_scale
, expand2
->wl_scale
, ftol
, abstol
);
709 cmp_real(fp
, "inputrec->expandedvals->init-wl-delta", -1, expand1
->init_wl_delta
, expand2
->init_wl_delta
, ftol
, abstol
);
710 cmp_real(fp
, "inputrec->expandedvals->wl-ratio", -1, expand1
->wl_ratio
, expand2
->wl_ratio
, ftol
, abstol
);
711 cmp_int(fp
, "inputrec->expandedvals->nstexpanded", -1, expand1
->nstexpanded
, expand2
->nstexpanded
);
712 cmp_int(fp
, "inputrec->expandedvals->lmc-seed", -1, expand1
->lmc_seed
, expand2
->lmc_seed
);
713 cmp_real(fp
, "inputrec->expandedvals->mc-temperature", -1, expand1
->mc_temp
, expand2
->mc_temp
, ftol
, abstol
);
716 static void cmp_fepvals(FILE *fp
, t_lambda
*fep1
, t_lambda
*fep2
, real ftol
, real abstol
)
719 cmp_int(fp
, "inputrec->nstdhdl", -1, fep1
->nstdhdl
, fep2
->nstdhdl
);
720 cmp_double(fp
, "inputrec->fepvals->init_fep_state", -1, fep1
->init_fep_state
, fep2
->init_fep_state
, ftol
, abstol
);
721 cmp_double(fp
, "inputrec->fepvals->delta_lambda", -1, fep1
->delta_lambda
, fep2
->delta_lambda
, ftol
, abstol
);
722 cmp_int(fp
, "inputrec->fepvals->n_lambda", -1, fep1
->n_lambda
, fep2
->n_lambda
);
723 for (i
= 0; i
< efptNR
; i
++)
725 for (j
= 0; j
< min(fep1
->n_lambda
, fep2
->n_lambda
); j
++)
727 cmp_double(fp
, "inputrec->fepvals->all_lambda", -1, fep1
->all_lambda
[i
][j
], fep2
->all_lambda
[i
][j
], ftol
, abstol
);
730 cmp_int(fp
, "inputrec->fepvals->lambda_neighbors", 1, fep1
->lambda_neighbors
,
731 fep2
->lambda_neighbors
);
732 cmp_real(fp
, "inputrec->fepvals->sc_alpha", -1, fep1
->sc_alpha
, fep2
->sc_alpha
, ftol
, abstol
);
733 cmp_int(fp
, "inputrec->fepvals->sc_power", -1, fep1
->sc_power
, fep2
->sc_power
);
734 cmp_real(fp
, "inputrec->fepvals->sc_r_power", -1, fep1
->sc_r_power
, fep2
->sc_r_power
, ftol
, abstol
);
735 cmp_real(fp
, "inputrec->fepvals->sc_sigma", -1, fep1
->sc_sigma
, fep2
->sc_sigma
, ftol
, abstol
);
736 cmp_int(fp
, "inputrec->fepvals->edHdLPrintEnergy", -1, fep1
->edHdLPrintEnergy
, fep1
->edHdLPrintEnergy
);
737 cmp_bool(fp
, "inputrec->fepvals->bScCoul", -1, fep1
->bScCoul
, fep1
->bScCoul
);
738 cmp_int(fp
, "inputrec->separate_dhdl_file", -1, fep1
->separate_dhdl_file
, fep2
->separate_dhdl_file
);
739 cmp_int(fp
, "inputrec->dhdl_derivatives", -1, fep1
->dhdl_derivatives
, fep2
->dhdl_derivatives
);
740 cmp_int(fp
, "inputrec->dh_hist_size", -1, fep1
->dh_hist_size
, fep2
->dh_hist_size
);
741 cmp_double(fp
, "inputrec->dh_hist_spacing", -1, fep1
->dh_hist_spacing
, fep2
->dh_hist_spacing
, ftol
, abstol
);
744 static void cmp_inputrec(FILE *fp
, t_inputrec
*ir1
, t_inputrec
*ir2
, real ftol
, real abstol
)
746 fprintf(fp
, "comparing inputrec\n");
748 /* gcc 2.96 doesnt like these defines at all, but issues a huge list
749 * of warnings. Maybe it will change in future versions, but for the
750 * moment I've spelled them out instead. /EL 000820
751 * #define CIB(s) cmp_int(fp,"inputrec->"#s,0,ir1->##s,ir2->##s)
752 * #define CII(s) cmp_int(fp,"inputrec->"#s,0,ir1->##s,ir2->##s)
753 * #define CIR(s) cmp_real(fp,"inputrec->"#s,0,ir1->##s,ir2->##s,ftol)
755 cmp_int(fp
, "inputrec->eI", -1, ir1
->eI
, ir2
->eI
);
756 cmp_int64(fp
, "inputrec->nsteps", ir1
->nsteps
, ir2
->nsteps
);
757 cmp_int64(fp
, "inputrec->init_step", ir1
->init_step
, ir2
->init_step
);
758 cmp_int(fp
, "inputrec->simulation_part", -1, ir1
->simulation_part
, ir2
->simulation_part
);
759 cmp_int(fp
, "inputrec->ePBC", -1, ir1
->ePBC
, ir2
->ePBC
);
760 cmp_int(fp
, "inputrec->bPeriodicMols", -1, ir1
->bPeriodicMols
, ir2
->bPeriodicMols
);
761 cmp_int(fp
, "inputrec->cutoff_scheme", -1, ir1
->cutoff_scheme
, ir2
->cutoff_scheme
);
762 cmp_int(fp
, "inputrec->ns_type", -1, ir1
->ns_type
, ir2
->ns_type
);
763 cmp_int(fp
, "inputrec->nstlist", -1, ir1
->nstlist
, ir2
->nstlist
);
764 cmp_int(fp
, "inputrec->nstcomm", -1, ir1
->nstcomm
, ir2
->nstcomm
);
765 cmp_int(fp
, "inputrec->comm_mode", -1, ir1
->comm_mode
, ir2
->comm_mode
);
766 cmp_int(fp
, "inputrec->nstlog", -1, ir1
->nstlog
, ir2
->nstlog
);
767 cmp_int(fp
, "inputrec->nstxout", -1, ir1
->nstxout
, ir2
->nstxout
);
768 cmp_int(fp
, "inputrec->nstvout", -1, ir1
->nstvout
, ir2
->nstvout
);
769 cmp_int(fp
, "inputrec->nstfout", -1, ir1
->nstfout
, ir2
->nstfout
);
770 cmp_int(fp
, "inputrec->nstcalcenergy", -1, ir1
->nstcalcenergy
, ir2
->nstcalcenergy
);
771 cmp_int(fp
, "inputrec->nstenergy", -1, ir1
->nstenergy
, ir2
->nstenergy
);
772 cmp_int(fp
, "inputrec->nstxout_compressed", -1, ir1
->nstxout_compressed
, ir2
->nstxout_compressed
);
773 cmp_double(fp
, "inputrec->init_t", -1, ir1
->init_t
, ir2
->init_t
, ftol
, abstol
);
774 cmp_double(fp
, "inputrec->delta_t", -1, ir1
->delta_t
, ir2
->delta_t
, ftol
, abstol
);
775 cmp_real(fp
, "inputrec->x_compression_precision", -1, ir1
->x_compression_precision
, ir2
->x_compression_precision
, ftol
, abstol
);
776 cmp_real(fp
, "inputrec->fourierspacing", -1, ir1
->fourier_spacing
, ir2
->fourier_spacing
, ftol
, abstol
);
777 cmp_int(fp
, "inputrec->nkx", -1, ir1
->nkx
, ir2
->nkx
);
778 cmp_int(fp
, "inputrec->nky", -1, ir1
->nky
, ir2
->nky
);
779 cmp_int(fp
, "inputrec->nkz", -1, ir1
->nkz
, ir2
->nkz
);
780 cmp_int(fp
, "inputrec->pme_order", -1, ir1
->pme_order
, ir2
->pme_order
);
781 cmp_real(fp
, "inputrec->ewald_rtol", -1, ir1
->ewald_rtol
, ir2
->ewald_rtol
, ftol
, abstol
);
782 cmp_int(fp
, "inputrec->ewald_geometry", -1, ir1
->ewald_geometry
, ir2
->ewald_geometry
);
783 cmp_real(fp
, "inputrec->epsilon_surface", -1, ir1
->epsilon_surface
, ir2
->epsilon_surface
, ftol
, abstol
);
784 cmp_int(fp
, "inputrec->bContinuation", -1, ir1
->bContinuation
, ir2
->bContinuation
);
785 cmp_int(fp
, "inputrec->bShakeSOR", -1, ir1
->bShakeSOR
, ir2
->bShakeSOR
);
786 cmp_int(fp
, "inputrec->etc", -1, ir1
->etc
, ir2
->etc
);
787 cmp_int(fp
, "inputrec->bPrintNHChains", -1, ir1
->bPrintNHChains
, ir2
->bPrintNHChains
);
788 cmp_int(fp
, "inputrec->epc", -1, ir1
->epc
, ir2
->epc
);
789 cmp_int(fp
, "inputrec->epct", -1, ir1
->epct
, ir2
->epct
);
790 cmp_real(fp
, "inputrec->tau_p", -1, ir1
->tau_p
, ir2
->tau_p
, ftol
, abstol
);
791 cmp_rvec(fp
, "inputrec->ref_p(x)", -1, ir1
->ref_p
[XX
], ir2
->ref_p
[XX
], ftol
, abstol
);
792 cmp_rvec(fp
, "inputrec->ref_p(y)", -1, ir1
->ref_p
[YY
], ir2
->ref_p
[YY
], ftol
, abstol
);
793 cmp_rvec(fp
, "inputrec->ref_p(z)", -1, ir1
->ref_p
[ZZ
], ir2
->ref_p
[ZZ
], ftol
, abstol
);
794 cmp_rvec(fp
, "inputrec->compress(x)", -1, ir1
->compress
[XX
], ir2
->compress
[XX
], ftol
, abstol
);
795 cmp_rvec(fp
, "inputrec->compress(y)", -1, ir1
->compress
[YY
], ir2
->compress
[YY
], ftol
, abstol
);
796 cmp_rvec(fp
, "inputrec->compress(z)", -1, ir1
->compress
[ZZ
], ir2
->compress
[ZZ
], ftol
, abstol
);
797 cmp_int(fp
, "refcoord_scaling", -1, ir1
->refcoord_scaling
, ir2
->refcoord_scaling
);
798 cmp_rvec(fp
, "inputrec->posres_com", -1, ir1
->posres_com
, ir2
->posres_com
, ftol
, abstol
);
799 cmp_rvec(fp
, "inputrec->posres_comB", -1, ir1
->posres_comB
, ir2
->posres_comB
, ftol
, abstol
);
800 cmp_real(fp
, "inputrec->verletbuf_tol", -1, ir1
->verletbuf_tol
, ir2
->verletbuf_tol
, ftol
, abstol
);
801 cmp_real(fp
, "inputrec->rlist", -1, ir1
->rlist
, ir2
->rlist
, ftol
, abstol
);
802 cmp_real(fp
, "inputrec->rlistlong", -1, ir1
->rlistlong
, ir2
->rlistlong
, ftol
, abstol
);
803 cmp_int(fp
, "inputrec->nstcalclr", -1, ir1
->nstcalclr
, ir2
->nstcalclr
);
804 cmp_real(fp
, "inputrec->rtpi", -1, ir1
->rtpi
, ir2
->rtpi
, ftol
, abstol
);
805 cmp_int(fp
, "inputrec->coulombtype", -1, ir1
->coulombtype
, ir2
->coulombtype
);
806 cmp_int(fp
, "inputrec->coulomb_modifier", -1, ir1
->coulomb_modifier
, ir2
->coulomb_modifier
);
807 cmp_real(fp
, "inputrec->rcoulomb_switch", -1, ir1
->rcoulomb_switch
, ir2
->rcoulomb_switch
, ftol
, abstol
);
808 cmp_real(fp
, "inputrec->rcoulomb", -1, ir1
->rcoulomb
, ir2
->rcoulomb
, ftol
, abstol
);
809 cmp_int(fp
, "inputrec->vdwtype", -1, ir1
->vdwtype
, ir2
->vdwtype
);
810 cmp_int(fp
, "inputrec->vdw_modifier", -1, ir1
->vdw_modifier
, ir2
->vdw_modifier
); cmp_real(fp
, "inputrec->rvdw_switch", -1, ir1
->rvdw_switch
, ir2
->rvdw_switch
, ftol
, abstol
);
811 cmp_real(fp
, "inputrec->rvdw", -1, ir1
->rvdw
, ir2
->rvdw
, ftol
, abstol
);
812 cmp_real(fp
, "inputrec->epsilon_r", -1, ir1
->epsilon_r
, ir2
->epsilon_r
, ftol
, abstol
);
813 cmp_real(fp
, "inputrec->epsilon_rf", -1, ir1
->epsilon_rf
, ir2
->epsilon_rf
, ftol
, abstol
);
814 cmp_real(fp
, "inputrec->tabext", -1, ir1
->tabext
, ir2
->tabext
, ftol
, abstol
);
815 cmp_int(fp
, "inputrec->implicit_solvent", -1, ir1
->implicit_solvent
, ir2
->implicit_solvent
);
816 cmp_int(fp
, "inputrec->gb_algorithm", -1, ir1
->gb_algorithm
, ir2
->gb_algorithm
);
817 cmp_int(fp
, "inputrec->nstgbradii", -1, ir1
->nstgbradii
, ir2
->nstgbradii
);
818 cmp_real(fp
, "inputrec->rgbradii", -1, ir1
->rgbradii
, ir2
->rgbradii
, ftol
, abstol
);
819 cmp_real(fp
, "inputrec->gb_saltconc", -1, ir1
->gb_saltconc
, ir2
->gb_saltconc
, ftol
, abstol
);
820 cmp_real(fp
, "inputrec->gb_epsilon_solvent", -1, ir1
->gb_epsilon_solvent
, ir2
->gb_epsilon_solvent
, ftol
, abstol
);
821 cmp_real(fp
, "inputrec->gb_obc_alpha", -1, ir1
->gb_obc_alpha
, ir2
->gb_obc_alpha
, ftol
, abstol
);
822 cmp_real(fp
, "inputrec->gb_obc_beta", -1, ir1
->gb_obc_beta
, ir2
->gb_obc_beta
, ftol
, abstol
);
823 cmp_real(fp
, "inputrec->gb_obc_gamma", -1, ir1
->gb_obc_gamma
, ir2
->gb_obc_gamma
, ftol
, abstol
);
824 cmp_real(fp
, "inputrec->gb_dielectric_offset", -1, ir1
->gb_dielectric_offset
, ir2
->gb_dielectric_offset
, ftol
, abstol
);
825 cmp_int(fp
, "inputrec->sa_algorithm", -1, ir1
->sa_algorithm
, ir2
->sa_algorithm
);
826 cmp_real(fp
, "inputrec->sa_surface_tension", -1, ir1
->sa_surface_tension
, ir2
->sa_surface_tension
, ftol
, abstol
);
828 cmp_int(fp
, "inputrec->eDispCorr", -1, ir1
->eDispCorr
, ir2
->eDispCorr
);
829 cmp_real(fp
, "inputrec->shake_tol", -1, ir1
->shake_tol
, ir2
->shake_tol
, ftol
, abstol
);
830 cmp_int(fp
, "inputrec->efep", -1, ir1
->efep
, ir2
->efep
);
831 cmp_fepvals(fp
, ir1
->fepvals
, ir2
->fepvals
, ftol
, abstol
);
832 cmp_int(fp
, "inputrec->bSimTemp", -1, ir1
->bSimTemp
, ir2
->bSimTemp
);
833 if ((ir1
->bSimTemp
== ir2
->bSimTemp
) && (ir1
->bSimTemp
))
835 cmp_simtempvals(fp
, ir1
->simtempvals
, ir2
->simtempvals
, min(ir1
->fepvals
->n_lambda
, ir2
->fepvals
->n_lambda
), ftol
, abstol
);
837 cmp_int(fp
, "inputrec->bExpanded", -1, ir1
->bExpanded
, ir2
->bExpanded
);
838 if ((ir1
->bExpanded
== ir2
->bExpanded
) && (ir1
->bExpanded
))
840 cmp_expandedvals(fp
, ir1
->expandedvals
, ir2
->expandedvals
, min(ir1
->fepvals
->n_lambda
, ir2
->fepvals
->n_lambda
), ftol
, abstol
);
842 cmp_int(fp
, "inputrec->nwall", -1, ir1
->nwall
, ir2
->nwall
);
843 cmp_int(fp
, "inputrec->wall_type", -1, ir1
->wall_type
, ir2
->wall_type
);
844 cmp_int(fp
, "inputrec->wall_atomtype[0]", -1, ir1
->wall_atomtype
[0], ir2
->wall_atomtype
[0]);
845 cmp_int(fp
, "inputrec->wall_atomtype[1]", -1, ir1
->wall_atomtype
[1], ir2
->wall_atomtype
[1]);
846 cmp_real(fp
, "inputrec->wall_density[0]", -1, ir1
->wall_density
[0], ir2
->wall_density
[0], ftol
, abstol
);
847 cmp_real(fp
, "inputrec->wall_density[1]", -1, ir1
->wall_density
[1], ir2
->wall_density
[1], ftol
, abstol
);
848 cmp_real(fp
, "inputrec->wall_ewald_zfac", -1, ir1
->wall_ewald_zfac
, ir2
->wall_ewald_zfac
, ftol
, abstol
);
850 cmp_bool(fp
, "inputrec->bPull", -1, ir1
->bPull
, ir2
->bPull
);
851 if (ir1
->bPull
&& ir2
->bPull
)
856 cmp_int(fp
, "inputrec->eDisre", -1, ir1
->eDisre
, ir2
->eDisre
);
857 cmp_real(fp
, "inputrec->dr_fc", -1, ir1
->dr_fc
, ir2
->dr_fc
, ftol
, abstol
);
858 cmp_int(fp
, "inputrec->eDisreWeighting", -1, ir1
->eDisreWeighting
, ir2
->eDisreWeighting
);
859 cmp_int(fp
, "inputrec->bDisreMixed", -1, ir1
->bDisreMixed
, ir2
->bDisreMixed
);
860 cmp_int(fp
, "inputrec->nstdisreout", -1, ir1
->nstdisreout
, ir2
->nstdisreout
);
861 cmp_real(fp
, "inputrec->dr_tau", -1, ir1
->dr_tau
, ir2
->dr_tau
, ftol
, abstol
);
862 cmp_real(fp
, "inputrec->orires_fc", -1, ir1
->orires_fc
, ir2
->orires_fc
, ftol
, abstol
);
863 cmp_real(fp
, "inputrec->orires_tau", -1, ir1
->orires_tau
, ir2
->orires_tau
, ftol
, abstol
);
864 cmp_int(fp
, "inputrec->nstorireout", -1, ir1
->nstorireout
, ir2
->nstorireout
);
865 cmp_real(fp
, "inputrec->em_stepsize", -1, ir1
->em_stepsize
, ir2
->em_stepsize
, ftol
, abstol
);
866 cmp_real(fp
, "inputrec->em_tol", -1, ir1
->em_tol
, ir2
->em_tol
, ftol
, abstol
);
867 cmp_int(fp
, "inputrec->niter", -1, ir1
->niter
, ir2
->niter
);
868 cmp_real(fp
, "inputrec->fc_stepsize", -1, ir1
->fc_stepsize
, ir2
->fc_stepsize
, ftol
, abstol
);
869 cmp_int(fp
, "inputrec->nstcgsteep", -1, ir1
->nstcgsteep
, ir2
->nstcgsteep
);
870 cmp_int(fp
, "inputrec->nbfgscorr", 0, ir1
->nbfgscorr
, ir2
->nbfgscorr
);
871 cmp_int(fp
, "inputrec->eConstrAlg", -1, ir1
->eConstrAlg
, ir2
->eConstrAlg
);
872 cmp_int(fp
, "inputrec->nProjOrder", -1, ir1
->nProjOrder
, ir2
->nProjOrder
);
873 cmp_real(fp
, "inputrec->LincsWarnAngle", -1, ir1
->LincsWarnAngle
, ir2
->LincsWarnAngle
, ftol
, abstol
);
874 cmp_int(fp
, "inputrec->nLincsIter", -1, ir1
->nLincsIter
, ir2
->nLincsIter
);
875 cmp_real(fp
, "inputrec->bd_fric", -1, ir1
->bd_fric
, ir2
->bd_fric
, ftol
, abstol
);
876 cmp_int64(fp
, "inputrec->ld_seed", ir1
->ld_seed
, ir2
->ld_seed
);
877 cmp_real(fp
, "inputrec->cos_accel", -1, ir1
->cos_accel
, ir2
->cos_accel
, ftol
, abstol
);
878 cmp_rvec(fp
, "inputrec->deform(a)", -1, ir1
->deform
[XX
], ir2
->deform
[XX
], ftol
, abstol
);
879 cmp_rvec(fp
, "inputrec->deform(b)", -1, ir1
->deform
[YY
], ir2
->deform
[YY
], ftol
, abstol
);
880 cmp_rvec(fp
, "inputrec->deform(c)", -1, ir1
->deform
[ZZ
], ir2
->deform
[ZZ
], ftol
, abstol
);
883 cmp_bool(fp
, "ir->bAdress->type", -1, ir1
->bAdress
, ir2
->bAdress
);
884 if (ir1
->bAdress
&& ir2
->bAdress
)
886 cmp_adress(fp
, ir1
->adress
, ir2
->adress
, ftol
, abstol
);
889 cmp_int(fp
, "inputrec->userint1", -1, ir1
->userint1
, ir2
->userint1
);
890 cmp_int(fp
, "inputrec->userint2", -1, ir1
->userint2
, ir2
->userint2
);
891 cmp_int(fp
, "inputrec->userint3", -1, ir1
->userint3
, ir2
->userint3
);
892 cmp_int(fp
, "inputrec->userint4", -1, ir1
->userint4
, ir2
->userint4
);
893 cmp_real(fp
, "inputrec->userreal1", -1, ir1
->userreal1
, ir2
->userreal1
, ftol
, abstol
);
894 cmp_real(fp
, "inputrec->userreal2", -1, ir1
->userreal2
, ir2
->userreal2
, ftol
, abstol
);
895 cmp_real(fp
, "inputrec->userreal3", -1, ir1
->userreal3
, ir2
->userreal3
, ftol
, abstol
);
896 cmp_real(fp
, "inputrec->userreal4", -1, ir1
->userreal4
, ir2
->userreal4
, ftol
, abstol
);
897 cmp_grpopts(fp
, &(ir1
->opts
), &(ir2
->opts
), ftol
, abstol
);
898 cmp_cosines(fp
, "ex", ir1
->ex
, ir2
->ex
, ftol
, abstol
);
899 cmp_cosines(fp
, "et", ir1
->et
, ir2
->et
, ftol
, abstol
);
902 static void comp_pull_AB(FILE *fp
, pull_params_t
*pull
, real ftol
, real abstol
)
906 for (i
= 0; i
< pull
->ncoord
; i
++)
908 fprintf(fp
, "comparing pull coord %d\n", i
);
909 cmp_real(fp
, "pull-coord->k", -1, pull
->coord
[i
].k
, pull
->coord
[i
].kB
, ftol
, abstol
);
913 static void comp_state(t_state
*st1
, t_state
*st2
,
914 gmx_bool bRMSD
, real ftol
, real abstol
)
918 fprintf(stdout
, "comparing flags\n");
919 cmp_int(stdout
, "flags", -1, st1
->flags
, st2
->flags
);
920 fprintf(stdout
, "comparing box\n");
921 cmp_rvecs(stdout
, "box", DIM
, st1
->box
, st2
->box
, FALSE
, ftol
, abstol
);
922 fprintf(stdout
, "comparing box_rel\n");
923 cmp_rvecs(stdout
, "box_rel", DIM
, st1
->box_rel
, st2
->box_rel
, FALSE
, ftol
, abstol
);
924 fprintf(stdout
, "comparing boxv\n");
925 cmp_rvecs(stdout
, "boxv", DIM
, st1
->boxv
, st2
->boxv
, FALSE
, ftol
, abstol
);
926 if (st1
->flags
& (1<<estSVIR_PREV
))
928 fprintf(stdout
, "comparing shake vir_prev\n");
929 cmp_rvecs(stdout
, "svir_prev", DIM
, st1
->svir_prev
, st2
->svir_prev
, FALSE
, ftol
, abstol
);
931 if (st1
->flags
& (1<<estFVIR_PREV
))
933 fprintf(stdout
, "comparing force vir_prev\n");
934 cmp_rvecs(stdout
, "fvir_prev", DIM
, st1
->fvir_prev
, st2
->fvir_prev
, FALSE
, ftol
, abstol
);
936 if (st1
->flags
& (1<<estPRES_PREV
))
938 fprintf(stdout
, "comparing prev_pres\n");
939 cmp_rvecs(stdout
, "pres_prev", DIM
, st1
->pres_prev
, st2
->pres_prev
, FALSE
, ftol
, abstol
);
941 cmp_int(stdout
, "ngtc", -1, st1
->ngtc
, st2
->ngtc
);
942 cmp_int(stdout
, "nhchainlength", -1, st1
->nhchainlength
, st2
->nhchainlength
);
943 if (st1
->ngtc
== st2
->ngtc
&& st1
->nhchainlength
== st2
->nhchainlength
)
945 for (i
= 0; i
< st1
->ngtc
; i
++)
947 nc
= i
*st1
->nhchainlength
;
948 for (j
= 0; j
< nc
; j
++)
950 cmp_real(stdout
, "nosehoover_xi",
951 i
, st1
->nosehoover_xi
[nc
+j
], st2
->nosehoover_xi
[nc
+j
], ftol
, abstol
);
955 cmp_int(stdout
, "nnhpres", -1, st1
->nnhpres
, st2
->nnhpres
);
956 if (st1
->nnhpres
== st2
->nnhpres
&& st1
->nhchainlength
== st2
->nhchainlength
)
958 for (i
= 0; i
< st1
->nnhpres
; i
++)
960 nc
= i
*st1
->nhchainlength
;
961 for (j
= 0; j
< nc
; j
++)
963 cmp_real(stdout
, "nosehoover_xi",
964 i
, st1
->nhpres_xi
[nc
+j
], st2
->nhpres_xi
[nc
+j
], ftol
, abstol
);
969 cmp_int(stdout
, "natoms", -1, st1
->natoms
, st2
->natoms
);
970 if (st1
->natoms
== st2
->natoms
)
972 if ((st1
->flags
& (1<<estX
)) && (st2
->flags
& (1<<estX
)))
974 fprintf(stdout
, "comparing x\n");
975 cmp_rvecs(stdout
, "x", st1
->natoms
, st1
->x
, st2
->x
, bRMSD
, ftol
, abstol
);
977 if ((st1
->flags
& (1<<estV
)) && (st2
->flags
& (1<<estV
)))
979 fprintf(stdout
, "comparing v\n");
980 cmp_rvecs(stdout
, "v", st1
->natoms
, st1
->v
, st2
->v
, bRMSD
, ftol
, abstol
);
985 void comp_tpx(const char *fn1
, const char *fn2
,
986 gmx_bool bRMSD
, real ftol
, real abstol
)
998 for (i
= 0; i
< (fn2
? 2 : 1); i
++)
1000 read_tpx_state(ff
[i
], &(ir
[i
]), &state
[i
], NULL
, &(mtop
[i
]));
1004 cmp_inputrec(stdout
, &ir
[0], &ir
[1], ftol
, abstol
);
1005 /* Convert gmx_mtop_t to t_topology.
1006 * We should implement direct mtop comparison,
1007 * but it might be useful to keep t_topology comparison as an option.
1009 top
[0] = gmx_mtop_t_to_t_topology(&mtop
[0]);
1010 top
[1] = gmx_mtop_t_to_t_topology(&mtop
[1]);
1011 cmp_top(stdout
, &top
[0], &top
[1], ftol
, abstol
);
1012 cmp_groups(stdout
, &mtop
[0].groups
, &mtop
[1].groups
,
1013 mtop
[0].natoms
, mtop
[1].natoms
);
1014 comp_state(&state
[0], &state
[1], bRMSD
, ftol
, abstol
);
1018 if (ir
[0].efep
== efepNO
)
1020 fprintf(stdout
, "inputrec->efep = %s\n", efep_names
[ir
[0].efep
]);
1026 comp_pull_AB(stdout
, ir
->pull
, ftol
, abstol
);
1028 /* Convert gmx_mtop_t to t_topology.
1029 * We should implement direct mtop comparison,
1030 * but it might be useful to keep t_topology comparison as an option.
1032 top
[0] = gmx_mtop_t_to_t_topology(&mtop
[0]);
1033 cmp_top(stdout
, &top
[0], NULL
, ftol
, abstol
);
1038 void comp_frame(FILE *fp
, t_trxframe
*fr1
, t_trxframe
*fr2
,
1039 gmx_bool bRMSD
, real ftol
, real abstol
)
1042 cmp_int(fp
, "flags", -1, fr1
->flags
, fr2
->flags
);
1043 cmp_int(fp
, "not_ok", -1, fr1
->not_ok
, fr2
->not_ok
);
1044 cmp_int(fp
, "natoms", -1, fr1
->natoms
, fr2
->natoms
);
1045 cmp_real(fp
, "t0", -1, fr1
->t0
, fr2
->t0
, ftol
, abstol
);
1046 if (cmp_bool(fp
, "bTitle", -1, fr1
->bTitle
, fr2
->bTitle
))
1048 cmp_str(fp
, "title", -1, fr1
->title
, fr2
->title
);
1050 if (cmp_bool(fp
, "bStep", -1, fr1
->bStep
, fr2
->bStep
))
1052 cmp_int(fp
, "step", -1, fr1
->step
, fr2
->step
);
1054 cmp_int(fp
, "step", -1, fr1
->step
, fr2
->step
);
1055 if (cmp_bool(fp
, "bTime", -1, fr1
->bTime
, fr2
->bTime
))
1057 cmp_real(fp
, "time", -1, fr1
->time
, fr2
->time
, ftol
, abstol
);
1059 if (cmp_bool(fp
, "bLambda", -1, fr1
->bLambda
, fr2
->bLambda
))
1061 cmp_real(fp
, "lambda", -1, fr1
->lambda
, fr2
->lambda
, ftol
, abstol
);
1063 if (cmp_bool(fp
, "bAtoms", -1, fr1
->bAtoms
, fr2
->bAtoms
))
1065 cmp_atoms(fp
, fr1
->atoms
, fr2
->atoms
, ftol
, abstol
);
1067 if (cmp_bool(fp
, "bPrec", -1, fr1
->bPrec
, fr2
->bPrec
))
1069 cmp_real(fp
, "prec", -1, fr1
->prec
, fr2
->prec
, ftol
, abstol
);
1071 if (cmp_bool(fp
, "bX", -1, fr1
->bX
, fr2
->bX
))
1073 cmp_rvecs(fp
, "x", min(fr1
->natoms
, fr2
->natoms
), fr1
->x
, fr2
->x
, bRMSD
, ftol
, abstol
);
1075 if (cmp_bool(fp
, "bV", -1, fr1
->bV
, fr2
->bV
))
1077 cmp_rvecs(fp
, "v", min(fr1
->natoms
, fr2
->natoms
), fr1
->v
, fr2
->v
, bRMSD
, ftol
, abstol
);
1079 if (cmp_bool(fp
, "bF", -1, fr1
->bF
, fr2
->bF
))
1081 cmp_rvecs(fp
, "f", min(fr1
->natoms
, fr2
->natoms
), fr1
->f
, fr2
->f
, bRMSD
, ftol
, abstol
);
1083 if (cmp_bool(fp
, "bBox", -1, fr1
->bBox
, fr2
->bBox
))
1085 cmp_rvecs(fp
, "box", 3, fr1
->box
, fr2
->box
, FALSE
, ftol
, abstol
);
1089 void comp_trx(const output_env_t oenv
, const char *fn1
, const char *fn2
,
1090 gmx_bool bRMSD
, real ftol
, real abstol
)
1095 t_trxstatus
*status
[2];
1100 fprintf(stderr
, "Comparing trajectory files %s and %s\n", fn1
, fn2
);
1101 for (i
= 0; i
< 2; i
++)
1103 b
[i
] = read_first_frame(oenv
, &status
[i
], fn
[i
], &fr
[i
], TRX_READ_X
|TRX_READ_V
|TRX_READ_F
);
1110 comp_frame(stdout
, &(fr
[0]), &(fr
[1]), bRMSD
, ftol
, abstol
);
1112 for (i
= 0; i
< 2; i
++)
1114 b
[i
] = read_next_frame(oenv
, status
[i
], &fr
[i
]);
1117 while (b
[0] && b
[1]);
1119 for (i
= 0; i
< 2; i
++)
1121 if (b
[i
] && !b
[1-i
])
1123 fprintf(stdout
, "\nEnd of file on %s but not on %s\n", fn
[1-i
], fn
[i
]);
1125 close_trj(status
[i
]);
1130 fprintf(stdout
, "\nBoth files read correctly\n");
1134 static real
ener_tensor_diag(int n
, int *ind1
, int *ind2
,
1137 t_energy e1
[], t_energy e2
[])
1146 d2
= tensi
[i
] - d1
*DIM
;
1148 /* Find the diagonal elements d1 and d2 */
1149 len
= strlen(enm1
[ind1
[i
]].name
);
1153 for (j
= 0; j
< n
; j
++)
1155 if (tensi
[j
] >= 0 &&
1156 strlen(enm1
[ind1
[j
]].name
) == len
&&
1157 strncmp(enm1
[ind1
[i
]].name
, enm1
[ind1
[j
]].name
, len
-2) == 0 &&
1158 (tensi
[j
] == d1
*DIM
+d1
|| tensi
[j
] == d2
*DIM
+d2
))
1160 prod1
*= fabs(e1
[ind1
[j
]].e
);
1161 prod2
*= fabs(e2
[ind2
[j
]].e
);
1168 return 0.5*(sqrt(prod1
) + sqrt(prod2
));
1176 static gmx_bool
enernm_equal(const char *nm1
, const char *nm2
)
1183 /* Remove " (bar)" at the end of a name */
1184 if (len1
> 6 && strcmp(nm1
+len1
-6, " (bar)") == 0)
1188 if (len2
> 6 && strcmp(nm2
+len2
-6, " (bar)") == 0)
1193 return (len1
== len2
&& gmx_strncasecmp(nm1
, nm2
, len1
) == 0);
1196 static void cmp_energies(FILE *fp
, int step1
, int step2
,
1197 t_energy e1
[], t_energy e2
[],
1199 real ftol
, real abstol
,
1200 int nre
, int *ind1
, int *ind2
, int maxener
)
1203 int *tensi
, len
, d1
, d2
;
1204 real ftol_i
, abstol_i
;
1206 snew(tensi
, maxener
);
1207 /* Check for tensor elements ending on "-XX", "-XY", ... , "-ZZ" */
1208 for (i
= 0; (i
< maxener
); i
++)
1212 len
= strlen(enm1
[ii
].name
);
1213 if (len
> 3 && enm1
[ii
].name
[len
-3] == '-')
1215 d1
= enm1
[ii
].name
[len
-2] - 'X';
1216 d2
= enm1
[ii
].name
[len
-1] - 'X';
1217 if (d1
>= 0 && d1
< DIM
&&
1218 d2
>= 0 && d2
< DIM
)
1220 tensi
[i
] = d1
*DIM
+ d2
;
1225 for (i
= 0; (i
< maxener
); i
++)
1227 /* Check if this is an off-diagonal tensor element */
1228 if (tensi
[i
] >= 0 && tensi
[i
] != 0 && tensi
[i
] != 4 && tensi
[i
] != 8)
1230 /* Turn on the relative tolerance check (4 is maximum relative diff.) */
1232 /* Do the relative tolerance through an absolute tolerance times
1233 * the size of diagonal components of the tensor.
1235 abstol_i
= ftol
*ener_tensor_diag(nre
, ind1
, ind2
, enm1
, tensi
, i
, e1
, e2
);
1238 fprintf(debug
, "tensor '%s' val %f diag %f\n",
1239 enm1
[i
].name
, e1
[i
].e
, abstol_i
/ftol
);
1243 /* We found a diagonal, we need to check with the minimum tolerance */
1244 abstol_i
= min(abstol_i
, abstol
);
1248 /* We did not find a diagonal, ignore the relative tolerance check */
1257 if (!equal_real(e1
[ind1
[i
]].e
, e2
[ind2
[i
]].e
, ftol_i
, abstol_i
))
1259 fprintf(fp
, "%-15s step %3d: %12g, step %3d: %12g\n",
1261 step1
, e1
[ind1
[i
]].e
,
1262 step2
, e2
[ind2
[i
]].e
);
1270 static void cmp_disres(t_enxframe
*fr1
, t_enxframe
*fr2
, real ftol
, real abstol
)
1273 char bav
[64], bt
[64], bs
[22];
1275 cmp_int(stdout
, "ndisre", -1, fr1
->ndisre
, fr2
->ndisre
);
1276 if ((fr1
->ndisre
== fr2
->ndisre
) && (fr1
->ndisre
> 0))
1278 sprintf(bav
, "step %s: disre rav", gmx_step_str(fr1
->step
, bs
));
1279 sprintf(bt
, "step %s: disre rt", gmx_step_str(fr1
->step
, bs
));
1280 for (i
= 0; (i
< fr1
->ndisre
); i
++)
1282 cmp_real(stdout
, bav
, i
, fr1
->disre_rm3tav
[i
], fr2
->disre_rm3tav
[i
], ftol
, abstol
);
1283 cmp_real(stdout
, bt
, i
, fr1
->disre_rt
[i
], fr2
->disre_rt
[i
], ftol
, abstol
);
1289 static void cmp_eblocks(t_enxframe
*fr1
, t_enxframe
*fr2
, real ftol
, real abstol
)
1292 char buf
[64], bs
[22];
1294 cmp_int(stdout
, "nblock", -1, fr1
->nblock
, fr2
->nblock
);
1295 if ((fr1
->nblock
== fr2
->nblock
) && (fr1
->nblock
> 0))
1297 for (j
= 0; (j
< fr1
->nblock
); j
++)
1299 t_enxblock
*b1
, *b2
; /* convenience vars */
1301 b1
= &(fr1
->block
[j
]);
1302 b2
= &(fr2
->block
[j
]);
1304 sprintf(buf
, "step %s: block[%d]", gmx_step_str(fr1
->step
, bs
), j
);
1305 cmp_int(stdout
, buf
, -1, b1
->nsub
, b2
->nsub
);
1306 cmp_int(stdout
, buf
, -1, b1
->id
, b2
->id
);
1308 if ( (b1
->nsub
== b2
->nsub
) && (b1
->id
== b2
->id
) )
1310 for (i
= 0; i
< b1
->nsub
; i
++)
1312 t_enxsubblock
*s1
, *s2
;
1317 cmp_int(stdout
, buf
, -1, (int)s1
->type
, (int)s2
->type
);
1318 cmp_int64(stdout
, buf
, s1
->nr
, s2
->nr
);
1320 if ((s1
->type
== s2
->type
) && (s1
->nr
== s2
->nr
))
1324 case xdr_datatype_float
:
1325 for (k
= 0; k
< s1
->nr
; k
++)
1327 cmp_float(stdout
, buf
, i
,
1328 s1
->fval
[k
], s2
->fval
[k
],
1332 case xdr_datatype_double
:
1333 for (k
= 0; k
< s1
->nr
; k
++)
1335 cmp_double(stdout
, buf
, i
,
1336 s1
->dval
[k
], s2
->dval
[k
],
1340 case xdr_datatype_int
:
1341 for (k
= 0; k
< s1
->nr
; k
++)
1343 cmp_int(stdout
, buf
, i
,
1344 s1
->ival
[k
], s2
->ival
[k
]);
1347 case xdr_datatype_int64
:
1348 for (k
= 0; k
< s1
->nr
; k
++)
1350 cmp_int64(stdout
, buf
,
1351 s1
->lval
[k
], s2
->lval
[k
]);
1354 case xdr_datatype_char
:
1355 for (k
= 0; k
< s1
->nr
; k
++)
1357 cmp_uc(stdout
, buf
, i
,
1358 s1
->cval
[k
], s2
->cval
[k
]);
1361 case xdr_datatype_string
:
1362 for (k
= 0; k
< s1
->nr
; k
++)
1364 cmp_str(stdout
, buf
, i
,
1365 s1
->sval
[k
], s2
->sval
[k
]);
1369 gmx_incons("Unknown data type!!");
1378 void comp_enx(const char *fn1
, const char *fn2
, real ftol
, real abstol
, const char *lastener
)
1380 int nre
, nre1
, nre2
, block
;
1381 ener_file_t in1
, in2
;
1382 int i
, j
, maxener
, *ind1
, *ind2
, *have
;
1384 gmx_enxnm_t
*enm1
= NULL
, *enm2
= NULL
;
1385 t_enxframe
*fr1
, *fr2
;
1388 fprintf(stdout
, "comparing energy file %s and %s\n\n", fn1
, fn2
);
1390 in1
= open_enx(fn1
, "r");
1391 in2
= open_enx(fn2
, "r");
1392 do_enxnms(in1
, &nre1
, &enm1
);
1393 do_enxnms(in2
, &nre2
, &enm2
);
1396 fprintf(stdout
, "There are %d and %d terms in the energy files\n\n",
1401 fprintf(stdout
, "There are %d terms in the energy files\n\n", nre1
);
1408 for (i
= 0; i
< nre1
; i
++)
1410 for (j
= 0; j
< nre2
; j
++)
1412 if (enernm_equal(enm1
[i
].name
, enm2
[j
].name
))
1421 if (nre
== 0 || ind1
[nre
-1] != i
)
1423 cmp_str(stdout
, "enm", i
, enm1
[i
].name
, "-");
1426 for (i
= 0; i
< nre2
; i
++)
1430 cmp_str(stdout
, "enm", i
, "-", enm2
[i
].name
);
1435 for (i
= 0; i
< nre
; i
++)
1437 if ((lastener
!= NULL
) && (strstr(enm1
[i
].name
, lastener
) != NULL
))
1444 fprintf(stdout
, "There are %d terms to compare in the energy files\n\n",
1447 for (i
= 0; i
< maxener
; i
++)
1449 cmp_str(stdout
, "unit", i
, enm1
[ind1
[i
]].unit
, enm2
[ind2
[i
]].unit
);
1456 b1
= do_enx(in1
, fr1
);
1457 b2
= do_enx(in2
, fr2
);
1460 fprintf(stdout
, "\nEnd of file on %s but not on %s\n", fn2
, fn1
);
1464 fprintf(stdout
, "\nEnd of file on %s but not on %s\n", fn1
, fn2
);
1466 else if (!b1
&& !b2
)
1468 fprintf(stdout
, "\nFiles read successfully\n");
1472 cmp_real(stdout
, "t", -1, fr1
->t
, fr2
->t
, ftol
, abstol
);
1473 cmp_int(stdout
, "step", -1, fr1
->step
, fr2
->step
);
1474 /* We don't want to print the nre mismatch for every frame */
1475 /* cmp_int(stdout,"nre",-1,fr1->nre,fr2->nre); */
1476 if ((fr1
->nre
>= nre
) && (fr2
->nre
>= nre
))
1478 cmp_energies(stdout
, fr1
->step
, fr1
->step
, fr1
->ener
, fr2
->ener
,
1479 enm1
, ftol
, abstol
, nre
, ind1
, ind2
, maxener
);
1481 /*cmp_disres(fr1,fr2,ftol,abstol);*/
1482 cmp_eblocks(fr1
, fr2
, ftol
, abstol
);