\secref{mod_nb_int}). One then has a buffer with the size equal to the
neighbor list cut-off less the longest interaction cut-off.
-With the
-group cut-off scheme, one can also choose to let {\tt mdrun} only
-update the neighbor list when required. That is when one or more
-particles have moved more than half the buffer size from the center of
-geometry of the \swapindex{charge}{group} to which they belong (see
-\secref{chargegroup}), as determined at the previous neighbor search.
-This option guarantees that there are no cut-off artifacts. {\bf
- Note} that for larger systems this comes at a high computational
-cost, since the neighbor list update frequency will be determined by
-just one or two particles moving slightly beyond the half buffer
-length (which not even necessarily implies that the neighbor list is
-invalid), while 99.99\% of the particles are fine.
-
} % Brace matches ifthenelse test for gmxlite
\subsubsection{Simple search\swapindexquiet{simple}{search}}
<dd>The neighbor list is only constructed once and never updated.
This is mainly useful for vacuum simulations in which all particles
see each other.</dd>
-<dt><b>-1</b></dt>
-<dd>Automated update frequency, only supported with <b>cutoff-scheme</b>=<b>group</b>.
-This can only be used with switched, shifted or user potentials where
-the cut-off can be smaller than <b>rlist</b>. One then has a buffer
-of size <b>rlist</b> minus the longest cut-off.
-The neighbor list is only updated when one or more particles have moved further
-than half the buffer size from the center of geometry of their charge group
-as determined at the previous neighbor search.
-Coordinate scaling due to pressure coupling or the <b>deform</b> option
-is taken into account.
-This option guarantees that their are no cut-off artifacts,
-but for larger systems this can come at a high computational cost,
-since the neighbor list update frequency will be determined
-by just one or two particles moving slightly beyond the half buffer length
-(which does not necessarily imply that the neighbor list is invalid),
-while 99.99% of the particles are fine.
-</dd>
+<dt><b><0</b></dt>
+<dd>Unused</dd>
</dl></dd>
<dt><b>nstcalclr: (-1) [steps]</b></dt>
{
warning_error(wi, "rlist should be >= 0");
}
+ sprintf(err_buf, "nstlist can not be smaller than 0. (If you were trying to use the heuristic neighbour-list update scheme for efficient buffering for improved energy conservation, please use the Verlet cut-off scheme instead.)");
+ CHECK(ir->nstlist < 0);
process_interaction_modifier(ir, &ir->coulomb_modifier);
process_interaction_modifier(ir, &ir->vdw_modifier);
{
warning_error(wi, "rlistlong can not be shorter than rlist");
}
- if (IR_TWINRANGE(*ir) && ir->nstlist <= 0)
+ if (IR_TWINRANGE(*ir) && ir->nstlist == 0)
{
- warning_error(wi, "Can not have nstlist<=0 with twin-range interactions");
+ warning_error(wi, "Can not have nstlist == 0 with twin-range interactions");
}
}
(ir->ePBC != epbcNONE) ||
(ir->rcoulomb != 0.0) || (ir->rvdw != 0.0));
- if (ir->nstlist < 0)
- {
- warning_error(wi, "Can not have heuristic neighborlist updates without cut-off");
- }
if (ir->nstlist > 0)
{
warning_note(wi, "Simulating without cut-offs can be (slightly) faster with nstlist=0, nstype=simple and only one MPI rank");
if (ir->cutoff_scheme == ecutsGROUP)
{
if (((ir->coulomb_modifier != eintmodNONE && ir->rcoulomb == ir->rlist) ||
- (ir->vdw_modifier != eintmodNONE && ir->rvdw == ir->rlist)) &&
- ir->nstlist != 1)
+ (ir->vdw_modifier != eintmodNONE && ir->rvdw == ir->rlist)))
{
warning_note(wi, "With exact cut-offs, rlist should be "
"larger than rcoulomb and rvdw, so that there "
warning_note(wi, "You have selected user tables with dispersion correction, the dispersion will be corrected to -C6/r^6 beyond rvdw_switch (the tabulated interaction between rvdw_switch and rvdw will not be double counted). Make sure that you really want dispersion correction to -C6/r^6.");
}
- if (ir->nstlist == -1)
- {
- sprintf(err_buf, "With nstlist=-1 rvdw and rcoulomb should be smaller than rlist to account for diffusion and possibly charge-group radii");
- CHECK(ir->rvdw >= ir->rlist || ir->rcoulomb >= ir->rlist);
- }
- sprintf(err_buf, "nstlist can not be smaller than -1");
- CHECK(ir->nstlist < -1);
-
if (ir->eI == eiLBFGS && (ir->coulombtype == eelCUT || ir->vdwtype == evdwCUT)
&& ir->rvdw != 0)
{
int read_nblist(FILE *in, FILE *out, int **mat, int natoms, gmx_bool bSymm);
/* Returns total number of neighbors. If bSymm the matrix is symmetrized. */
-int natoms_beyond_ns_buffer(t_inputrec *ir, t_forcerec *fr, t_block *cgs,
- matrix scale_tot, rvec *x);
-/* Returns the number of atoms that moved beyond the ns buffer */
-
void reallocate_nblist(t_nblist *nl);
/* List reallocation, only exported for Verlet scheme use with FEP */
*
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
- * Copyright (c) 2011, by the GROMACS development team, led by
+ * Copyright (c) 2011,2014, 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.
data). This means that the only meaningful values are positive,
negative or zero. */
enum {
- eglsNABNSB, eglsCHKPT, eglsSTOPCOND, eglsRESETCOUNTERS, eglsNR
+ eglsCHKPT, eglsSTOPCOND, eglsRESETCOUNTERS, eglsNR
};
typedef struct {
+++ /dev/null
-/*
- * This file is part of the GROMACS molecular simulation package.
- *
- * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
- * Copyright (c) 2001-2004, The GROMACS development team.
- * Copyright (c) 2011,2013, 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.
- *
- * GROMACS is free software; you can redistribute it and/or
- * modify it under the terms of the GNU Lesser General Public License
- * as published by the Free Software Foundation; either version 2.1
- * of the License, or (at your option) any later version.
- *
- * GROMACS is distributed in the hope that it will be useful,
- * but WITHOUT ANY WARRANTY; without even the implied warranty of
- * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
- * Lesser General Public License for more details.
- *
- * You should have received a copy of the GNU Lesser General Public
- * License along with GROMACS; if not, see
- * http://www.gnu.org/licenses, or write to the Free Software Foundation,
- * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
- *
- * If you want to redistribute modifications to GROMACS, please
- * consider that scientific software is very special. Version
- * control is crucial - bugs must be traceable. We will be happy to
- * consider code for inclusion in the official distribution, but
- * derived work must not be called official GROMACS. Details are found
- * in the README & COPYING files - if they are missing, get the
- * official version at http://www.gromacs.org.
- *
- * To help us fund GROMACS development, we humbly ask that you cite
- * the research papers on the package. Check out http://www.gromacs.org.
- */
-#ifndef _nlistheuristics_h
-#define _nlistheuristics_h
-
-#ifdef __cplusplus
-extern "C" {
-#endif
-
-#if 0
-}
-/* Hack to make automatic indenting work */
-#endif
-
-typedef struct {
- gmx_bool bGStatEveryStep;
- gmx_int64_t step_ns;
- gmx_int64_t step_nscheck;
- gmx_int64_t nns;
- matrix scale_tot;
- int nabnsb;
- double s1;
- double s2;
- double ab;
- double lt_runav;
- double lt_runav2;
-} gmx_nlheur_t;
-
-void reset_nlistheuristics(gmx_nlheur_t *nlh, gmx_int64_t step);
-
-void init_nlistheuristics(gmx_nlheur_t *nlh,
- gmx_bool bGStatEveryStep, gmx_int64_t step);
-
-void update_nliststatistics(gmx_nlheur_t *nlh, gmx_int64_t step);
-
-void set_nlistheuristics(gmx_nlheur_t *nlh, gmx_bool bReset, gmx_int64_t step);
-
-#ifdef __cplusplus
-}
-#endif
-
-#endif
t_mdatoms *md,
t_state *state,
rvec force[], /* forces on home particles */
- matrix *scale_tot,
matrix pcoupl_mu,
t_nrnb *nrnb,
gmx_update_t upd);
#include "gromacs/utility/smalloc.h"
/* Is the signal in one simulation independent of other simulations? */
-gmx_bool gs_simlocal[eglsNR] = { TRUE, FALSE, FALSE, TRUE };
+gmx_bool gs_simlocal[eglsNR] = { FALSE, FALSE, TRUE };
/* check which of the multisim simulations has the shortest number of
steps and return that number of nsteps */
{
real t;
gmx_bool bNS;
- int nabnsb;
tensor force_vir, shake_vir, ekin;
real dvdl_constr, prescorr, enercorr, dvdlcorr;
real terminate = 0;
{
bNS = TRUE;
}
- else if (inputrec->nstlist == -1)
- {
- nabnsb = natoms_beyond_ns_buffer(inputrec, fr, &top->cgs, NULL, ems->s.x);
- if (PAR(cr))
- {
- gmx_sumi(1, &nabnsb, cr);
- }
- bNS = (nabnsb > 0);
- }
}
if (vsite)
+++ /dev/null
-/*
- * This file is part of the GROMACS molecular simulation package.
- *
- * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
- * Copyright (c) 2001-2004, The GROMACS development team.
- * Copyright (c) 2012,2013,2014, 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.
- *
- * GROMACS is free software; you can redistribute it and/or
- * modify it under the terms of the GNU Lesser General Public License
- * as published by the Free Software Foundation; either version 2.1
- * of the License, or (at your option) any later version.
- *
- * GROMACS is distributed in the hope that it will be useful,
- * but WITHOUT ANY WARRANTY; without even the implied warranty of
- * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
- * Lesser General Public License for more details.
- *
- * You should have received a copy of the GNU Lesser General Public
- * License along with GROMACS; if not, see
- * http://www.gnu.org/licenses, or write to the Free Software Foundation,
- * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
- *
- * If you want to redistribute modifications to GROMACS, please
- * consider that scientific software is very special. Version
- * control is crucial - bugs must be traceable. We will be happy to
- * consider code for inclusion in the official distribution, but
- * derived work must not be called official GROMACS. Details are found
- * in the README & COPYING files - if they are missing, get the
- * official version at http://www.gromacs.org.
- *
- * To help us fund GROMACS development, we humbly ask that you cite
- * the research papers on the package. Check out http://www.gromacs.org.
- */
-#include "gmxpre.h"
-
-#include "gromacs/legacyheaders/typedefs.h"
-#include "gromacs/legacyheaders/types/nlistheuristics.h"
-#include "gromacs/math/vec.h"
-#include "gromacs/utility/cstringutil.h"
-#include "gromacs/utility/fatalerror.h"
-
-void reset_nlistheuristics(gmx_nlheur_t *nlh, gmx_int64_t step)
-{
- nlh->lt_runav = 0;
- nlh->lt_runav2 = 0;
- nlh->step_nscheck = step;
-}
-
-void init_nlistheuristics(gmx_nlheur_t *nlh,
- gmx_bool bGStatEveryStep, gmx_int64_t step)
-{
- nlh->bGStatEveryStep = bGStatEveryStep;
- nlh->nns = 0;
- nlh->nabnsb = 0;
- nlh->s1 = 0;
- nlh->s2 = 0;
- nlh->ab = 0;
-
- reset_nlistheuristics(nlh, step);
-}
-
-void update_nliststatistics(gmx_nlheur_t *nlh, gmx_int64_t step)
-{
- gmx_int64_t nl_lt;
- char sbuf[STEPSTRSIZE], sbuf2[STEPSTRSIZE];
-
- /* Determine the neighbor list life time */
- nl_lt = step - nlh->step_ns;
- if (debug)
- {
- fprintf(debug, "%d atoms beyond ns buffer, updating neighbor list after %s steps\n", nlh->nabnsb, gmx_step_str(nl_lt, sbuf));
- }
- nlh->nns++;
- nlh->s1 += nl_lt;
- nlh->s2 += nl_lt*nl_lt;
- nlh->ab += nlh->nabnsb;
- if (nlh->lt_runav == 0)
- {
- nlh->lt_runav = nl_lt;
- /* Initialize the fluctuation average
- * such that at startup we check after 0 steps.
- */
- nlh->lt_runav2 = sqr(nl_lt/2.0);
- }
- /* Running average with 0.9 gives an exp. history of 9.5 */
- nlh->lt_runav2 = 0.9*nlh->lt_runav2 + 0.1*sqr(nlh->lt_runav - nl_lt);
- nlh->lt_runav = 0.9*nlh->lt_runav + 0.1*nl_lt;
- if (nlh->bGStatEveryStep)
- {
- /* Always check the nlist validity */
- nlh->step_nscheck = step;
- }
- else
- {
- /* We check after: <life time> - 2*sigma
- * The factor 2 is quite conservative,
- * but we assume that with nstlist=-1 the user
- * prefers exact integration over performance.
- */
- nlh->step_nscheck = step
- + (int)(nlh->lt_runav - 2.0*sqrt(nlh->lt_runav2)) - 1;
- }
- if (debug)
- {
- fprintf(debug, "nlist life time %s run av. %4.1f sig %3.1f check %s check with -gcom %d\n",
- gmx_step_str(nl_lt, sbuf), nlh->lt_runav, sqrt(nlh->lt_runav2),
- gmx_step_str(nlh->step_nscheck-step+1, sbuf2),
- (int)(nlh->lt_runav - 2.0*sqrt(nlh->lt_runav2)));
- }
-}
-
-void set_nlistheuristics(gmx_nlheur_t *nlh, gmx_bool bReset, gmx_int64_t step)
-{
- int d;
-
- if (bReset)
- {
- reset_nlistheuristics(nlh, step);
- }
- else
- {
- update_nliststatistics(nlh, step);
- }
-
- nlh->step_ns = step;
- /* Initialize the cumulative coordinate scaling matrix */
- clear_mat(nlh->scale_tot);
- for (d = 0; d < DIM; d++)
- {
- nlh->scale_tot[d][d] = 1.0;
- }
-}
return nsearch;
}
-
-int natoms_beyond_ns_buffer(t_inputrec *ir, t_forcerec *fr, t_block *cgs,
- matrix scale_tot, rvec *x)
-{
- int cg0, cg1, cg, a0, a1, a, i, j;
- real rint, hbuf2, scale;
- rvec *cg_cm, cgsc;
- gmx_bool bIsotropic;
- int nBeyond;
-
- nBeyond = 0;
-
- rint = max(ir->rcoulomb, ir->rvdw);
- if (ir->rlist < rint)
- {
- gmx_fatal(FARGS, "The neighbor search buffer has negative size: %f nm",
- ir->rlist - rint);
- }
- cg_cm = fr->cg_cm;
-
- cg0 = fr->cg0;
- cg1 = fr->hcg;
-
- if (!EI_DYNAMICS(ir->eI) || !DYNAMIC_BOX(*ir))
- {
- hbuf2 = sqr(0.5*(ir->rlist - rint));
- for (cg = cg0; cg < cg1; cg++)
- {
- a0 = cgs->index[cg];
- a1 = cgs->index[cg+1];
- for (a = a0; a < a1; a++)
- {
- if (distance2(cg_cm[cg], x[a]) > hbuf2)
- {
- nBeyond++;
- }
- }
- }
- }
- else
- {
- bIsotropic = TRUE;
- scale = scale_tot[0][0];
- for (i = 1; i < DIM; i++)
- {
- /* With anisotropic scaling, the original spherical ns volumes become
- * ellipsoids. To avoid costly transformations we use the minimum
- * eigenvalue of the scaling matrix for determining the buffer size.
- * Since the lower half is 0, the eigenvalues are the diagonal elements.
- */
- scale = min(scale, scale_tot[i][i]);
- if (scale_tot[i][i] != scale_tot[i-1][i-1])
- {
- bIsotropic = FALSE;
- }
- for (j = 0; j < i; j++)
- {
- if (scale_tot[i][j] != 0)
- {
- bIsotropic = FALSE;
- }
- }
- }
- hbuf2 = sqr(0.5*(scale*ir->rlist - rint));
- if (bIsotropic)
- {
- for (cg = cg0; cg < cg1; cg++)
- {
- svmul(scale, cg_cm[cg], cgsc);
- a0 = cgs->index[cg];
- a1 = cgs->index[cg+1];
- for (a = a0; a < a1; a++)
- {
- if (distance2(cgsc, x[a]) > hbuf2)
- {
- nBeyond++;
- }
- }
- }
- }
- else
- {
- /* Anistropic scaling */
- for (cg = cg0; cg < cg1; cg++)
- {
- /* Since scale_tot contains the transpose of the scaling matrix,
- * we need to multiply with the transpose.
- */
- tmvmul_ur0(scale_tot, cg_cm[cg], cgsc);
- a0 = cgs->index[cg];
- a1 = cgs->index[cg+1];
- for (a = a0; a < a1; a++)
- {
- if (distance2(cgsc, x[a]) > hbuf2)
- {
- nBeyond++;
- }
- }
- }
- }
- }
-
- return nBeyond;
-}
}
static void deform(gmx_update_t upd,
- int start, int homenr, rvec x[], matrix box, matrix *scale_tot,
+ int start, int homenr, rvec x[], matrix box,
const t_inputrec *ir, gmx_int64_t step)
{
matrix bnew, invbox, mu;
x[i][YY] = mu[YY][YY]*x[i][YY]+mu[ZZ][YY]*x[i][ZZ];
x[i][ZZ] = mu[ZZ][ZZ]*x[i][ZZ];
}
- if (scale_tot != NULL)
- {
- /* The transposes of the scaling matrices are stored,
- * so we need to do matrix multiplication in the inverse order.
- */
- mmul_ur0(*scale_tot, mu, *scale_tot);
- }
}
void update_tcouple(gmx_int64_t step,
t_mdatoms *md,
t_state *state,
rvec force[], /* forces on home particles */
- matrix *scale_tot,
matrix pcoupl_mu,
t_nrnb *nrnb,
gmx_update_t upd)
break;
}
- if ((!(IR_NPT_TROTTER(inputrec) || IR_NPH_TROTTER(inputrec))) && scale_tot)
- {
- /* The transposes of the scaling matrices are stored,
- * therefore we need to reverse the order in the multiplication.
- */
- mmul_ur0(*scale_tot, pcoupl_mu, *scale_tot);
- }
-
if (DEFORM(*inputrec))
{
- deform(upd, start, homenr, state->x, state->box, scale_tot, inputrec, step);
+ deform(upd, start, homenr, state->x, state->box, inputrec, step);
}
where();
dump_it_all(fplog, "After update",
#include "gromacs/legacyheaders/types/iteratedconstraints.h"
#include "gromacs/legacyheaders/types/mdatom.h"
#include "gromacs/legacyheaders/types/membedt.h"
-#include "gromacs/legacyheaders/types/nlistheuristics.h"
#include "gromacs/legacyheaders/types/nrnb.h"
#include "gromacs/legacyheaders/types/oenv.h"
#include "gromacs/legacyheaders/types/shellfc.h"
t_vcm *vcm;
t_state *bufstate = NULL;
matrix pcoupl_mu, M;
- gmx_nlheur_t nlh;
t_trxframe rerun_fr;
gmx_repl_ex_t repl_ex = NULL;
int nchkpt = 1;
step = ir->init_step;
step_rel = 0;
- init_nlistheuristics(&nlh, bGStatEveryStep, step);
-
if (MULTISIM(cr) && (repl_ex_nst <= 0 ))
{
/* check how many steps are left in other sims */
/* Determine whether or not to do Neighbour Searching and LR */
bNStList = (ir->nstlist > 0 && step % ir->nstlist == 0);
- bNS = (bFirstStep || bExchanged || bNeedRepartition || bNStList || bDoFEP ||
- (ir->nstlist == -1 && nlh.nabnsb > 0));
-
- if (bNS && ir->nstlist == -1)
- {
- set_nlistheuristics(&nlh, bFirstStep || bExchanged || bNeedRepartition || bDoFEP, step);
- }
+ bNS = (bFirstStep || bExchanged || bNeedRepartition || bNStList || bDoFEP);
}
/* check whether we should stop because another simulation has
/* Do we need global communication ? */
bGStat = (bCalcVir || bCalcEner || bStopCM ||
- do_per_step(step, nstglobalcomm) || (bVV && IR_NVT_TROTTER(ir) && do_per_step(step-1, nstglobalcomm)) ||
- (ir->nstlist == -1 && !bRerunMD && step >= nlh.step_nscheck));
+ do_per_step(step, nstglobalcomm) || (bVV && IR_NVT_TROTTER(ir) && do_per_step(step-1, nstglobalcomm)));
do_ene = (do_per_step(step, ir->nstenergy) || bLastStep);
gs.sig[eglsRESETCOUNTERS] = 1;
}
- if (ir->nstlist == -1 && !bRerunMD)
- {
- /* When bGStatEveryStep=FALSE, global_stat is only called
- * when we check the atom displacements, not at NS steps.
- * This means that also the bonded interaction count check is not
- * performed immediately after NS. Therefore a few MD steps could
- * be performed with missing interactions.
- * But wrong energies are never written to file,
- * since energies are only written after global_stat
- * has been called.
- */
- if (step >= nlh.step_nscheck)
- {
- nlh.nabnsb = natoms_beyond_ns_buffer(ir, fr, &top->cgs,
- nlh.scale_tot, state->x);
- }
- else
- {
- /* This is not necessarily true,
- * but step_nscheck is determined quite conservatively.
- */
- nlh.nabnsb = 0;
- }
- }
-
/* In parallel we only have to check for checkpointing in steps
* where we do global communication,
* otherwise the other nodes don't know.
*/
if (bGStat || (!EI_VV(ir->eI) && do_per_step(step+1, nstglobalcomm)))
{
- if (ir->nstlist == -1 && bFirstIterate)
- {
- gs.sig[eglsNABNSB] = nlh.nabnsb;
- }
compute_globals(fplog, gstat, cr, ir, fr, ekind, state, state_global, mdatoms, nrnb, vcm,
wcycle, enerd, force_vir, shake_vir, total_vir, pres, mu_tot,
constr,
| (bFirstIterate ? CGLO_FIRSTITERATE : 0)
| CGLO_CONSTRAINT
);
- if (ir->nstlist == -1 && bFirstIterate)
- {
- nlh.nabnsb = gs.set[eglsNABNSB];
- gs.set[eglsNABNSB] = 0;
- }
}
/* bIterate is set to keep it from eliminating the old ekin kinetic energy terms */
/* ############# END CALC EKIN AND PRESSURE ################# */
sum_dhdl(enerd, state->lambda, ir->fepvals);
}
update_box(fplog, step, ir, mdatoms, state, f,
- ir->nstlist == -1 ? &nlh.scale_tot : NULL, pcoupl_mu, nrnb, upd);
+ pcoupl_mu, nrnb, upd);
/* ################# END UPDATE STEP 2 ################# */
/* #### We now have r(t+dt) and v(t+dt/2) ############# */
done_mdoutf(outf);
debug_gmx();
- if (ir->nstlist == -1 && nlh.nns > 0 && fplog)
- {
- fprintf(fplog, "Average neighborlist lifetime: %.1f steps, std.dev.: %.1f steps\n", nlh.s1/nlh.nns, sqrt(nlh.s2/nlh.nns - sqr(nlh.s1/nlh.nns)));
- fprintf(fplog, "Average number of atoms that crossed the half buffer length: %.1f\n\n", nlh.ab/nlh.nns);
- }
-
if (pme_loadbal != NULL)
{
pme_loadbal_done(pme_loadbal, cr, fplog,