Initializing the System State

Overview

Questions

  • How do I place particles in the initial condition?

  • What units does HOOMD-blue use?

Objectives

  • Describe the basic properties of the initial condition, including particle position, orientation, type and the periodic box.

  • Briefly describe HOOMD-blues system of units.

  • Demonstrate writing a system to a GSD file.

  • Show how to initialize Simulation state from a GSD file.

Boilerplate code

[1]:
import hoomd
import math
import itertools
import numpy

The render function in the next (hidden) cell will render the initial condition using fresnel.

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

Components of the system state

You need to initialize the system state before you can run a simulation. The initial condition describes the the position and orientation of every particle in the system and the periodic box at the start of the simulation.

The hard regular octahedra system self-assembles the the bcc structure. Self-assembly is a process where particles will organize themselves into an ordered structure at equilibrium. Most self-assembly studies run simulations of many thousands of particles for tens of hours. To keep this tutorial short, it simulations a small number of particles commensurate with the bcc structure (2 * m**3, where m is an integer).

[3]:
m = 5
N_particles = 2 * m**3

Placing particles

In hard particle Monte Carlo, valid particle configurations have no overlaps. The octahedron particle in this tutorial sits inside a sphere of diameter 1, so place particles a little bit further than that apart on a KxKxK simple cubic lattice of width L. The following sections of this tutorial will demonstrate how to randomize this ordered configuration and compress it to a higher density.

[4]:
spacing = 1.1
K = math.ceil(N_particles**(1/3))
L = K * spacing

In HOOMD, particle positions must be placed inside a period box. Cubic boxes range from -L/2 to L/2. Use itertools and numpy to place the positions on the lattice in this range.

[5]:
x = numpy.linspace(-L / 2, L / 2, K, endpoint=False)
position = list(itertools.product(x, repeat=3))
print(position[0:4])
[(-3.8500000000000005, -3.8500000000000005, -3.8500000000000005), (-3.8500000000000005, -3.8500000000000005, -2.7500000000000004), (-3.8500000000000005, -3.8500000000000005, -1.6500000000000004), (-3.8500000000000005, -3.8500000000000005, -0.5500000000000003)]

Filter this list down to N particles because K*K*K >= N_particles .

[6]:
position = position[0:N_particles]

You also need to set the orientation for each particle. HOOMD represents orientations with quaternions. The quaternion (1, 0, 0, 0) represents no rotation. Set that for all particles:

[7]:
orientation = [(1, 0, 0, 0)]*N_particles

Here is what the system looks like now:

[8]:
render(position, orientation, L)
[8]:
../../_images/tutorial_00-Introducing-HOOMD-blue_03-Initializing-the-System-State_15_0.png

Units

HOOMD-blue does not adopt any particular real system of units. Instead, HOOMD-blue uses an internally self-consistent system of units and is compatible with many systems of units. For example: if you select the units of meter, Joule, and kilogram for length, energy and mass then the units of force will be Newtons and velocity will be meters/second. A popular system of units for nano-scale systems is nanometers, kilojoules/mol, and atomic mass units.

In HPMC, the primary unit is that of length. Mass exists, but is factored out of the partition function and does not enter into the simulation. The scale of energy is irrelevant in athermal HPMC systems where overlapping energies are infinite and valid configurations have zero potential energy. However, energy does appear implicitly in derived units like \([\mathrm{pressure}] = \left(\frac{\mathrm{[energy]}}{\mathrm{[length]}^3}\right)\). In HPMC, \(kT\) is assumed to be 1 energy. You can convert a pressure \(P\) in HPMC’s units to a non-dimensional value with \(P^* = \frac{P\sigma^3}{kT}\) where \(\sigma\) is a representative length in your system, such as the particle’s diameter or edge length.

Writing the configuration to the file system

GSD files store the periodic box, particle positions, orientations, and other properties of the state. Use the GSD Python package to write this file.

[9]:
import gsd.hoomd

The Snapshot object stores the state of the system.

[10]:
snapshot = gsd.hoomd.Snapshot()
snapshot.particles.N = N_particles
snapshot.particles.position = position
snapshot.particles.orientation = orientation

Each particle also has a type. In HPMC, each type has its own shape parameters. There is a single type in this tutorial, so set a type id of 0 for all particles (type ids are 0-indexed):

[11]:
snapshot.particles.typeid = [0]*N_particles

Every particle type needs a name. Names can be any string. HOOMD-blue uses the type names to match the parameters you specify in operations with the type names present in the state. Give the name 'octahedron' to the type id 0:

[12]:
snapshot.particles.types = ['octahedron']

GSD represents boxes with a 6-element array. Three box lengths L_x, L_y, L_z, and 3 tilt factors. Set equal box lengths 0 tilt factors to define a cubic box.

[13]:
snapshot.configuration.box = [L, L, L, 0, 0, 0]

Write this snapshot to lattice.gsd:

[14]:
with gsd.hoomd.open(name='lattice.gsd', mode='wb') as f:
    f.append(snapshot)

Initializing a Simulation

You can use the file to initialize the Simulation state. First, create a Simulation:

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

Then use the create_state_from_gsd factory function to read in the GSD file and populate the state with the particles and box from the file:

[16]:
sim.create_state_from_gsd(filename='lattice.gsd')

The next section in this tutorial will show you how to use HPMC to randomize the particle positions and orientations.