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-2013, The GROMACS development team.
6 * Copyright (c) 2013,2014,2015, 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/commandline/pargs.h"
44 #include "gromacs/fileio/confio.h"
45 #include "gromacs/fileio/enxio.h"
46 #include "gromacs/fileio/gmxfio.h"
47 #include "gromacs/fileio/tpxio.h"
48 #include "gromacs/fileio/trx.h"
49 #include "gromacs/fileio/trxio.h"
50 #include "gromacs/fileio/xtcio.h"
51 #include "gromacs/legacyheaders/names.h"
52 #include "gromacs/legacyheaders/txtdump.h"
53 #include "gromacs/legacyheaders/types/ifunc.h"
54 #include "gromacs/math/units.h"
55 #include "gromacs/math/vec.h"
56 #include "gromacs/pbcutil/pbc.h"
57 #include "gromacs/tools/compare.h"
58 #include "gromacs/topology/atomprop.h"
59 #include "gromacs/topology/block.h"
60 #include "gromacs/topology/index.h"
61 #include "gromacs/topology/mtop_util.h"
62 #include "gromacs/utility/arraysize.h"
63 #include "gromacs/utility/fatalerror.h"
64 #include "gromacs/utility/futil.h"
65 #include "gromacs/utility/smalloc.h"
87 static void tpx2system(FILE *fp
, gmx_mtop_t
*mtop
)
90 gmx_mtop_atomloop_block_t aloop
;
93 fprintf(fp
, "\\subsection{Simulation system}\n");
94 aloop
= gmx_mtop_atomloop_block_init(mtop
);
95 while (gmx_mtop_atomloop_block_next(aloop
, &atom
, &nmol
))
97 if (atom
->ptype
== eptVSite
)
102 fprintf(fp
, "A system of %d molecules (%d atoms) was simulated.\n",
103 mtop
->mols
.nr
, mtop
->natoms
-nvsite
);
106 fprintf(fp
, "Virtual sites were used in some of the molecules.\n");
111 static void tpx2params(FILE *fp
, t_inputrec
*ir
)
113 fprintf(fp
, "\\subsection{Simulation settings}\n");
114 fprintf(fp
, "A total of %g ns were simulated with a time step of %g fs.\n",
115 ir
->nsteps
*ir
->delta_t
*0.001, 1000*ir
->delta_t
);
116 fprintf(fp
, "Neighbor searching was performed every %d steps.\n", ir
->nstlist
);
117 fprintf(fp
, "The %s algorithm was used for electrostatic interactions.\n",
118 EELTYPE(ir
->coulombtype
));
119 fprintf(fp
, "with a cut-off of %g nm.\n", ir
->rcoulomb
);
120 if (ir
->coulombtype
== eelPME
)
122 fprintf(fp
, "A reciprocal grid of %d x %d x %d cells was used with %dth order B-spline interpolation.\n", ir
->nkx
, ir
->nky
, ir
->nkz
, ir
->pme_order
);
124 if (ir
->rvdw
> ir
->rlist
)
126 fprintf(fp
, "A twin-range Van der Waals cut-off (%g/%g nm) was used, where the long range forces were updated during neighborsearching.\n", ir
->rlist
, ir
->rvdw
);
130 fprintf(fp
, "A single cut-off of %g was used for Van der Waals interactions.\n", ir
->rlist
);
134 fprintf(fp
, "Temperature coupling was done with the %s algorithm.\n",
135 etcoupl_names
[ir
->etc
]);
139 fprintf(fp
, "Pressure coupling was done with the %s algorithm.\n",
140 epcoupl_names
[ir
->epc
]);
145 static void tpx2methods(const char *tpx
, const char *tex
)
152 read_tpx_state(tpx
, &ir
, &state
, &mtop
);
153 fp
= gmx_fio_fopen(tex
, "w");
154 fprintf(fp
, "\\section{Methods}\n");
155 tpx2system(fp
, &mtop
);
160 static void chk_coords(int frame
, int natoms
, rvec
*x
, matrix box
, real fac
, real tol
)
166 for (i
= 0; (i
< natoms
); i
++)
168 for (j
= 0; (j
< DIM
); j
++)
170 if ((vol
> 0) && (fabs(x
[i
][j
]) > fac
*box
[j
][j
]))
172 printf("Warning at frame %d: coordinates for atom %d are large (%g)\n",
176 if ((fabs(x
[i
][XX
]) < tol
) &&
177 (fabs(x
[i
][YY
]) < tol
) &&
178 (fabs(x
[i
][ZZ
]) < tol
))
185 printf("Warning at frame %d: there are %d particles with all coordinates zero\n",
190 static void chk_vels(int frame
, int natoms
, rvec
*v
)
194 for (i
= 0; (i
< natoms
); i
++)
196 for (j
= 0; (j
< DIM
); j
++)
198 if (fabs(v
[i
][j
]) > 500)
200 printf("Warning at frame %d. Velocities for atom %d are large (%g)\n",
207 static void chk_forces(int frame
, int natoms
, rvec
*f
)
211 for (i
= 0; (i
< natoms
); i
++)
213 for (j
= 0; (j
< DIM
); j
++)
215 if (fabs(f
[i
][j
]) > 10000)
217 printf("Warning at frame %d. Forces for atom %d are large (%g)\n",
224 static void chk_bonds(t_idef
*idef
, int ePBC
, rvec
*x
, matrix box
, real tol
)
226 int ftype
, k
, ai
, aj
, type
;
227 real b0
, blen
, deviation
;
231 set_pbc(&pbc
, ePBC
, box
);
232 for (ftype
= 0; (ftype
< F_NRE
); ftype
++)
234 if ((interaction_function
[ftype
].flags
& IF_CHEMBOND
) == IF_CHEMBOND
)
236 for (k
= 0; (k
< idef
->il
[ftype
].nr
); )
238 type
= idef
->il
[ftype
].iatoms
[k
++];
239 ai
= idef
->il
[ftype
].iatoms
[k
++];
240 aj
= idef
->il
[ftype
].iatoms
[k
++];
245 b0
= idef
->iparams
[type
].harmonic
.rA
;
248 b0
= std::sqrt(idef
->iparams
[type
].harmonic
.rA
);
251 b0
= idef
->iparams
[type
].morse
.b0A
;
254 b0
= idef
->iparams
[type
].cubic
.b0
;
257 b0
= idef
->iparams
[type
].constr
.dA
;
264 pbc_dx(&pbc
, x
[ai
], x
[aj
], dx
);
266 deviation
= sqr(blen
-b0
);
267 if (std::sqrt(deviation
/sqr(b0
)) > tol
)
269 fprintf(stderr
, "Distance between atoms %d and %d is %.3f, should be %.3f\n", ai
+1, aj
+1, blen
, b0
);
277 void chk_trj(const gmx_output_env_t
*oenv
, const char *fn
, const char *tpr
, real tol
)
281 t_fr_time first
, last
;
282 int j
= -1, new_natoms
, natoms
;
284 gmx_bool bShowTimestep
= TRUE
, newline
= FALSE
;
287 gmx_localtop_t
*top
= NULL
;
293 read_tpx_state(tpr
, &ir
, &state
, &mtop
);
294 top
= gmx_mtop_generate_local_top(&mtop
, &ir
);
299 printf("Checking file %s\n", fn
);
329 read_first_frame(oenv
, &status
, fn
, &fr
, TRX_READ_X
| TRX_READ_V
| TRX_READ_F
);
335 fprintf(stderr
, "\n# Atoms %d\n", fr
.natoms
);
338 fprintf(stderr
, "Precision %g (nm)\n", 1/fr
.prec
);
342 if ((natoms
> 0) && (new_natoms
!= natoms
))
344 fprintf(stderr
, "\nNumber of atoms at t=%g don't match (%d, %d)\n",
345 old_t1
, natoms
, new_natoms
);
350 if (fabs((fr
.time
-old_t1
)-(old_t1
-old_t2
)) >
351 0.1*(fabs(fr
.time
-old_t1
)+fabs(old_t1
-old_t2
)) )
353 bShowTimestep
= FALSE
;
354 fprintf(stderr
, "%sTimesteps at t=%g don't match (%g, %g)\n",
355 newline
? "\n" : "", old_t1
, old_t1
-old_t2
, fr
.time
-old_t1
);
361 chk_bonds(&top
->idef
, ir
.ePBC
, fr
.x
, fr
.box
, tol
);
365 chk_coords(j
, natoms
, fr
.x
, fr
.box
, 1e5
, tol
);
369 chk_vels(j
, natoms
, fr
.v
);
373 chk_forces(j
, natoms
, fr
.f
);
379 new_natoms
= fr
.natoms
;
380 #define INC(s, n, f, l, item) if (s.item != 0) { if (n.item == 0) { first.item = fr.time; } last.item = fr.time; n.item++; \
382 INC(fr
, count
, first
, last
, bStep
);
383 INC(fr
, count
, first
, last
, bTime
);
384 INC(fr
, count
, first
, last
, bLambda
);
385 INC(fr
, count
, first
, last
, bX
);
386 INC(fr
, count
, first
, last
, bV
);
387 INC(fr
, count
, first
, last
, bF
);
388 INC(fr
, count
, first
, last
, bBox
);
391 while (read_next_frame(oenv
, status
, &fr
));
393 fprintf(stderr
, "\n");
397 fprintf(stderr
, "\nItem #frames");
400 fprintf(stderr
, " Timestep (ps)");
402 fprintf(stderr
, "\n");
403 #define PRINTITEM(label, item) fprintf(stderr, "%-10s %6d", label, count.item); if ((bShowTimestep) && (count.item > 1)) {fprintf(stderr, " %g\n", (last.item-first.item)/(count.item-1)); }else fprintf(stderr, "\n")
404 PRINTITEM ( "Step", bStep
);
405 PRINTITEM ( "Time", bTime
);
406 PRINTITEM ( "Lambda", bLambda
);
407 PRINTITEM ( "Coords", bX
);
408 PRINTITEM ( "Velocities", bV
);
409 PRINTITEM ( "Forces", bF
);
410 PRINTITEM ( "Box", bBox
);
413 void chk_tps(const char *fn
, real vdw_fac
, real bon_lo
, real bon_hi
)
423 gmx_bool bV
, bX
, bB
, bFirst
, bOut
;
424 real r2
, ekin
, temp1
, temp2
, dist2
, vdwfac2
, bonlo2
, bonhi2
;
428 fprintf(stderr
, "Checking coordinate file %s\n", fn
);
429 read_tps_conf(fn
, &top
, &ePBC
, &x
, &v
, box
, TRUE
);
432 fprintf(stderr
, "%d atoms in file\n", atoms
->nr
);
434 /* check coordinates and box */
437 for (i
= 0; (i
< natom
) && !(bV
&& bX
); i
++)
439 for (j
= 0; (j
< DIM
) && !(bV
&& bX
); j
++)
441 bV
= bV
|| (v
[i
][j
] != 0);
442 bX
= bX
|| (x
[i
][j
] != 0);
446 for (i
= 0; (i
< DIM
) && !bB
; i
++)
448 for (j
= 0; (j
< DIM
) && !bB
; j
++)
450 bB
= bB
|| (box
[i
][j
] != 0);
454 fprintf(stderr
, "coordinates %s\n", bX
? "found" : "absent");
455 fprintf(stderr
, "box %s\n", bB
? "found" : "absent");
456 fprintf(stderr
, "velocities %s\n", bV
? "found" : "absent");
457 fprintf(stderr
, "\n");
459 /* check velocities */
463 for (i
= 0; (i
< natom
); i
++)
465 for (j
= 0; (j
< DIM
); j
++)
467 ekin
+= 0.5*atoms
->atom
[i
].m
*v
[i
][j
]*v
[i
][j
];
470 temp1
= (2.0*ekin
)/(natom
*DIM
*BOLTZ
);
471 temp2
= (2.0*ekin
)/(natom
*(DIM
-1)*BOLTZ
);
472 fprintf(stderr
, "Kinetic energy: %g (kJ/mol)\n", ekin
);
473 fprintf(stderr
, "Assuming the number of degrees of freedom to be "
474 "Natoms * %d or Natoms * %d,\n"
475 "the velocities correspond to a temperature of the system\n"
476 "of %g K or %g K respectively.\n\n", DIM
, DIM
-1, temp1
, temp2
);
479 /* check coordinates */
482 vdwfac2
= sqr(vdw_fac
);
483 bonlo2
= sqr(bon_lo
);
484 bonhi2
= sqr(bon_hi
);
487 "Checking for atoms closer than %g and not between %g and %g,\n"
488 "relative to sum of Van der Waals distance:\n",
489 vdw_fac
, bon_lo
, bon_hi
);
490 snew(atom_vdw
, natom
);
491 aps
= gmx_atomprop_init();
492 for (i
= 0; (i
< natom
); i
++)
494 gmx_atomprop_query(aps
, epropVDW
,
495 *(atoms
->resinfo
[atoms
->atom
[i
].resind
].name
),
496 *(atoms
->atomname
[i
]), &(atom_vdw
[i
]));
499 fprintf(debug
, "%5d %4s %4s %7g\n", i
+1,
500 *(atoms
->resinfo
[atoms
->atom
[i
].resind
].name
),
501 *(atoms
->atomname
[i
]),
505 gmx_atomprop_destroy(aps
);
508 set_pbc(&pbc
, ePBC
, box
);
512 for (i
= 0; (i
< natom
); i
++)
516 fprintf(stderr
, "\r%5d", i
+1);
518 for (j
= i
+1; (j
< natom
); j
++)
522 pbc_dx(&pbc
, x
[i
], x
[j
], dx
);
526 rvec_sub(x
[i
], x
[j
], dx
);
529 dist2
= sqr(atom_vdw
[i
]+atom_vdw
[j
]);
530 if ( (r2
<= dist2
*bonlo2
) ||
531 ( (r2
>= dist2
*bonhi2
) && (r2
<= dist2
*vdwfac2
) ) )
535 fprintf(stderr
, "\r%5s %4s %8s %5s %5s %4s %8s %5s %6s\n",
536 "atom#", "name", "residue", "r_vdw",
537 "atom#", "name", "residue", "r_vdw", "distance");
541 "\r%5d %4s %4s%4d %-5.3g %5d %4s %4s%4d %-5.3g %-6.4g\n",
542 i
+1, *(atoms
->atomname
[i
]),
543 *(atoms
->resinfo
[atoms
->atom
[i
].resind
].name
),
544 atoms
->resinfo
[atoms
->atom
[i
].resind
].nr
,
546 j
+1, *(atoms
->atomname
[j
]),
547 *(atoms
->resinfo
[atoms
->atom
[j
].resind
].name
),
548 atoms
->resinfo
[atoms
->atom
[j
].resind
].nr
,
556 fprintf(stderr
, "\rno close atoms found\n");
558 fprintf(stderr
, "\r \n");
565 for (i
= 0; (i
< natom
) && (k
< 10); i
++)
568 for (j
= 0; (j
< DIM
) && !bOut
; j
++)
570 bOut
= bOut
|| (x
[i
][j
] < 0) || (x
[i
][j
] > box
[j
][j
]);
577 fprintf(stderr
, "Atoms outside box ( ");
578 for (j
= 0; (j
< DIM
); j
++)
580 fprintf(stderr
, "%g ", box
[j
][j
]);
582 fprintf(stderr
, "):\n"
583 "(These may occur often and are normally not a problem)\n"
584 "%5s %4s %8s %5s %s\n",
585 "atom#", "name", "residue", "r_vdw", "coordinate");
589 "%5d %4s %4s%4d %-5.3g",
590 i
, *(atoms
->atomname
[i
]),
591 *(atoms
->resinfo
[atoms
->atom
[i
].resind
].name
),
592 atoms
->resinfo
[atoms
->atom
[i
].resind
].nr
, atom_vdw
[i
]);
593 for (j
= 0; (j
< DIM
); j
++)
595 fprintf(stderr
, " %6.3g", x
[i
][j
]);
597 fprintf(stderr
, "\n");
602 fprintf(stderr
, "(maybe more)\n");
606 fprintf(stderr
, "no atoms found outside box\n");
608 fprintf(stderr
, "\n");
613 void chk_ndx(const char *fn
)
619 grps
= init_index(fn
, &grpname
);
622 pr_blocka(debug
, 0, fn
, grps
, FALSE
);
626 printf("Contents of index file %s\n", fn
);
627 printf("--------------------------------------------------\n");
628 printf("Nr. Group #Entries First Last\n");
629 for (i
= 0; (i
< grps
->nr
); i
++)
631 printf("%4d %-20s%8d%8d%8d\n", i
, grpname
[i
],
632 grps
->index
[i
+1]-grps
->index
[i
],
633 grps
->a
[grps
->index
[i
]]+1,
634 grps
->a
[grps
->index
[i
+1]-1]+1);
637 for (i
= 0; (i
< grps
->nr
); i
++)
645 void chk_enx(const char *fn
)
649 gmx_enxnm_t
*enm
= NULL
;
653 real t0
, old_t1
, old_t2
;
656 fprintf(stderr
, "Checking energy file %s\n\n", fn
);
658 in
= open_enx(fn
, "r");
659 do_enxnms(in
, &nre
, &enm
);
660 fprintf(stderr
, "%d groups in energy file", nre
);
669 while (do_enx(in
, fr
))
673 if (fabs((fr
->t
-old_t1
)-(old_t1
-old_t2
)) >
674 0.1*(fabs(fr
->t
-old_t1
)+fabs(old_t1
-old_t2
)) )
677 fprintf(stderr
, "\nTimesteps at t=%g don't match (%g, %g)\n",
678 old_t1
, old_t1
-old_t2
, fr
->t
-old_t1
);
690 fprintf(stderr
, "\rframe: %6s (index %6d), t: %10.3f\n",
691 gmx_step_str(fr
->step
, buf
), fnr
, fr
->t
);
695 fprintf(stderr
, "\n\nFound %d frames", fnr
);
696 if (bShowTStep
&& fnr
> 1)
698 fprintf(stderr
, " with a timestep of %g ps", (old_t1
-t0
)/(fnr
-1));
700 fprintf(stderr
, ".\n");
703 free_enxnms(nre
, enm
);
707 int gmx_check(int argc
, char *argv
[])
709 const char *desc
[] = {
710 "[THISMODULE] reads a trajectory ([REF].tng[ref], [REF].trr[ref] or ",
711 "[REF].xtc[ref]), an energy file ([REF].edr[ref])",
712 "or an index file ([REF].ndx[ref])",
713 "and prints out useful information about them.[PAR]",
714 "Option [TT]-c[tt] checks for presence of coordinates,",
715 "velocities and box in the file, for close contacts (smaller than",
716 "[TT]-vdwfac[tt] and not bonded, i.e. not between [TT]-bonlo[tt]",
717 "and [TT]-bonhi[tt], all relative to the sum of both Van der Waals",
718 "radii) and atoms outside the box (these may occur often and are",
719 "no problem). If velocities are present, an estimated temperature",
720 "will be calculated from them.[PAR]",
721 "If an index file, is given its contents will be summarized.[PAR]",
722 "If both a trajectory and a [REF].tpr[ref] file are given (with [TT]-s1[tt])",
723 "the program will check whether the bond lengths defined in the tpr",
724 "file are indeed correct in the trajectory. If not you may have",
725 "non-matching files due to e.g. deshuffling or due to problems with",
726 "virtual sites. With these flags, [TT]gmx check[tt] provides a quick check for such problems.[PAR]",
727 "The program can compare two run input ([REF].tpr[ref])",
729 "when both [TT]-s1[tt] and [TT]-s2[tt] are supplied.",
730 "Similarly a pair of trajectory files can be compared (using the [TT]-f2[tt]",
731 "option), or a pair of energy files (using the [TT]-e2[tt] option).[PAR]",
732 "For free energy simulations the A and B state topology from one",
733 "run input file can be compared with options [TT]-s1[tt] and [TT]-ab[tt].[PAR]",
734 "In case the [TT]-m[tt] flag is given a LaTeX file will be written",
735 "consisting of a rough outline for a methods section for a paper."
738 { efTRX
, "-f", NULL
, ffOPTRD
},
739 { efTRX
, "-f2", NULL
, ffOPTRD
},
740 { efTPR
, "-s1", "top1", ffOPTRD
},
741 { efTPR
, "-s2", "top2", ffOPTRD
},
742 { efTPS
, "-c", NULL
, ffOPTRD
},
743 { efEDR
, "-e", NULL
, ffOPTRD
},
744 { efEDR
, "-e2", "ener2", ffOPTRD
},
745 { efNDX
, "-n", NULL
, ffOPTRD
},
746 { efTEX
, "-m", NULL
, ffOPTWR
}
748 #define NFILE asize(fnm)
749 const char *fn1
= NULL
, *fn2
= NULL
, *tex
= NULL
;
751 gmx_output_env_t
*oenv
;
752 static real vdw_fac
= 0.8;
753 static real bon_lo
= 0.4;
754 static real bon_hi
= 0.7;
755 static gmx_bool bRMSD
= FALSE
;
756 static real ftol
= 0.001;
757 static real abstol
= 0.001;
758 static gmx_bool bCompAB
= FALSE
;
759 static char *lastener
= NULL
;
760 static t_pargs pa
[] = {
761 { "-vdwfac", FALSE
, etREAL
, {&vdw_fac
},
762 "Fraction of sum of VdW radii used as warning cutoff" },
763 { "-bonlo", FALSE
, etREAL
, {&bon_lo
},
764 "Min. fract. of sum of VdW radii for bonded atoms" },
765 { "-bonhi", FALSE
, etREAL
, {&bon_hi
},
766 "Max. fract. of sum of VdW radii for bonded atoms" },
767 { "-rmsd", FALSE
, etBOOL
, {&bRMSD
},
768 "Print RMSD for x, v and f" },
769 { "-tol", FALSE
, etREAL
, {&ftol
},
770 "Relative tolerance for comparing real values defined as [MATH]2*(a-b)/([MAG]a[mag]+[MAG]b[mag])[math]" },
771 { "-abstol", FALSE
, etREAL
, {&abstol
},
772 "Absolute tolerance, useful when sums are close to zero." },
773 { "-ab", FALSE
, etBOOL
, {&bCompAB
},
774 "Compare the A and B topology from one file" },
775 { "-lastener", FALSE
, etSTR
, {&lastener
},
776 "Last energy term to compare (if not given all are tested). It makes sense to go up until the Pressure." }
779 if (!parse_common_args(&argc
, argv
, 0, NFILE
, fnm
, asize(pa
), pa
,
780 asize(desc
), desc
, 0, NULL
, &oenv
))
785 fn1
= opt2fn_null("-f", NFILE
, fnm
);
786 fn2
= opt2fn_null("-f2", NFILE
, fnm
);
787 tex
= opt2fn_null("-m", NFILE
, fnm
);
790 comp_trx(oenv
, fn1
, fn2
, bRMSD
, ftol
, abstol
);
794 chk_trj(oenv
, fn1
, opt2fn_null("-s1", NFILE
, fnm
), ftol
);
798 fprintf(stderr
, "Please give me TWO trajectory (.xtc/.trr/.tng) files!\n");
801 fn1
= opt2fn_null("-s1", NFILE
, fnm
);
802 fn2
= opt2fn_null("-s2", NFILE
, fnm
);
803 if ((fn1
&& fn2
) || bCompAB
)
809 gmx_fatal(FARGS
, "With -ab you need to set the -s1 option");
813 comp_tpx(fn1
, fn2
, bRMSD
, ftol
, abstol
);
817 tpx2methods(fn1
, tex
);
819 else if ((fn1
&& !opt2fn_null("-f", NFILE
, fnm
)) || (!fn1
&& fn2
))
821 fprintf(stderr
, "Please give me TWO run input (.tpr) files\n"
822 "or specify the -m flag to generate a methods.tex file\n");
825 fn1
= opt2fn_null("-e", NFILE
, fnm
);
826 fn2
= opt2fn_null("-e2", NFILE
, fnm
);
829 comp_enx(fn1
, fn2
, ftol
, abstol
, lastener
);
833 chk_enx(ftp2fn(efEDR
, NFILE
, fnm
));
837 fprintf(stderr
, "Please give me TWO energy (.edr/.ene) files!\n");
840 if (ftp2bSet(efTPS
, NFILE
, fnm
))
842 chk_tps(ftp2fn(efTPS
, NFILE
, fnm
), vdw_fac
, bon_lo
, bon_hi
);
845 if (ftp2bSet(efNDX
, NFILE
, fnm
))
847 chk_ndx(ftp2fn(efNDX
, NFILE
, fnm
));