hoomd.update

Overview

BoxResize

Resizes the box between an initial and final box.

CustomUpdater

User-defined updater.

FilterUpdater

Update sets of particles associated with a filter.

RemoveDrift

Remove the average drift from a restrained system.

Details

Updaters.

Updater operations modify the simulation state when they act.

class hoomd.update.BoxResize(trigger, box1, box2, variant, filter=hoomd._hoomd.ParticleFilter)

Bases: Updater

Resizes the box between an initial and final box.

BoxResize resizes the box between gradually from the initial box to the final box. The simulation box follows the linear interpolation between the initial and final boxes where the minimum of the variant gives box1 and the maximum gives box2:

\[\begin{split}\begin{align*} L_{x}' &= \lambda L_{2x} + (1 - \lambda) L_{1x} \\ L_{y}' &= \lambda L_{2y} + (1 - \lambda) L_{1y} \\ L_{z}' &= \lambda L_{2z} + (1 - \lambda) L_{1z} \\ xy' &= \lambda xy_{2} + (1 - \lambda) xy_{2} \\ xz' &= \lambda xz_{2} + (1 - \lambda) xz_{2} \\ yz' &= \lambda yz_{2} + (1 - \lambda) yz_{2} \\ \end{align*}\end{split}\]

Where box1 is \((L_{1x}, L_{1y}, L_{1z}, xy_1, xz_1, yz_1)\), box2 is \((L_{2x}, L_{2y}, L_{2z}, xy_2, xz_2, yz_2)\), \(\lambda = \frac{f(t) - \min f}{\max f - \min f}\), \(t\) is the timestep, and \(f(t)\) is given by variant.

For each particle \(i\) matched by filter, BoxResize scales the particle to fit in the new box:

\[\vec{r}_i \leftarrow s_x \vec{a}_1' + s_y \vec{a}_2' + s_z \vec{a}_3' - \frac{\vec{a}_1' + \vec{a}_2' + \vec{a}_3'}{2}\]

where \(\vec{a}_k'\) are the new box vectors determined by \((L_x', L_y', L_z', xy', xz', yz')\) and the scale factors are determined by the current particle position \(\vec{r}_i\) and the old box vectors \(\vec{a}_k\):

\[\vec{r}_i = s_x \vec{a}_1 + s_y \vec{a}_2 + s_z \vec{a}_3 - \frac{\vec{a}_1 + \vec{a}_2 + \vec{a}_3}{2}\]

After scaling particles that match the filter, BoxResize wraps all particles \(j\) back into the new box:

\[\vec{r_j} \leftarrow \mathrm{minimum\_image}_{\vec{a}_k}' (\vec{r}_j)\]

Important

The passed Variant must be bounded on the interval \(t \in [0,\infty)\) or the behavior of the updater is undefined.

Warning

Rescaling particles fails in HPMC simulations with more than one MPI rank.

Note

When using rigid bodies, ensure that the BoxResize updater is last in the operations updater list. Immediately after the BoxResize updater triggers, rigid bodies (hoomd.md.constrain.Rigid) will be temporarily deformed. hoomd.md.Integrator will run after the last updater and resets the constituent particle positions before computing forces.

Parameters
box1

The box associated with the minimum of the passed variant.

Type

hoomd.Box

box2

The box associated with the maximum of the passed variant.

Type

hoomd.Box

variant

A variant used to interpolate between the two boxes.

Type

hoomd.variant.Variant

trigger

The trigger to activate this updater.

Type

hoomd.trigger.Trigger

filter

The subset of particles to update.

Type

hoomd.filter.filter_like

get_box(timestep)

Get the box for a given timestep.

Parameters

timestep (int) – The timestep to use for determining the resized box.

Returns

The box used at the given timestep. None before the first call to Simulation.run.

Return type

Box

static update(state, box, filter=hoomd._hoomd.ParticleFilter)

Immediately scale the particle in the system state to the given box.

Parameters
class hoomd.update.CustomUpdater(trigger, action)

Bases: CustomOperation, Updater

User-defined updater.

Parameters

CustomUpdater is a hoomd.operation.Updater that wraps a user-defined hoomd.custom.Action object so the action can be added to a hoomd.Operations instance for use with hoomd.Simulation objects.

Updaters modify the system state.

class hoomd.update.FilterUpdater(trigger, filters)

Bases: Updater

Update sets of particles associated with a filter.

hoomd.Simulation caches the particles selected by hoomd.filter.filter_like objects to avoid the cost of re-running the filter on every time step. The particles selected by a filter will remain static until recomputed. This class provides a mechanism to update the cached list of particles. For example, periodically update a MD integration method’s group so that the integration method applies to particles in a given region of space.

Tip

To improve performance, use a hoomd.trigger.Trigger subclass, to update only when there is a known change to the particles that a filter would select.

Note

Some actions automatically recompute all filter particles such as adding or removing particles to the hoomd.Simulation.state.

Parameters
trigger

The trigger associated with the updater.

Type

hoomd.trigger.Trigger

__eq__(other)

Return whether two objects are equivalent.

property filters

filters to update select particles.

Type

list[hoomd.filter.filter_like]

class hoomd.update.RemoveDrift(reference_positions, trigger=1)

Bases: Updater

Remove the average drift from a restrained system.

Parameters

RemoveDrift computes the mean drift \(\vec{D}\) from the given reference_positions (\(\vec{r}_{ref, i}\)):

\[\vec{D} = \frac{1}{\mathrm{N_{particles}}} \sum_{i=0}^\mathrm{N_{particles-1}} \mathrm{minimum\_image}(\vec{r}_i - \vec{r}_{ref,i})\]

RemoveDrift then shifts all particles in the system by \(-\vec{D}\):

\[\vec{r}_i \leftarrow \mathrm{minimum\_image}(\vec{r}_i - \vec{D})\]

Tip

Use RemoveDrift with hoomd.hpmc.external.field.Harmonic to improve the accuracy of Frenkel-Ladd calculations.

reference_positions

the reference positions \([\mathrm{length}]\).

Type

(N_particles, 3) numpy.ndarray of float

trigger

The timesteps to remove drift.

Type

hoomd.trigger.Trigger