# How does verlet neighbour searching list set rlist for dynamic simulations **URL:** \ **Category:** User discussions\ **Created:** [June 21, 2022, 9:55am UTC](https://gromacs.bioexcel.eu/t/how-does-verlet-neighbour-searching-list-set-rlist-for-dynamic-simulations/4328 "2022-06-21T09:55:28Z")\ **Posts on this page:** 2\ **Page:** 1
**Author:** ![Audreyl](https://avatars.discourse-cdn.com/v4/letter/a/e19b73/32.png) [@Audreyl](https://gromacs.bioexcel.eu/u/Audreyl)\ **Post date:** [June 21, 2022, 9:55am UTC](https://gromacs.bioexcel.eu/t/how-does-verlet-neighbour-searching-list-set-rlist-for-dynamic-simulations/4328/1 "2022-06-21T09:55:28Z")
I am running dynamics simulations using Gromacs 2021.5. In my MDP file, I use Verlet for Neighbour searching. I understand from the manual that the buffer is set automatically and rlist is ignored. However, when I run my simulations the rlist is always set to 1nm as per my log files unless I assign the rlist a different value. **This is what my log file looks like without setting rlist** Input Parameters: integrator = md tinit = 0 dt = 0.001 nsteps = 3000000 init-step = 0 simulation-part = 1 mts = false comm-mode = Linear nstcomm = 100 bd-fric = 0 ld-seed = -8388813 emtol = 10 emstep = 0.01 niter = 20 fcstep = 0 nstcgsteep = 1000 nbfgscorr = 10 rtpi = 0.05 nstxout = 2000 nstvout = 2000 nstfout = 0 nstlog = 2000 nstcalcenergy = 100 nstenergy = 2000 nstxout-compressed = 0 compressed-x-precision = 1000 cutoff-scheme = Verlet nstlist = 10 pbc = xyz periodic-molecules = false verlet-buffer-tolerance = 0.005 rlist = 1 coulombtype = PME coulomb-modifier = Potential-shift rcoulomb-switch = 0 rcoulomb = 1 epsilon-r = 1 epsilon-rf = inf vdw-type = PME vdw-modifier = Potential-shift rvdw-switch = 0 rvdw = 1 DispCorr = No table-extension = 1 fourierspacing = 0.12 fourier-nx = 96 fourier-ny = 96 fourier-nz = 96 pme-order = 4 ewald-rtol = 1e-05 ewald-rtol-lj = 0.001 lj-pme-comb-rule = Geometric ewald-geometry = 0 epsilon-surface = 0 tcoupl = V-rescale nsttcouple = 10 nh-chain-length = 0 print-nose-hoover-chain-variables = false pcoupl = Berendsen pcoupltype = Isotropic nstpcouple = 10 tau-p = 1 compressibility (3x3): compressibility[0]={ 4.50000e-05, 0.00000e+00, 0.00000e+00} compressibility[1]={ 0.00000e+00, 4.50000e-05, 0.00000e+00} compressibility[2]={ 0.00000e+00, 0.00000e+00, 4.50000e-05} ref-p (3x3): ref-p[0]={ 1.00000e+02, 0.00000e+00, 0.00000e+00} ref-p[1]={ 0.00000e+00, 1.00000e+02, 0.00000e+00} ref-p[2]={ 0.00000e+00, 0.00000e+00, 1.00000e+02} refcoord-scaling = No posres-com (3): posres-com[0]= 0.00000e+00 posres-com[1]= 0.00000e+00 posres-com[2]= 0.00000e+00 posres-comB (3): posres-comB[0]= 0.00000e+00 posres-comB[1]= 0.00000e+00 posres-comB[2]= 0.00000e+00 QMMM = false qm-opts: ngQM = 0 constraint-algorithm = Lincs continuation = false Shake-SOR = false shake-tol = 0.0001 lincs-order = 4 lincs-iter = 1 lincs-warnangle = 30 nwall = 0 wall-type = 9-3 wall-r-linpot = -1 wall-atomtype[0] = -1 wall-atomtype[1] = -1 wall-density[0] = 0 wall-density[1] = 0 wall-ewald-zfac = 3 pull = false awh = false rotation = false interactiveMD = false disre = No disre-weighting = Conservative disre-mixed = false dr-fc = 1000 dr-tau = 0 nstdisreout = 100 orire-fc = 0 orire-tau = 0 nstorireout = 100 free-energy = no cos-acceleration = 0 deform (3x3): deform[0]={ 0.00000e+00, 0.00000e+00, 0.00000e+00} deform[1]={ 0.00000e+00, 0.00000e+00, 0.00000e+00} deform[2]={ 0.00000e+00, 0.00000e+00, 0.00000e+00} simulated-tempering = false swapcoords = no userint1 = 0 userint2 = 0 userint3 = 0 userint4 = 0 userreal1 = 0 userreal2 = 0 userreal3 = 0 userreal4 = 0 applied-forces: electric-field: x: E0 = 0 omega = 0 t0 = 0 sigma = 0 y: E0 = 0 omega = 0 t0 = 0 sigma = 0 z: E0 = 0 omega = 0 t0 = 0 sigma = 0 density-guided-simulation: active = false group = protein similarity-measure = inner-product atom-spreading-weight = unity force-constant = 1e+09 gaussian-transform-spreading-width = 0.2 gaussian-transform-spreading-range-in-multiples-of-width = 4 reference-density-filename = reference.mrc nst = 1 normalize-densities = true adaptive-force-scaling = false adaptive-force-scaling-time-constant = 4 shift-vector = transformation-matrix = grpopts: nrdf: 26797 ref-t: 300 tau-t: 0.1 annealing: No annealing-npoints: 0 acc: 0 0 0 nfreeze: N N N energygrp-flags[0]: 0 Changing nstlist from 10 to 100, rlist from 1 to 1 1 GPU selected for this run. Mapping of GPU IDs to the 1 GPU task in the 1 rank on this node: PP:0 PP tasks will do (non-perturbed) short-ranged and most bonded interactions on the GPU PP task will update and constrain coordinates on the CPU Using 1 MPI thread Non-default thread affinity set, disabling internal thread affinity Using 20 OpenMP threads System total charge: -0.000 Will do PME sum in reciprocal space for LJ dispersion interactions. ++++ PLEASE READ AND CITE THE FOLLOWING REFERENCE ++++ U. Essmann, L. Perera, M. L. Berkowitz, T. Darden, H. Lee and L. G. Pedersen A smooth particle mesh Ewald method J. Chem. Phys. 103 (1995) pp. 8577-8592 -------- -------- — Thank You — -------- -------- Using a Gaussian width (1/beta) of 0.298423 nm for LJ Ewald Will do PME sum in reciprocal space for electrostatic interactions. ++++ PLEASE READ AND CITE THE FOLLOWING REFERENCE ++++ U. Essmann, L. Perera, M. L. Berkowitz, T. Darden, H. Lee and L. G. Pedersen A smooth particle mesh Ewald method J. Chem. Phys. 103 (1995) pp. 8577-8592 -------- -------- — Thank You — -------- -------- Using a Gaussian width (1/beta) of 0.320163 nm for Ewald Potential shift: LJ r^-12: -1.000e+00 r^-6: -1.000e-03, Ewald -1.000e-05 Initialized non-bonded Ewald tables, spacing: 9.33e-04 size: 1073 Using shifted Lennard-Jones, switch between 0 and 1 nm Generated table with 1000 data points for 1-4 COUL. Tabscale = 500 points/nm Generated table with 1000 data points for 1-4 LJ6. Tabscale = 500 points/nm Generated table with 1000 data points for 1-4 LJ12. Tabscale = 500 points/nm Using GPU 8x8 nonbonded short-range kernels Using a 8x8 pair-list setup: updated every 100 steps, buffer 0.000 nm, rlist 1.000 nm At tolerance 0.005 kJ/mol/ps per atom, equivalent classical 1x1 list would be: updated every 100 steps, buffer 0.000 nm, rlist 1.000 nm Removing pbc first time \*_This is what my log file with setting rlist=1.5 look like_ Input Parameters: integrator = md tinit = 0 dt = 0.001 nsteps = 3500000 init-step = 0 simulation-part = 1 mts = false comm-mode = Linear nstcomm = 100 bd-fric = 0 ld-seed = 1945558523 emtol = 10 emstep = 0.01 niter = 20 fcstep = 0 nstcgsteep = 1000 nbfgscorr = 10 rtpi = 0.05 nstxout = 5000 nstvout = 5000 nstfout = 0 nstlog = 5000 nstcalcenergy = 100 nstenergy = 5000 nstxout-compressed = 0 compressed-x-precision = 1000 cutoff-scheme = Verlet nstlist = 10 pbc = xyz periodic-molecules = false verlet-buffer-tolerance = 0.005 rlist = 1.5 coulombtype = PME coulomb-modifier = Potential-shift rcoulomb-switch = 0 rcoulomb = 1.5 epsilon-r = 1 epsilon-rf = inf vdw-type = PME vdw-modifier = Potential-shift rvdw-switch = 0 rvdw = 1.5 DispCorr = EnerPres table-extension = 1 fourierspacing = 0.12 fourier-nx = 84 fourier-ny = 84 fourier-nz = 84 pme-order = 4 ewald-rtol = 1e-05 ewald-rtol-lj = 0.001 lj-pme-comb-rule = Geometric ewald-geometry = 0 epsilon-surface = 0 tcoupl = V-rescale nsttcouple = 10 nh-chain-length = 0 print-nose-hoover-chain-variables = false pcoupl = Berendsen pcoupltype = Isotropic nstpcouple = 10 tau-p = 1 compressibility (3x3): compressibility[0]={ 4.50000e-05, 0.00000e+00, 0.00000e+00} compressibility[1]={ 0.00000e+00, 4.50000e-05, 0.00000e+00} compressibility[2]={ 0.00000e+00, 0.00000e+00, 4.50000e-05} ref-p (3x3): ref-p[0]={ 1.00000e+02, 0.00000e+00, 0.00000e+00} ref-p[1]={ 0.00000e+00, 1.00000e+02, 0.00000e+00} ref-p[2]={ 0.00000e+00, 0.00000e+00, 1.00000e+02} refcoord-scaling = No posres-com (3): posres-com[0]= 0.00000e+00 posres-com[1]= 0.00000e+00 posres-com[2]= 0.00000e+00 posres-comB (3): posres-comB[0]= 0.00000e+00 posres-comB[1]= 0.00000e+00 posres-comB[2]= 0.00000e+00 QMMM = false qm-opts: ngQM = 0 constraint-algorithm = Lincs continuation = false Shake-SOR = false shake-tol = 0.0001 lincs-order = 4 lincs-iter = 1 lincs-warnangle = 30 nwall = 0 wall-type = 9-3 wall-r-linpot = -1 wall-atomtype[0] = -1 wall-atomtype[1] = -1 wall-density[0] = 0 wall-density[1] = 0 wall-ewald-zfac = 3 pull = false awh = false rotation = false interactiveMD = false disre = No disre-weighting = Conservative disre-mixed = false dr-fc = 1000 dr-tau = 0 nstdisreout = 100 orire-fc = 0 orire-tau = 0 nstorireout = 100 free-energy = no cos-acceleration = 0 deform (3x3): deform[0]={ 0.00000e+00, 0.00000e+00, 0.00000e+00} deform[1]={ 0.00000e+00, 0.00000e+00, 0.00000e+00} deform[2]={ 0.00000e+00, 0.00000e+00, 0.00000e+00} simulated-tempering = false swapcoords = no userint1 = 0 userint2 = 0 userint3 = 0 userint4 = 0 userreal1 = 0 userreal2 = 0 userreal3 = 0 userreal4 = 0 applied-forces: electric-field: x: E0 = 0 omega = 0 t0 = 0 sigma = 0 y: E0 = 0 omega = 0 t0 = 0 sigma = 0 z: E0 = 0 omega = 0 t0 = 0 sigma = 0 density-guided-simulation: active = false group = protein similarity-measure = inner-product atom-spreading-weight = unity force-constant = 1e+09 gaussian-transform-spreading-width = 0.2 gaussian-transform-spreading-range-in-multiples-of-width = 4 reference-density-filename = reference.mrc nst = 1 normalize-densities = true adaptive-force-scaling = false adaptive-force-scaling-time-constant = 4 shift-vector = transformation-matrix = grpopts: nrdf: 23112 ref-t: 300 tau-t: 0.1 annealing: No annealing-npoints: 0 acc: 0 0 0 nfreeze: N N N energygrp-flags[0]: 0 Changing nstlist from 10 to 100, rlist from 1.5 to 1.5 Initializing Domain Decomposition on 4 ranks Dynamic load balancing: auto Minimum cell size due to atom displacement: 0.482 nm Initial maximum distances in bonded interactions: two-body bonded interactions: 0.404 nm, LJ-14, atoms 5258 5261 multi-body bonded interactions: 0.404 nm, Ryckaert-Bell., atoms 5258 5261 Minimum cell size due to bonded interactions: 0.444 nm Maximum distance for 5 constraints, at 120 deg. angles, all-trans: 0.218 nm Estimated maximum distance required for P-LINCS: 0.218 nm Scaling the initial minimum size with 1/0.8 (option -dds) = 1.25 Using 0 separate PME ranks Optimizing the DD grid for 4 cells with a minimum initial size of 0.603 nm The maximum allowed number of cells is: X 16 Y 16 Z 16 Domain decomposition grid 4 x 1 x 1, separate PME ranks 0 PME domain decomposition: 4 x 1 x 1 Domain decomposition rank 0, coordinates 0 0 0 The initial number of communication pulses is: X 1 The initial domain decomposition cell size is: X 2.50 nm The maximum allowed distance for atoms involved in interactions is: non-bonded interactions 1.500 nm (the following are initial values, they could change due to box deformation) two-body bonded interactions (-rdd) 1.500 nm multi-body bonded interactions (-rdd) 1.500 nm atoms separated by up to 5 constraints (-rcon) 2.500 nm When dynamic load balancing gets turned on, these settings will change to: The maximum number of communication pulses is: X 1 The minimum size for domain decomposition cells is 1.500 nm The requested allowed shrink of DD cells (option -dds) is: 0.80 The allowed shrink of domain decomposition cells is: X 0.60 The maximum allowed distance for atoms involved in interactions is: non-bonded interactions 1.500 nm two-body bonded interactions (-rdd) 1.500 nm multi-body bonded interactions (-rdd) 1.500 nm atoms separated by up to 5 constraints (-rcon) 1.500 nm On host r2i7n8 1 GPU selected for this run. Mapping of GPU IDs to the 4 GPU tasks in the 4 ranks on this node: PP:0,PP:0,PP:0,PP:0 PP tasks will do (non-perturbed) short-ranged and most bonded interactions on the GPU PP task will update and constrain coordinates on the CPU Using 4 MPI processes Non-default thread affinity set, disabling internal thread affinity Using 5 OpenMP threads per MPI process System total charge: -0.000 Will do PME sum in reciprocal space for LJ dispersion interactions. ++++ PLEASE READ AND CITE THE FOLLOWING REFERENCE ++++ U. Essmann, L. Perera, M. L. Berkowitz, T. Darden, H. Lee and L. G. Pedersen A smooth particle mesh Ewald method J. Chem. Phys. 103 (1995) pp. 8577-8592 -------- -------- — Thank You — -------- -------- Using a Gaussian width (1/beta) of 0.447634 nm for LJ Ewald Will do PME sum in reciprocal space for electrostatic interactions. ++++ PLEASE READ AND CITE THE FOLLOWING REFERENCE ++++ U. Essmann, L. Perera, M. L. Berkowitz, T. Darden, H. Lee and L. G. Pedersen A smooth particle mesh Ewald method J. Chem. Phys. 103 (1995) pp. 8577-8592 -------- -------- — Thank You — -------- -------- Using a Gaussian width (1/beta) of 0.480244 nm for Ewald Potential shift: LJ r^-12: -7.707e-03 r^-6: -8.779e-05, Ewald -6.667e-06 Initialized non-bonded Ewald tables, spacing: 1.14e-03 size: 1314 Using shifted Lennard-Jones, switch between 0 and 1.5 nm Generated table with 1250 data points for 1-4 COUL. Tabscale = 500 points/nm Generated table with 1250 data points for 1-4 LJ6. Tabscale = 500 points/nm Generated table with 1250 data points for 1-4 LJ12. Tabscale = 500 points/nm Long Range LJ corr.: 3.0097e-15 Using GPU 8x8 nonbonded short-range kernels Using a 8x8 pair-list setup: updated every 100 steps, buffer 0.000 nm, rlist 1.500 nm At tolerance 0.005 kJ/mol/ps per atom, equivalent classical 1x1 list would be: updated every 100 steps, buffer 0.000 nm, rlist 1.500 nm Removing pbc first time Initializing Parallel LINear Constraint Solver ++++ PLEASE READ AND CITE THE FOLLOWING REFERENCE ++++ B. Hess P-LINCS: A Parallel Linear Constraint Solver for molecular simulation J. Chem. Theory Comput. 4 (2008) pp. 116-122 -------- -------- — Thank You — -------- -------- The number of constraints is 6900 There are constraints between atoms in different decomposition domains, will communicate selected coordinates each lincs iteration Linking all bonded interactions to atoms Intra-simulation communication will occur every 10 steps. **I would like to understand how this works**. How is the rlist determined ?? Thank you. ---
**Author:** ![ebriand](https://dub1.discourse-cdn.com/flex017/user_avatar/gromacs.bioexcel.eu/ebriand/32/2314_2.png) [@ebriand](https://gromacs.bioexcel.eu/u/ebriand)\ **Post date:** [June 21, 2022, 11:54am UTC](https://gromacs.bioexcel.eu/t/how-does-verlet-neighbour-searching-list-set-rlist-for-dynamic-simulations/4328/2 "2022-06-21T11:54:34Z")
Under the verlet buffer scheme, the important parameter is `verlet-buffer-tolerance`. This parameter is the acceptable level of error/drift due to particles moving such that they should be in another particle short-range interaction list, but are not (because neighbor search happen every `nstlist` step, not every step). These borderline particles are included in the neighbor list, by taking all pairs within r\_{list} = r\_{cutoff} + r\_{buffer}, with r\_{buffer} some positive value. Under NVT condition, one can estimate the diffusion of the borderline particles , so that r\_{buffer} can be set to include them and thus satisfy the error tolerance. Additionally, at a given allowed error level, you can increase `nstlist` if you increase the buffer (and conversely). That gives a free parameter that can be tuned for performance. This is what it done when the `rlist`/`nstlist` are set automatically. If you really want to set nstlist/rlist exactly (which is not recommended in the general case), you need to set `verlet-buffer-tolerance` to -1. See [this part of the reference manual](https://manual.gromacs.org/documentation/current/reference-manual/algorithms/molecular-dynamics.html#neighbor-searching), [this part of the user manual](https://manual.gromacs.org/documentation/current/user-guide/mdp-options.html#neighbor-searching), and [this publication](https://aip.scitation.org/doi/full/10.1063/5.0018516).