hpmc.update

Overview

hpmc.update.boxmc Apply box updates to sample isobaric and related ensembles.
hpmc.update.clusters Equilibrate the system according to the geometric cluster algorithm (GCA).
hpmc.update.muvt Insert and remove particles in the muVT ensemble.
hpmc.update.remove_drift Remove the center of mass drift from a system restrained on a lattice.
hpmc.update.wall Apply wall updates with a user-provided python callback.

Details

HPMC updaters.

class hoomd.hpmc.update.boxmc(mc, betaP, seed)

Apply box updates to sample isobaric and related ensembles.

Parameters:
  • mc (hoomd.hpmc.integrate) – HPMC integrator object for system on which to apply box updates
  • betaP (float or hoomd.variant) – \(\frac{p}{k_{\mathrm{B}}T}\). (units of inverse area in 2D or inverse volume in 3D) Apply your chosen reduced pressure convention externally.
  • seed (int) – random number seed for MC box changes

One or more Monte Carlo move types are applied to evolve the simulation box. By default, no moves are applied. Activate desired move types using the following methods with a non-zero weight:

  • aspect() - box aspect ratio moves
  • length() - change box lengths independently
  • shear() - shear the box
  • volume() - scale the box lengths uniformly
  • ln_volume() - scale the box lengths uniformly with logarithmic increments

Pressure inputs to update.boxmc are defined as \(\beta P\). Conversions from a specific definition of reduced pressure \(P^*\) are left for the user to perform.

Note

All delta and weight values for all move types default to 0.

Example:

mc = hpmc.integrate.sphere(seed=415236, d=0.3)
boxMC = hpmc.update.boxmc(mc, betaP=1.0, seed=9876)
boxMC.set_betap(2.0)
boxMC.ln_volume(delta=0.01, weight=2.0)
boxMC.length(delta=(0.1,0.1,0.1), weight=4.0)
run(30) # perform approximately 10 volume moves and 20 length moves
aspect(delta=None, weight=None)

Enable/disable aspect ratio move and set parameters.

Parameters:
  • delta (float) – maximum relative change of aspect ratio
  • weight (float) – relative weight of this box move type relative to other box move types. 0 disables this move type.

Rescale aspect ratio along a randomly chosen dimension.

Note

When an argument is None, the value is left unchanged from its current state.

Example:

box_update.aspect(delta=0.01)
box_update.aspect(delta=0.01, weight=2)
box_update.aspect(delta=0.01, weight=0.15)
Returns:A dict with the current values of delta, and weight.
disable()

Disables the updater.

Examples:

updater.disable()

Executing the disable command will remove the updater from the system. Any hoomd.run() command executed after disabling an updater will not use that updater during the simulation. A disabled updater can be re-enabled with enable()

enable()

Enables the updater.

Example:

box_updater.set_params(isotropic=True)
run(1e5)
box_updater.disable()
update.box_resize(dLy = 10)
box_updater.enable()
run(1e5)

See updater base class documentation for more information

get_aspect_acceptance()

Get the average acceptance ratio for aspect changing moves.

Returns:The average aspect change acceptance for the last run

Example:

mc = hpmc.integrate.shape(..);
mc_shape_param[name].set(....);
box_update = hpmc.update.boxmc(mc, betaP=10, seed=1)
run(100)
a_accept = box_update.get_aspect_acceptance()
get_ln_volume_acceptance()

Get the average acceptance ratio for log(V) changing moves.

Returns:The average volume change acceptance for the last run

Example:

mc = hpmc.integrate.shape(..);
mc.shape_param[name].set(....);
box_update = hpmc.update.boxmc(mc, betaP=10, seed=1)
run(100)
v_accept = box_update.get_ln_volume_acceptance()
get_shear_acceptance()

Get the average acceptance ratio for shear changing moves.

Returns:The average shear change acceptance for the last run

Example:

mc = hpmc.integrate.shape(..);
mc.shape_param[name].set(....);
box_update = hpmc.update.boxmc(mc, betaP=10, seed=1)
run(100)
s_accept = box_update.get_shear_acceptance()
get_volume_acceptance()

Get the average acceptance ratio for volume changing moves.

Returns:The average volume change acceptance for the last run

Example:

mc = hpmc.integrate.shape(..);
mc.shape_param[name].set(....);
box_update = hpmc.update.boxmc(mc, betaP=10, seed=1)
run(100)
v_accept = box_update.get_volume_acceptance()
length(delta=None, weight=None)

Enable/disable isobaric box dimension move and set parameters.

Parameters:
  • delta (float or tuple) – maximum change of the box thickness for each pair of parallel planes connected by the corresponding box edges. I.e. maximum change of HOOMD-blue box parameters Lx, Ly, Lz. A single float x is equivalent to (x, x, x).
  • weight (float) – relative weight of this box move type relative to other box move types. 0 disables this move type.

Sample the isobaric distribution of box dimensions by rescaling the plane-to-plane distance of box faces, Lx, Ly, Lz (see Periodic boundary conditions).

Note

When an argument is None, the value is left unchanged from its current state.

Example:

box_update.length(delta=(0.01, 0.01, 0.0)) # 2D box changes
box_update.length(delta=(0.01, 0.01, 0.01), weight=2)
box_update.length(delta=0.01, weight=2)
box_update.length(delta=(0.10, 0.01, 0.01), weight=0.15) # sample Lx more aggressively
Returns:A dict with the current values of delta and weight.
ln_volume(delta=None, weight=None)

Enable/disable isobaric volume move and set parameters.

Parameters:
  • delta (float) – maximum change of ln(V) (where V is box area (2D) or volume (3D)).
  • weight (float) – relative weight of this box move type relative to other box move types. 0 disables this move type.

Sample the isobaric distribution of box volumes by rescaling the box.

Note

When an argument is None, the value is left unchanged from its current state.

Example:

box_update.ln_volume(delta=0.001)
box_update.ln_volume(delta=0.001, weight=2)
box_update.ln_volume(delta=0.001, weight=0.15)
Returns:A dict with the current values of delta and weight.
restore_state()

Restore the state information from the file used to initialize the simulations

set_betap(betaP)

Update the pressure set point for Metropolis Monte Carlo volume updates.

Parameters:betaP (float) or (hoomd.variant) – \(\frac{p}{k_{\mathrm{B}}T}\). (units of inverse area in 2D or inverse volume in 3D) Apply your chosen reduced pressure convention externally.
set_period(period)

Changes the updater period.

Parameters:period (int) – New period to set.

Examples:

updater.set_period(100);
updater.set_period(1);

While the simulation is running, the action of each updater is executed every period time steps. Changing the period does not change the phase set when the analyzer was first created.

shear(delta=None, weight=None, reduce=None)

Enable/disable box shear moves and set parameters.

Parameters:
  • delta (tuple) – maximum change of the box tilt factor xy, xz, yz.
  • reduce (float) – Maximum number of lattice vectors of shear to allow before applying lattice reduction. Shear of +/- 0.5 cannot be lattice reduced, so set to a value < 0.5 to disable (default 0) Note that due to precision errors, lattice reduction may introduce small overlaps which can be resolved, but which temporarily break detailed balance.
  • weight (float) – relative weight of this box move type relative to other box move types. 0 disables this move type.

Sample the distribution of box shear by adjusting the HOOMD-blue tilt factor parameters xy, xz, and yz. (see Periodic boundary conditions).

Note

When an argument is None, the value is left unchanged from its current state.

Example:

box_update.shear(delta=(0.01, 0.00, 0.0)) # 2D box changes
box_update.shear(delta=(0.01, 0.01, 0.01), weight=2)
box_update.shear(delta=(0.10, 0.01, 0.01), weight=0.15) # sample xy more aggressively
Returns:A dict with the current values of delta, weight, and reduce.
volume(delta=None, weight=None)

Enable/disable isobaric volume move and set parameters.

Parameters:
  • delta (float) – maximum change of the box area (2D) or volume (3D).
  • weight (float) – relative weight of this box move type relative to other box move types. 0 disables this move type.

Sample the isobaric distribution of box volumes by rescaling the box.

Note

When an argument is None, the value is left unchanged from its current state.

Example:

box_update.volume(delta=0.01)
box_update.volume(delta=0.01, weight=2)
box_update.volume(delta=0.01, weight=0.15)
Returns:A dict with the current values of delta and weight.
class hoomd.hpmc.update.clusters(mc, seed, period=1)

Equilibrate the system according to the geometric cluster algorithm (GCA).

The GCA as described in Liu and Lujten (2004), http://doi.org/10.1103/PhysRevLett.92.035504 is used for hard shape, patch interactions and depletants.

With depletants, Clusters are defined by a simple distance cut-off criterion. Two particles belong to the same cluster if the circumspheres of the depletant-excluded volumes overlap.

Supported moves include pivot moves (point reflection), line reflections (pi rotation around an axis), and type swaps. Only the pivot move is rejection free. With anisotropic particles, the pivot move cannot be used because it would create a chiral mirror image of the particle, and only line reflections are employed. Line reflections are not rejection free because of periodic boundary conditions, as discussed in Sinkovits et al. (2012), http://doi.org/10.1063/1.3694271 .

The type swap move works between two types of spherical particles and exchanges their identities.

The clusters updater support TBB execution on multiple CPU cores. See Installing binaries for more information on how to compile HOOMD with TBB support.

Parameters:
  • mc (hoomd.hpmc.integrate) – MC integrator.
  • seed (int) – The seed of the pseudo-random number generator (Needs to be the same across partitions of the same Gibbs ensemble)
  • period (int) – Number of timesteps between histogram evaluations.

Example:

mc = hpmc.integrate.sphere(seed=415236)
hpmc.update.clusters(mc=mc, seed=123)
disable()

Disables the updater.

Examples:

updater.disable()

Executing the disable command will remove the updater from the system. Any hoomd.run() command executed after disabling an updater will not use that updater during the simulation. A disabled updater can be re-enabled with enable()

enable()

Enables the updater.

Examples:

updater.enable()

See also

disable()

get_pivot_acceptance()

Get the average acceptance ratio for pivot moves

Returns:The average acceptance rate for pivot moves during the last run
get_reflection_acceptance()

Get the average acceptance ratio for reflection moves

Returns:The average acceptance rate for reflection moves during the last run
get_swap_acceptance()

Get the average acceptance ratio for swap moves

Returns:The average acceptance rate for type swap moves during the last run
restore_state()

Restore the state information from the file used to initialize the simulations

set_params(move_ratio=None, flip_probability=None, swap_move_ratio=None, delta_mu=None, swap_types=None)

Set options for the clusters moves.

Parameters:
  • move_ratio (float) – Set the ratio between pivot and reflection moves (default 0.5)
  • flip_probability (float) – Set the probability for transforming an individual cluster (default 0.5)
  • swap_move_ratio (float) – Set the ratio between type swap moves and geometric moves (default 0.5)
  • delta_mu (float) – The chemical potential difference between types to be swapped
  • swap_types (list) – A pair of two types whose identities are swapped

Note

When an argument is None, the value is left unchanged from its current state.

Example:

clusters = hpmc.update.clusters(mc, seed=123)
clusters.set_params(move_ratio = 1.0)
clusters.set_params(swap_types=['A','B'], delta_mu = -0.001)
set_period(period)

Changes the updater period.

Parameters:period (int) – New period to set.

Examples:

updater.set_period(100);
updater.set_period(1);

While the simulation is running, the action of each updater is executed every period time steps. Changing the period does not change the phase set when the analyzer was first created.

class hoomd.hpmc.update.muvt(mc, seed, period=1, transfer_types=None, ngibbs=1)

Insert and remove particles in the muVT ensemble.

Parameters:
  • mc (hoomd.hpmc.integrate) – MC integrator.
  • seed (int) – The seed of the pseudo-random number generator (Needs to be the same across partitions of the same Gibbs ensemble)
  • period (int) – Number of timesteps between histogram evaluations.
  • transfer_types (list) – List of type names that are being transferred from/to the reservoir or between boxes (if None, all types)
  • ngibbs (int) – The number of partitions to use in Gibbs ensemble simulations (if == 1, perform grand canonical muVT)

The muVT (or grand-canonical) ensemble simulates a system at constant fugacity.

Gibbs ensemble simulations are also supported, where particles and volume are swapped between two or more boxes. Every box correspond to one MPI partition, and can therefore run on multiple ranks. See hoomd.comm and the –nrank command line option for how to split a MPI task into partitions.

Note

Multiple Gibbs ensembles are also supported in a single parallel job, with the ngibbs option to update.muvt(), where the number of partitions can be a multiple of ngibbs.

Example:

mc = hpmc.integrate.sphere(seed=415236)
update.muvt(mc=mc, period)
disable()

Disables the updater.

Examples:

updater.disable()

Executing the disable command will remove the updater from the system. Any hoomd.run() command executed after disabling an updater will not use that updater during the simulation. A disabled updater can be re-enabled with enable()

enable()

Enables the updater.

Examples:

updater.enable()

See also

disable()

restore_state()

Restore the state information from the file used to initialize the simulations

set_fugacity(type, fugacity)

Change muVT fugacities.

Parameters:
  • type (str) – Particle type to set parameters for
  • fugacity (float) – Fugacity of this particle type (dimension of volume^-1)

Example:

muvt = hpmc.update.muvt(mc, period=10)
muvt.set_fugacity(type='A', fugacity=1.23)
variant = hoomd.variant.linear_interp(points=[(0,1e1), (1e5, 4.56)])
muvt.set_fugacity(type='A', fugacity=variant)
set_params(dV=None, move_ratio=None, transfer_ratio=None)

Set muVT parameters.

Parameters:
  • dV (float) – (if set) Set volume rescaling factor (dimensionless)
  • move_ratio (float) – (if set) Set the ratio between volume and exchange/transfer moves (applies to Gibbs ensemble)
  • transfer_ratio (float) – (if set) Set the ratio between transfer and exchange moves

Example:

muvt = hpmc.update.muvt(mc, period = 10)
muvt.set_params(dV=0.1)
muvt.set_params(n_trial=2)
muvt.set_params(move_ratio=0.05)
set_period(period)

Changes the updater period.

Parameters:period (int) – New period to set.

Examples:

updater.set_period(100);
updater.set_period(1);

While the simulation is running, the action of each updater is executed every period time steps. Changing the period does not change the phase set when the analyzer was first created.

class hoomd.hpmc.update.remove_drift(mc, external_lattice, period=1)

Remove the center of mass drift from a system restrained on a lattice.

Parameters:

The command hpmc.update.remove_drift sets up an updater that removes the center of mass drift of a system every period timesteps,

Example:

mc = hpmc.integrate.convex_polyhedron(seed=seed);
mc.shape_param.set("A", vertices=verts)
mc.set_params(d=0.005, a=0.005)
lattice = hpmc.compute.lattice_field(mc=mc, position=fcc_lattice, k=1000.0);
remove_drift = update.remove_drift(mc=mc, external_lattice=lattice, period=1000);
disable()

Disables the updater.

Examples:

updater.disable()

Executing the disable command will remove the updater from the system. Any hoomd.run() command executed after disabling an updater will not use that updater during the simulation. A disabled updater can be re-enabled with enable()

enable()

Enables the updater.

Examples:

updater.enable()

See also

disable()

restore_state()

Restore the state information from the file used to initialize the simulations

set_period(period)

Changes the updater period.

Parameters:period (int) – New period to set.

Examples:

updater.set_period(100);
updater.set_period(1);

While the simulation is running, the action of each updater is executed every period time steps. Changing the period does not change the phase set when the analyzer was first created.

class hoomd.hpmc.update.wall(mc, walls, py_updater, move_ratio, seed, period=1)

Apply wall updates with a user-provided python callback.

Parameters:
  • mc (hoomd.hpmc.integrate) – MC integrator.
  • walls (hoomd.hpmc.field.wall) – the wall class instance to be updated
  • py_updater (callable) – the python callback that performs the update moves. This must be a python method that is a function of the timestep of the simulation. It must actually update the hoomd.hpmc.field.wall) managed object.
  • move_ratio (float) – the probability with which an update move is attempted
  • seed (int) – the seed of the pseudo-random number generator that determines whether or not an update move is attempted
  • period (int) – the number of timesteps between update move attempt attempts Every period steps, a walls update move is tried with probability move_ratio. This update move is provided by the py_updater callback. Then, update.wall only accepts an update move provided by the python callback if it maintains confinement conditions associated with all walls. Otherwise, it reverts back to a non-updated copy of the walls.

Once initialized, the update provides the following log quantities that can be logged via hoomd.analyze.log:

  • hpmc_wall_acceptance_ratio - the acceptance ratio for wall update moves

Example:

mc = hpmc.integrate.sphere(seed = 415236);
ext_wall = hpmc.compute.wall(mc);
ext_wall.add_sphere_wall(radius = 1.0, origin = [0, 0, 0], inside = True);
def perturb(timestep):
  r = np.sqrt(ext_wall.get_sphere_wall_param(index = 0, param = "rsq"));
  ext_wall.set_sphere_wall(index = 0, radius = 1.5*r, origin = [0, 0, 0], inside = True);
wall_updater = hpmc.update.wall(mc, ext_wall, perturb, move_ratio = 0.5, seed = 27, period = 50);
log = analyze.log(quantities=['hpmc_wall_acceptance_ratio'], period=100, filename='log.dat', overwrite=True);

Example:

mc = hpmc.integrate.sphere(seed = 415236);
ext_wall = hpmc.compute.wall(mc);
ext_wall.add_sphere_wall(radius = 1.0, origin = [0, 0, 0], inside = True);
def perturb(timestep):
  r = np.sqrt(ext_wall.get_sphere_wall_param(index = 0, param = "rsq"));
  ext_wall.set_sphere_wall(index = 0, radius = 1.5*r, origin = [0, 0, 0], inside = True);
wall_updater = hpmc.update.wall(mc, ext_wall, perturb, move_ratio = 0.5, seed = 27, period = 50);
disable()

Disables the updater.

Examples:

updater.disable()

Executing the disable command will remove the updater from the system. Any hoomd.run() command executed after disabling an updater will not use that updater during the simulation. A disabled updater can be re-enabled with enable()

enable()

Enables the updater.

Examples:

updater.enable()

See also

disable()

get_accepted_count(mode=0)

Get the number of accepted wall update moves.

Parameters:mode (int) – specify the type of count to return. If mode!=0, return absolute quantities. If mode=0, return quantities relative to the start of the run. DEFAULTS to 0.
Returns:the number of accepted wall update moves

Example:

mc = hpmc.integrate.sphere(seed = 415236);
ext_wall = hpmc.compute.wall(mc);
ext_wall.add_sphere_wall(radius = 1.0, origin = [0, 0, 0], inside = True);
def perturb(timestep):
  r = np.sqrt(ext_wall.get_sphere_wall_param(index = 0, param = "rsq"));
  ext_wall.set_sphere_wall(index = 0, radius = 1.5*r, origin = [0, 0, 0], inside = True);
wall_updater = hpmc.update.wall(mc, ext_wall, perturb, move_ratio = 0.5, seed = 27, period = 50);
run(100);
acc_count = wall_updater.get_accepted_count(mode = 0);
get_total_count(mode=0)

Get the number of attempted wall update moves.

Parameters:mode (int) – specify the type of count to return. If mode!=0, return absolute quantities. If mode=0, return quantities relative to the start of the run. DEFAULTS to 0.
Returns:the number of attempted wall update moves

Example:

mc = hpmc.integrate.sphere(seed = 415236);
ext_wall = hpmc.compute.wall(mc);
ext_wall.add_sphere_wall(radius = 1.0, origin = [0, 0, 0], inside = True);
def perturb(timestep):
  r = np.sqrt(ext_wall.get_sphere_wall_param(index = 0, param = "rsq"));
  ext_wall.set_sphere_wall(index = 0, radius = 1.5*r, origin = [0, 0, 0], inside = True);
wall_updater = hpmc.update.wall(mc, ext_wall, perturb, move_ratio = 0.5, seed = 27, period = 50);
run(100);
tot_count = wall_updater.get_total_count(mode = 0);
restore_state()

Restore the state information from the file used to initialize the simulations

set_period(period)

Changes the updater period.

Parameters:period (int) – New period to set.

Examples:

updater.set_period(100);
updater.set_period(1);

While the simulation is running, the action of each updater is executed every period time steps. Changing the period does not change the phase set when the analyzer was first created.