How to choose the neighbor list buffer distance¶
Set the neighbor list buffer distance (hoomd.md.nlist.NeighborList.buffer
) to maximize the
performance of your simulation. The neighbor list recomputes itself more often when buffer
is
small, and the pair force computation takes more time when buffer
is large. There
is an optimal value between the extremes that strongly depends on system size, density, the hardware
device, the pair potential, step size, and more. Test your specific model, change buffer
and
measure the performance to find the optimal. For example:
import hoomd
kT = 1.5
sample_steps = 1000
# Prepare a MD simulation.
device = hoomd.device.auto_select()
simulation = hoomd.Simulation(device=device)
simulation.create_state_from_gsd(filename='spheres.gsd')
simulation.state.thermalize_particle_momenta(filter=hoomd.filter.All(), kT=kT)
neighbor_list = hoomd.md.nlist.Cell(buffer=0.4)
lj = hoomd.md.pair.LJ(nlist=neighbor_list)
lj.params[('A', 'A')] = dict(epsilon=1.0, sigma=1.0)
lj.r_cut[('A', 'A')] = 2.5
bussi = hoomd.md.methods.thermostats.Bussi(kT=kT)
constant_volume = hoomd.md.methods.ConstantVolume(filter=hoomd.filter.All(),
thermostat=bussi)
simulation.operations.integrator = hoomd.md.Integrator(
dt=0.001, methods=[constant_volume], forces=[lj])
# Complete GPU kernel autotuning before making sensitive timing measurements.
if isinstance(device, hoomd.device.GPU):
simulation.run(0)
while not simulation.operations.is_tuning_complete:
simulation.run(1000)
for buffer in [0, 0.05, 0.1, 0.2, 0.3]:
neighbor_list.buffer = buffer
simulation.run(sample_steps)
device.notice(f'buffer={buffer}: TPS={simulation.tps:0.3g}, '
f'num_builds={neighbor_list.num_builds}')
Example output:
buffer=0: TPS=212, num_builds=1001
buffer=0.05: TPS=584, num_builds=197
buffer=0.1: TPS=839, num_builds=102
buffer=0.2: TPS=954, num_builds=53
buffer=0.3: TPS=880, num_builds=37
Important
Ensure that you run sufficient steps to sample many neighbor list builds to properly sample the amortized time spent per build.
Tip
Measure the optimal value of buffer
in testing, then apply that fixed value in production
runs to avoid wasting time at non-optimal values.
See also
hoomd.md.tune.NeighborListBuffer
can automate this process.