From: Berk Hess Date: Fri, 20 Jun 2014 07:01:33 +0000 (+0200) Subject: Fix g_mindist with pbc=XY X-Git-Url: http://biod.pnpi.spb.ru/gitweb/?a=commitdiff_plain;h=b72bb55d47f6cb569aa35b7c3dfcab0c40a31d3d;p=alexxy%2Fgromacs.git Fix g_mindist with pbc=XY Fixes #1189. Change-Id: I7c8242df6f2644bebc49c13c86412f8bc1d37428 --- diff --git a/src/tools/gmx_mindist.c b/src/tools/gmx_mindist.c index 94718f2010..2116cf270e 100644 --- a/src/tools/gmx_mindist.c +++ b/src/tools/gmx_mindist.c @@ -58,20 +58,37 @@ #include "rmpbc.h" #include "xtcio.h" #include "gmx_ana.h" +#include "names.h" -static void periodic_dist(matrix box, rvec x[], int n, atom_id index[], +static void periodic_dist(int ePBC, + matrix box, rvec x[], int n, atom_id index[], real *rmin, real *rmax, int *min_ind) { -#define NSHIFT 26 - int sx, sy, sz, i, j, s; +#define NSHIFT_MAX 26 + int nsz, nshift, sx, sy, sz, i, j, s; real sqr_box, r2min, r2max, r2; - rvec shift[NSHIFT], d0, d; + rvec shift[NSHIFT_MAX], d0, d; - sqr_box = sqr(min(norm(box[XX]), min(norm(box[YY]), norm(box[ZZ])))); + sqr_box = min(norm2(box[XX]), norm2(box[YY])); + if (ePBC == epbcXYZ) + { + sqr_box = min(sqr_box, norm2(box[ZZ])); + nsz = 1; + } + else if (ePBC == epbcXY) + { + nsz = 0; + } + else + { + gmx_fatal(FARGS, "pbc = %s is not supported by g_mindist", + epbc_names[ePBC]); + nsz = 0; /* Keep compilers quiet */ + } - s = 0; - for (sz = -1; sz <= 1; sz++) + nshift = 0; + for (sz = -nsz; sz <= nsz; sz++) { for (sy = -1; sy <= 1; sy++) { @@ -81,9 +98,10 @@ static void periodic_dist(matrix box, rvec x[], int n, atom_id index[], { for (i = 0; i < DIM; i++) { - shift[s][i] = sx*box[XX][i]+sy*box[YY][i]+sz*box[ZZ][i]; + shift[nshift][i] = + sx*box[XX][i] + sy*box[YY][i] + sz*box[ZZ][i]; } - s++; + nshift++; } } } @@ -102,7 +120,7 @@ static void periodic_dist(matrix box, rvec x[], int n, atom_id index[], { r2max = r2; } - for (s = 0; s < NSHIFT; s++) + for (s = 0; s < nshift; s++) { rvec_add(d0, shift[s], d); r2 = norm2(d); @@ -164,7 +182,7 @@ static void periodic_mindist_plot(const char *trxfn, const char *outfn, gmx_rmpbc(gpbc, natoms, box, x); } - periodic_dist(box, x, n, index, &rmin, &rmax, ind_min); + periodic_dist(ePBC, box, x, n, index, &rmin, &rmax, ind_min); if (rmin < rmint) { rmint = rmin;