hoomd.update
Overview
Resizes the box between an initial and final box. |
|
User-defined updater. |
|
Update sets of particles associated with a filter. |
|
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 givesbox1
and the maximum givesbox2
:\[\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 byvariant
.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 theBoxResize
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
trigger (hoomd.trigger.trigger_like) – The trigger to activate this updater.
box1 (hoomd.box.box_like) – The box associated with the minimum of the passed variant.
box2 (hoomd.box.box_like) – The box associated with the maximum of the passed variant.
variant (hoomd.variant.variant_like) – A variant used to interpolate between the two boxes.
filter (hoomd.filter.filter_like) – The subset of particle positions to update.
- variant
A variant used to interpolate between the two boxes.
- trigger
The trigger to activate this updater.
- filter
The subset of particles to update.
- 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 toSimulation.run
.- Return type
- static update(state, box, filter=hoomd._hoomd.ParticleFilter)
Immediately scale the particle in the system state to the given box.
- Parameters
state (State) – System state to scale.
box (hoomd.box.box_like) – New box.
filter (hoomd.filter.filter_like) – The subset of particles to update.
- class hoomd.update.CustomUpdater(trigger, action)
Bases:
CustomOperation
,Updater
User-defined updater.
- Parameters
action (hoomd.custom.Action) – The action to call.
trigger (hoomd.trigger.trigger_like) – Select the timesteps to call the action.
CustomUpdater
is ahoomd.operation.Updater
that wraps a user-definedhoomd.custom.Action
object so the action can be added to ahoomd.Operations
instance for use withhoomd.Simulation
objects.Updaters modify the system state.
See also
The base class
hoomd.custom.CustomOperation
.
- class hoomd.update.FilterUpdater(trigger, filters)
Bases:
Updater
Update sets of particles associated with a filter.
hoomd.Simulation
caches the particles selected byhoomd.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 (hoomd.trigger.trigger_like) – A trigger to use for determining when to update particles associated with a filter.
filters (list[hoomd.filter.filter_like]) – A list of
hoomd.filter.filter_like
objects to update.
- trigger
The trigger associated with the updater.
- __eq__(other)
Return whether two objects are equivalent.
- property filters
filters to update select particles.
- Type
- class hoomd.update.RemoveDrift(reference_positions, trigger=1)
Bases:
Updater
Remove the average drift from a restrained system.
- Parameters
reference_positions ((N_particles, 3)
numpy.ndarray
offloat
) – the reference positions \([\mathrm{length}]\).trigger (hoomd.trigger.trigger_like) – Select the timesteps to remove drift.
RemoveDrift
computes the mean drift \(\vec{D}\) from the givenreference_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
withhoomd.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
offloat
- trigger
The timesteps to remove drift.