Compare shapes¶
It is sometimes of interest to compare nanoparticles of different structural motifs. This example illustrates how truncated octahedral, decahedral, and icosahedral particles constructed based on the same surface energies can be compared with WulffPack.
Preparations¶
Before creating the particles, some parameters should be chosen (the values here are arbitrary). For simplicity, we will only consider {111} and {100} facets, but WulffPack can handle an arbitrary number of facets (although the code may take a fairly long time to run). When doing quantitative analyses of the particles, it is recommended to use electronvolts as unit for energy and Ångström as unit for distance.
from wulffpack import (SingleCrystal,
Decahedron,
Icosahedron)
from ase.build import bulk
# Specify parameters
prim = bulk('Au',
crystalstructure='fcc',
a=4.08)
surface_energies = {(1, 0, 0): 0.066,
(1, 1, 1): 0.06}
twin_energy = 1e-3
natoms = 1000
We then create three particles, one of every kind.
oh = SingleCrystal(surface_energies,
primitive_structure=prim,
natoms=natoms)
dh = Decahedron(surface_energies,
twin_energy=twin_energy,
primitive_structure=prim,
natoms=natoms)
ih = Icosahedron(surface_energies,
twin_energy=twin_energy,
primitive_structure=prim,
natoms=natoms)
Volume¶
We may now compare the three particles. For peace of mind, we first check that they have the same volume. In WulffPack, the volume is determined from the specified number of atoms and the primitive structure together. It does not take into account that the specified number of atoms may be inadequate to create the specified shape.
print('Volume (sanity check):')
for name, particle in zip(['Oh', 'Dh', 'Ih'], [oh, dh, ih]):
volume = particle.volume
print('{}: {:.2f}'.format(name, volume))
print()
This snippet produces the following output:
Volume (sanity check):
Oh: 16979.33
Dh: 16979.33
Ih: 16979.33
Area¶
We may now turn to something more interesting: the surface area. The syntax is the same as for the volume:
print('Total surface area (excluding twin boundaries):')
for name, particle in zip(['Oh', 'Dh', 'Ih'], [oh, dh, ih]):
area = particle.area
print('{}: {:.2f}'.format(name, area))
print()
Total surface area (excluding twin boundaries):
Oh: 3488.01
Dh: 3468.81
Ih: 3400.92
As one might have expected, the surface area of the nearly spherical icosahedron is the smallest.
Fraction of a specified facet¶
It can also be of interest to compare how large fraction of the area that is, say, {100}:
print('Fraction of {100} facets:')
for name, particle in zip(['Oh', 'Dh', 'Ih'], [oh, dh, ih]):
fraction = particle.facet_fractions.get((1, 0, 0), 0.)
print('{}: {:.4f}'.format(name, fraction))
print()
Fraction of {100} facets:
Oh: 0.2775
Dh: 0.2373
Ih: 0.0000
The truncated octahedron has the largest fraction of {100} facets, which is the primary reason decahedra and icosahedra are often energetically favorable for small nanoparticles made of FCC metals, for which the {111} energy is typically lower than {100}.
Surface energy¶
We may compare the total surface energy of the three particles, this time including the energy of the twin boundaries in the decahedral and icosahedral particles:
print('Total surface energy (including twin energy):')
for name, particle in zip(['Oh', 'Dh', 'Ih'], [oh, dh, ih]):
energy = particle.surface_energy
print('{}: {:.2f}'.format(name, energy))
print()
Total surface energy (including twin energy):
Oh: 215.09
Dh: 215.07
Ih: 208.82
The output indicates that the icosahedron is the most stable shape, followed by the decahedron. Note, however, that this model only includes surface energies (including twin boundary energies) and no contributions from, e.g., strain (see below).
Average surface energy¶
It is sometimes of interest to extract the average surface energy of the particle, i.e. the average of the surface energy per facet \(E_f\) weighted with the area fractions of each facet \(f\),
where \(A_f\) is the area of facet \(f\) and \(A_\text{tot}\) is the total area of the particle.
print('Average surface energy (excluding twin energy):')
for name, particle in zip(['Oh', 'Dh', 'Ih'], [oh, dh, ih]):
energy = particle.average_surface_energy
print('{}: {:.4f}'.format(name, energy))
print()
Average surface energy (excluding twin energy):
Oh: 0.0617
Dh: 0.0614
Ih: 0.0600
This quantity is of particular interest in the single crystalline (Oh) case, since it provides an estimate for the surface energy that would be measured experimentally.
Strain energy¶
Decahedral and icosahedral particles can be built from FCC crystalline grain only if they are strained. This strain energy can be estimated based on continuum mechanics if we know the shear modulus and Poisson’s ratio of the material [HowMar84].
shear_modulus = 0.17 # eV / Angstrom^3
poissons_ratio = 0.42 # unitless
print('Estimated strain energies:')
for name, particle in zip(['Dh', 'Ih'], [dh, ih]):
energy = particle.get_strain_energy(shear_modulus=shear_modulus,
poissons_ratio=poissons_ratio)
print('{}: {:.2f}'.format(name, energy))
print()
Estimated strain energies:
Dh: 0.52
Ih: 5.94
We note that for this size, the strain energy is sufficient to make the decahedron less stable than the truncated octahedron (the sum of total surface energy and strain energy is larger for the decahedron). Note that these strain energies should only be taken as fairly crude approximations.
Number of corners and total edge length¶
Finally, it can be interesting to compare how many corners (vertices) the particles have, as well as the total length of their edges. Note that these only includes corners and edges on the surface of the particle, not the 1D or 2D defects inside decahedra and icosahedra.
print('Number of corners and total edge length:')
for name, particle in zip(['Oh', 'Dh', 'Ih'], [oh, dh, ih]):
ncorners = particle.number_of_corners
edge_length = particle.edge_length
print('{}: {} corners, edge length: {:.2f}'.format(name,
ncorners,
edge_length))
Number of corners and total edge length:
Oh: 24 corners, edge length: 409.48
Dh: 42 corners, edge length: 610.03
Ih: 72 corners, edge length: 614.70
Note that for both the decahedron and the icosahedron, there are very small “reentrant” surfaces close to where the grains meet at the fivefold axis. These reentrant surfaces increase the number of corners significantly, but are unlikely to manifest in an atomistic representation of these shapes. One can in principle remove these facets by specifying a very small value for the twin energy.
Source code¶
The complete source code is available in
examples/compare_shapes.py
from wulffpack import (SingleCrystal,
Decahedron,
Icosahedron)
from ase.build import bulk
# Specify parameters
prim = bulk('Au',
crystalstructure='fcc',
a=4.08)
surface_energies = {(1, 0, 0): 0.066,
(1, 1, 1): 0.06}
twin_energy = 1e-3
natoms = 1000
# Create one particle per type
oh = SingleCrystal(surface_energies,
primitive_structure=prim,
natoms=natoms)
dh = Decahedron(surface_energies,
twin_energy=twin_energy,
primitive_structure=prim,
natoms=natoms)
ih = Icosahedron(surface_energies,
twin_energy=twin_energy,
primitive_structure=prim,
natoms=natoms)
# Compare volumes
print('Volume (sanity check):')
for name, particle in zip(['Oh', 'Dh', 'Ih'], [oh, dh, ih]):
volume = particle.volume
print('{}: {:.2f}'.format(name, volume))
print()
# Compare areas
print('Total surface area (excluding twin boundaries):')
for name, particle in zip(['Oh', 'Dh', 'Ih'], [oh, dh, ih]):
area = particle.area
print('{}: {:.2f}'.format(name, area))
print()
# Compare fractions of {100} facets
print('Fraction of {100} facets:')
for name, particle in zip(['Oh', 'Dh', 'Ih'], [oh, dh, ih]):
fraction = particle.facet_fractions.get((1, 0, 0), 0.)
print('{}: {:.4f}'.format(name, fraction))
print()
# Compare total surface energies
print('Total surface energy (including twin energy):')
for name, particle in zip(['Oh', 'Dh', 'Ih'], [oh, dh, ih]):
energy = particle.surface_energy
print('{}: {:.2f}'.format(name, energy))
print()
# Compare average surface energies
print('Average surface energy (excluding twin energy):')
for name, particle in zip(['Oh', 'Dh', 'Ih'], [oh, dh, ih]):
energy = particle.average_surface_energy
print('{}: {:.4f}'.format(name, energy))
print()
# Compare strain energies
shear_modulus = 0.17 # eV / Angstrom^3
poissons_ratio = 0.42 # unitless
print('Estimated strain energies:')
for name, particle in zip(['Dh', 'Ih'], [dh, ih]):
energy = particle.get_strain_energy(shear_modulus=shear_modulus,
poissons_ratio=poissons_ratio)
print('{}: {:.2f}'.format(name, energy))
print()
# Compare edge length and number of corners
print('Number of corners and total edge length:')
for name, particle in zip(['Oh', 'Dh', 'Ih'], [oh, dh, ih]):
ncorners = particle.number_of_corners
edge_length = particle.edge_length
print('{}: {} corners, edge length: {:.2f}'.format(name,
ncorners,
edge_length))