Basic usage

The simplest usages of WulffPack are illustrated here. Please note that the surface energies used in the examples are arbitrary.

Modules to import

To run the following examples, three WulffPack classes are imported (SingleCrystal, Decahedron, and Icosahedron). In addition, two functions from ASE will be used, ase.build.bulk() and ase.io.write().

from wulffpack import (SingleCrystal,
                       Decahedron,
                       Icosahedron)
from ase.build import bulk
from ase.io import write

Regular Wulff construction

In the simplest case, the only input WulffPack needs to make a Wulff construction is one or several surface energies. Surface energies are specified as a dictionary, where the keys are the Miller indices as a tuple and the values are surface energies as floats.

WulffPack does not impose any units on surface energies but uses only their ratio to determine the shape. This means that the surface energies can be specified in any units of choice, such as J/m^2 or eV/Å^2. Note, however, that ASE, which is used to define the primitive structure in WulffPack, uses Ångström, meaning that eV/Å^2 is recommended to obtain a sensible total surface energy from the Wulff construction.

surface_energies = {(1, 0, 0): 1.1,
                    (1, 1, 1): 1.,
                    (1, 1, 0): 1.15,
                    (2, 1, 1): 1.14}
particle = SingleCrystal(surface_energies)
particle.view()
write('single_crystal.xyz', particle.atoms)

The last two lines show a rendition of the particle with the aid of matplotlib and save an xyz file with an atomistic version of the structure. By default, the particle will be carved out from an FCC crystal and have roughly 1000 atoms. The crystal structure and the number of atoms can both be modified, as demonstrated below.

Regular Wulff construction for an HCP crystal

WulffPack allows Wulff constructions based on arbitrary crystals. To use a crystal other than FCC, one specifies a primitive structure in the form of an ASE Atoms object. In the below example, a Wulff construction based on an hexagonal close-packed (HCP) crystal is created. Bravais-Miller indices are used to denote facets.

prim = bulk('Co', crystalstructure='hcp')
surface_energies = {(1, 0, -1, 0): 1.1,
                    (0, 0, 0, 1): 1.0,
                    (1, 1, -2, 0): 1.0}
particle = SingleCrystal(surface_energies,
                         primitive_structure=prim,
                         natoms=5000)
particle.view()
write('hcp_single_crystal.xyz', particle.atoms)

Note

WulffPack uses Spglib to standardize the input primitive structure. The Miller (or Bravais-Miller) indices should always refer to this standardized version of the input primitive structure. In the case of FCC, for example, the Miller indices should refer to the conventional, cubic cell even if WulffPack is provided with a primitive FCC cell. The conventions for the standardization are detailed in the Spglib documentation. If unsure, the standardized cell is a property of the SingleCrystal object:

from ase.visualize import view
view(particle.standardized_structure)

The Miller indices should thus always refer to the reciprocal cell corresponding to the cell of that standardized structure.

If a larger atoms object with the same shape is wanted, the quickest solution is to simply set the natoms attribute of the same particle.

particle.natoms = 10000
write('hcp_single_crystal_larger.xyz', particle.atoms)

By doing so, all the coordinates are scaled such that the volume matches the new number of atoms.

Wulff construction of twinned particles

WulffPack can also construct decahedral and icosahedral particles using the so-called modified Wulff construction [Mar83]. In this case a cubic crystal is required, and FCC is probably the only crystal structure that makes physical sense. For this construction, it is also mandatory to specify a twin boundary energy.

surface_energies = {(1, 0, 0): 1.1,
                    (1, 1, 1): 1.,
                    (1, 1, 0): 1.15}
prim = bulk('Pd', a=3.9)
particle = Decahedron(surface_energies,
                      twin_energy=0.04,
                      primitive_structure=prim)
particle.view()
write('decahedron.xyz', particle.atoms)

The interface to construct an icosahedral particle is identical.

particle = Icosahedron(surface_energies,
                       twin_energy=0.04,
                       primitive_structure=prim)
particle.view()
write('icosahedron.xyz', particle.atoms)

Source code

The complete source code is available in examples/basic_usage.py

from wulffpack import (SingleCrystal,
                       Decahedron,
                       Icosahedron)
from ase.build import bulk
from ase.io import write

# Show a regular Wulff construction, cubic crystal
surface_energies = {(1, 0, 0): 1.1,
                    (1, 1, 1): 1.,
                    (1, 1, 0): 1.15,
                    (2, 1, 1): 1.14}
particle = SingleCrystal(surface_energies)
particle.view()
write('single_crystal.xyz', particle.atoms)

# Wulff construction for hcp crystal
prim = bulk('Co', crystalstructure='hcp')
surface_energies = {(1, 0, -1, 0): 1.1,
                    (0, 0, 0, 1): 1.0,
                    (1, 1, -2, 0): 1.0}
particle = SingleCrystal(surface_energies,
                         primitive_structure=prim,
                         natoms=5000)
particle.view()
write('hcp_single_crystal.xyz', particle.atoms)

# Change number of atoms
particle.natoms = 10000
write('hcp_single_crystal_larger.xyz', particle.atoms)

# Wulff construction for decahedron
surface_energies = {(1, 0, 0): 1.1,
                    (1, 1, 1): 1.,
                    (1, 1, 0): 1.15}
prim = bulk('Pd', a=3.9)
particle = Decahedron(surface_energies,
                      twin_energy=0.04,
                      primitive_structure=prim)
particle.view()
write('decahedron.xyz', particle.atoms)

# Wulff construction for icosahedron
particle = Icosahedron(surface_energies,
                       twin_energy=0.04,
                       primitive_structure=prim)
particle.view()
write('icosahedron.xyz', particle.atoms)