Overview
When performing surface science calculations or studying crystal structures with specific crystallographic orientations, creating properly oriented slab and bulk models is essential. This tutorial covers two complementary approaches: manual transformation using VESTA (covered in a previous article Crystal Structure Modeling with VESTA) and automated generation using ASE's powerful structure-building tools.
Creating Slab Structures with ASE
Basic Surface Generation
ASE provides convenient functions to directly generate slab structures exposing specific crystallographic facets. The following example creates an Al(111) surface slab:
import os
from ase.build import fcc111
from ase.io import write
# Create Al(111) surface slab
slab = fcc111('Al', size=(2, 2, 7), a=4.05, vacuum=15.0)
# Save structure
output_dir = '../1_structures/initial_structures'
os.makedirs(output_dir, exist_ok=True)
# Save as VASP format
write(f'{output_dir}/al_slab_111_2x2x7.vasp', slab, format='vasp')Sorting Atoms by z-Coordinate
Automatically generated structures often have atoms ordered in the opposite direction along the z-axis compared to conventional expectations. The following enhanced script sorts atoms from bottom to top (increasing z-coordinates):
"""Create aluminium (111) slab structure along z-axis with atoms sorted by z-coordinate"""
import os
import numpy as np
from ase import Atoms
from ase.build import fcc111
from ase.io import write
def sort_atoms_by_z(atoms, reverse=True):
"""
Sort atoms by z-coordinate
Parameters:
atoms: ASE Atoms object
reverse: True for descending (top to bottom), False for ascending
Returns:
Sorted Atoms object
"""
positions = atoms.get_positions()
z_coords = positions[:, 2]
# Get sorted indices
if reverse:
sorted_indices = np.argsort(z_coords)[::-1] # Descending
else:
sorted_indices = np.argsort(z_coords) # Ascending
# Rearrange atoms according to sorted indices
sorted_atoms = atoms[sorted_indices]
print(f"Atoms sorted by z-coordinate ({'descending' if reverse else 'ascending'})")
print(f"z range: {positions[:, 2].min():.3f} ~ {positions[:, 2].max():.3f} Å")
return sorted_atoms
# Create periodic bulk structure
slab_periodic = fcc111('Al', size=(2, 2, 7), a=4.05, vacuum=15.0)
# Sort by z-coordinate in descending order
slab_sorted = sort_atoms_by_z(slab_periodic, reverse=True)
# Save in different formats
output_dir = os.path.join(os.path.dirname(__file__), '..', '1_structures', 'initial_structures')
os.makedirs(output_dir, exist_ok=True)
# Save as VASP format
write(os.path.join(output_dir, 'al_slab_111_2x2x7.vasp'), slab_sorted, format='vasp')
# Save as CIF format
# write(os.path.join(output_dir, 'al_bulk_111_2x2x7.cif'), slab_sorted)Generating Bulk Structures with Specific Orientations
Creating Oriented Unit Cells
When creating bulk structures with a specific crystallographic plane oriented along the z-axis, ASE's slab-building functions can be adapted by setting vacuum=0.0. However, this approach generates duplicate atoms at the top and bottom surfaces due to periodic boundary conditions.
import os
import numpy as np
from ase import Atoms
from ase.build import fcc111
from ase.io import write
slab = fcc111('Al', size=(2, 2, 7), a=4.05, vacuum=0.0)
# Save in different formats
output_dir = os.path.join(os.path.dirname(__file__), '..', '1_structures', 'initial_structures')
os.makedirs(output_dir, exist_ok=True)
# Save as VASP format
write(os.path.join(output_dir, 'al_bulk_111_2x2x7.vasp'), slab, format='vasp')
# Save as CIF format
# write(os.path.join(output_dir, 'al_bulk_111_2x2x7.cif'), slab)Removing Periodic Duplicate Atoms
To create a proper bulk structure without redundant atoms, we can implement a function to detect and remove duplicates at periodic boundaries:
"""Create aluminium bulk structure with specific crystallographic plane along z-axis,
removing duplicate atoms caused by periodic boundary conditions"""
import os
import numpy as np
from ase import Atoms
from ase.build import fcc111
from ase.io import write
def remove_periodic_duplicates(atoms, tolerance=1e-5, keep='top'):
"""
Remove duplicate atoms caused by periodic boundary conditions
Parameters:
atoms: ASE Atoms object
tolerance: Distance threshold for identifying duplicates (Å)
keep: 'top' to retain top layer atoms, 'bottom' to retain bottom layer
Returns:
Atoms object with duplicates removed
"""
cell = atoms.get_cell()
positions = atoms.get_positions()
scaled_positions = atoms.get_scaled_positions()
# Track indices to remove
remove_indices = set()
for i in range(len(positions)):
if i in remove_indices:
continue
frac_pos_i = scaled_positions[i]
for j in range(i + 1, len(positions)):
if j in remove_indices:
continue
frac_pos_j = scaled_positions[j]
# Calculate minimum distance under periodic boundaries
diff = frac_pos_i - frac_pos_j
# Map difference to [-0.5, 0.5] range
diff = diff - np.round(diff)
# Convert back to Cartesian coordinates
diff_cart = np.dot(diff, cell)
distance = np.linalg.norm(diff_cart)
if distance < tolerance:
# Found duplicate, decide which to keep based on 'keep' parameter
if keep == 'top':
# Retain higher z-coordinate (top), remove lower (bottom)
if positions[i][2] > positions[j][2]:
remove_indices.add(j)
else:
remove_indices.add(i)
else: # keep == 'bottom'
# Retain lower z-coordinate (bottom), remove higher (top)
if positions[i][2] < positions[j][2]:
remove_indices.add(j)
else:
remove_indices.add(i)
# Create index list of atoms to retain
keep_indices = [i for i in range(len(positions)) if i not in remove_indices]
# Create new structure without duplicate atoms
new_atoms = atoms[keep_indices]
print(f"Original structure: {len(atoms)} atoms")
print(f"After removing duplicates: {len(new_atoms)} atoms")
print(f"Removed {len(atoms) - len(new_atoms)} duplicate atoms (retained {keep} layer)")
return new_atoms
# Create periodic bulk structure
slab_periodic = fcc111('Al', size=(1, 3, 7), a=4.05, vacuum=0.0)
# Remove periodic boundary duplicates (retain bottom, remove top)
slab = remove_periodic_duplicates(slab_periodic, keep='bottom')
# Save in different formats
output_dir = os.path.join(os.path.dirname(__file__), '..', '1_structures', 'initial_structures')
os.makedirs(output_dir, exist_ok=True)
# Save as VASP format
write(os.path.join(output_dir, 'al_bulk_111_1x3x7.vasp'), slab, format='vasp')
# Save as CIF format
write(os.path.join(output_dir, 'al_bulk_111_1x3x7.cif'), slab)Reference Aluminium Structures in VASP Format
FCC Unit Cell
Aluminium conventional FCC unit cell (filename: al4b-1.vasp):
al4b-1.vasp
1.0
4.0500000000 0.0000000000 0.0000000000
0.0000000000 4.0500000000 0.0000000000
0.0000000000 0.0000000000 4.0500000000
Al
4
Cartesian
0.000000000 0.000000000 0.000000000
0.000000000 2.025000000 2.025000000
2.025000000 0.000000000 2.025000000
2.025000000 2.025000000 0.000000000-Oriented Unit Cell
Aluminium unit cell with (001) plane along z-axis (filename: al2b-1.vasp):
al2b-1.vasp
1.0
2.8637826443 0.0000000000 0.0000000000
0.0000000000 2.8637826443 0.0000000000
0.0000000000 0.0000000000 4.0500000000
Al
2
Cartesian
0.000000000 0.000000000 0.000000000
1.431891322 1.431891322 2.025000000-Oriented Unit Cell
Aluminium unit cell with (110) plane along z-axis (filename: al2b-2.vasp):
al2b-2.vasp
1.0
2.8637826443 0.0000000000 0.0000000000
0.0000000000 4.0500000000 0.0000000000
0.0000000000 0.0000000000 2.8637826443
Al
2
Cartesian
0.000000000 0.000000000 0.000000000
1.431891322 2.025000000 1.431891322-Oriented Structures
Orthorhombic Unit Cell
Aluminium orthorhombic unit cell with (111) plane along z-axis (filename: al6b-1.vasp):
al6b-1.vasp
1.0
4.9602169991 0.0000000000 0.0000000000
0.0000000000 2.8637826443 0.0000000000
0.0000000000 0.0000000000 7.0148062706
Al
6
Cartesian
0.000000000 0.000000000 0.000000000
2.480108500 1.431891322 0.000000000
1.653405716 0.000000000 2.338268827
4.133514067 1.431891322 2.338268827
0.826702858 1.431891322 4.676537653
3.306811431 0.000000000 4.676537653Hexagonal Unit Cell (3 atoms)
Hexagonal unit cell generated via ASE (filename: al3b-1.vasp):
al3b-1
1.0
2.8637824059 0.0000000000 0.0000000000
1.4318912029 2.4801083144 0.0000000000
0.0000000000 0.0000000000 7.0148057938
Al
3
Cartesian
0.000000000 0.000000000 0.000000000
1.431891118 0.826702722 2.338268459
2.863782406 1.653405444 4.676536917Hexagonal Unit Cell (6 atoms)
Extended hexagonal unit cell generated via ASE (filename: al6b-3.vasp):
Al
1.0
2.8637824059 0.0000000000 0.0000000000
2.8637824059 4.9602166288 0.0000000000
0.0000000000 0.0000000000 7.0148057938
Al
6
Cartesian
0.000000000 0.000000000 0.000000000
1.431891118 2.480108167 0.000000000
1.431891118 0.826702722 2.338268459
2.863782235 3.306810889 2.338268459
2.863782406 1.653405444 4.676536917
4.295673438 4.133513759 4.676536917