Randomizing the System

Overview

Questions

  • How can I generate a random initial condition?

Objectives

  • Show how to use HPMC to randomize the initial condition.

  • Demonstrate how to run a simulation.

  • Show how to use HPMC integrator properties to examine the acceptance ratio.

  • Explain that short simulations at low density effectively randomize the system.

Boilerplate code

[1]:
import hoomd
import math

The render function in the next (hidden) cell will render a snapshot using fresnel.

This is not intended as a full tutorial on fresnel - see the fresnel user documentation if you would like to learn more.

Method

The previous section of this tutorial placed all the particles on a simple cubic lattice. This is a convenient way to place non-overlapping particles, but it starts the simulation in a highly ordered state. You should randomize the the system enough so that it forgets this initial state and self-assembly can proceed without influence by the initial condition.

You cannot draw random numbers trivially for the particle positions, as that will result in overlaps between particles. Instead, start from the lattice and use HPMC to move particles randomly while ensuring that they do not overlap. In low density configurations, like the lattice generated in the previous section, a short simulation will quickly randomize the system.

Set up the simulation

The following code block creates the Simulation, configures the HPMC integrator, and initializes the system state from lattice.gsd as has been discussed in previous sections in this tutorial:

[3]:
cpu = hoomd.device.CPU()
sim = hoomd.Simulation(device=cpu)

mc = hoomd.hpmc.integrate.ConvexPolyhedron(seed=2)
mc.shape['octahedron'] = dict(vertices=[(-0.5, 0, 0),
                                         (0.5, 0, 0),
                                         (0, -0.5, 0),
                                         (0, 0.5, 0),
                                         (0, 0, -0.5),
                                         (0, 0, 0.5)])

sim.operations.integrator = mc
sim.create_state_from_gsd(filename='lattice.gsd')

Run the simulation

Save a snapshot of the current state of the system. This tutorial uses this later to see how far particles have moved.

[4]:
initial_snapshot = sim.state.snapshot
render(initial_snapshot)
[4]:
../../_images/tutorial_00-Introducing-HOOMD-blue_04-Randomizing-the-System_8_0.png

Run the simulation to randomize the particle positions and orientations. The run method takes the number of steps to run as an argument. 10,000 steps is enough to randomize a low density system:

[5]:
sim.run(10e3)

You can query properties of the HPMC integrator to see what it did. translate_moves is a tuple with the number of accepted and rejected translation moves. The acceptance ratio, the fraction of attempted moves which is accepted, is very high at this low density.

[6]:
mc.translate_moves
[6]:
(4487025, 510892)
[7]:
mc.translate_moves[0] / sum(mc.translate_moves)
[7]:
0.8977790147375396

rotate_moves similarly provides the number of accepted and rejected rotation moves.

[8]:
mc.rotate_moves
[8]:
(4753854, 248229)
[9]:
mc.rotate_moves[0] / sum(mc.rotate_moves)
[9]:
0.9503748738275635

overlaps reports the number of overlapping particle pairs in the state. There are no overlaps in the final configuration:

[10]:
mc.overlaps
[10]:
0

The final configuration

Look at the final particle positions and orientations and see how they have changed:

[11]:
final_snapshot = sim.state.snapshot
render(final_snapshot)
[11]:
../../_images/tutorial_00-Introducing-HOOMD-blue_04-Randomizing-the-System_20_0.png
[12]:
initial_snapshot.particles.position[0:4]
[12]:
array([[-3.8499999 , -3.8499999 , -3.8499999 ],
       [-3.8499999 , -3.8499999 , -2.75      ],
       [-3.8499999 , -3.8499999 , -1.64999998],
       [-3.8499999 , -3.8499999 , -0.55000001]])
[13]:
final_snapshot.particles.position[0:4]
[13]:
array([[ 3.5075011 ,  3.45149507, -0.23362428],
       [-0.34716783, -1.7078039 ,  1.89944438],
       [ 2.12477682, -0.39671483,  0.2136613 ],
       [-0.19598204,  1.5719255 ,  3.72206555]])
[14]:
initial_snapshot.particles.orientation[0:4]
[14]:
array([[1., 0., 0., 0.],
       [1., 0., 0., 0.],
       [1., 0., 0., 0.],
       [1., 0., 0., 0.]])
[15]:
final_snapshot.particles.orientation[0:4]
[15]:
array([[-0.91302374, -0.19004781, -0.36089446,  0.00496663],
       [ 0.49797178,  0.75066879,  0.39672083, -0.17644562],
       [ 0.23839869,  0.41443363, -0.87414225,  0.08535899],
       [ 0.1602925 , -0.45188123,  0.76177025, -0.43567851]])

The particle positions and orientations have indeed changed significantly, telling us that the system is well randomized.

Save the final configuration to a GSD file for use in the next stage of the simulation:

[16]:
hoomd.write.GSD.write(state=sim.state, filename='random.gsd')

The next section of the tutorial takes random.gsd and compresses it down to a higher density.