Crystal Builder Tutorial¶

This tutorial demonstrates how to use the molpy.builder.crystal module to generate crystal structures in a LAMMPS-like style.

"What you'll learn" By the end of this tutorial, you'll be able to: - Create crystal lattices with predefined or custom structures - Define regions in lattice or Cartesian coordinates - Generate large crystal structures efficiently using vectorized operations - Work with generated structures in the MolPy ecosystem

Overview¶

The crystal builder module provides:

Lattice definition: Define Bravais lattices with basis sites
Predefined structures: SC, BCC, FCC, rocksalt (NaCl)
Region control: Define regions in lattice or Cartesian coordinates
Vectorized generation: Efficient NumPy-based atom creation (no Python loops)

"LAMMPS-like workflow" The workflow is similar to LAMMPS: 1. Define a lattice → 2. Define a region → 3. Create atoms in the region

In [1]:

Copied!

from collections import Counter

import numpy as np

from molpy.builder.crystal import BlockRegion, CrystalBuilder, Lattice, Site
from collections import Counter

import numpy as np

from molpy.builder.crystal import BlockRegion, CrystalBuilder, Lattice, Site

1. Simple Cubic (SC) Structure¶

The simplest crystal structure with one atom per unit cell.

Simple Cubic Simple cubic (SC) is the most basic crystal structure where atoms are located only at the corners of a cube. Each unit cell contains exactly one atom.

In [2]:

Copied!





# Create a simple cubic lattice with lattice constant a = 2.0 Å
lat_sc = Lattice.cubic_sc(a=2.0, species="Cu")

print("Lattice vectors (3x3 matrix):")
print(lat_sc.cell)
print(f"\nNumber of basis sites: {len(lat_sc.basis)}")
print(f"Basis site fractional coordinates: {lat_sc.basis[0].frac}")
# Create a simple cubic lattice with lattice constant a = 2.0 Å
lat_sc = Lattice.cubic_sc(a=2.0, species="Cu")

print("Lattice vectors (3x3 matrix):")
print(lat_sc.cell)
print(f"\nNumber of basis sites: {len(lat_sc.basis)}")
print(f"Basis site fractional coordinates: {lat_sc.basis[0].frac}")

Lattice vectors (3x3 matrix):
[[2. 0. 0.]
 [0. 2. 0.]
 [0. 0. 2.]]

Number of basis sites: 1
Basis site fractional coordinates: (0.0, 0.0, 0.0)

In [3]:

Copied!





# Define a region in lattice coordinates
# This creates a 5×5×5 block of unit cells
region_sc = BlockRegion(0, 5, 0, 5, 0, 5, coord_system="lattice")

# Build the crystal structure
builder_sc = CrystalBuilder(lat_sc)
structure_sc = builder_sc.build_block(region_sc)

print(f"Generated {len(list(structure_sc.atoms))} atoms")
print(f"Expected: 5x5x5 = {5 * 5 * 5} atoms")
print("\nFirst 3 atom positions (in Å):")
print(structure_sc.xyz[:3])
# Define a region in lattice coordinates
# This creates a 5×5×5 block of unit cells
region_sc = BlockRegion(0, 5, 0, 5, 0, 5, coord_system="lattice")

# Build the crystal structure
builder_sc = CrystalBuilder(lat_sc)
structure_sc = builder_sc.build_block(region_sc)

print(f"Generated {len(list(structure_sc.atoms))} atoms")
print(f"Expected: 5x5x5 = {5 * 5 * 5} atoms")
print("\nFirst 3 atom positions (in Å):")
print(structure_sc.xyz[:3])

Generated 125 atoms
Expected: 5x5x5 = 125 atoms

First 3 atom positions (in Å):
[[0. 0. 0.]
 [0. 0. 2.]
 [0. 0. 4.]]

2. Face-Centered Cubic (FCC) Structure¶

FCC has four atoms per unit cell. Many metals like Ni, Cu, Al crystallize in FCC structure.

In [4]:

Copied!





# Create FCC Ni lattice (lattice constant for Ni is 3.52 Å)
lat_fcc = Lattice.cubic_fcc(a=3.52, species="Ni")

print(f"Number of basis sites: {len(lat_fcc.basis)}")
print("Basis site positions (fractional coordinates):")
for i, site in enumerate(lat_fcc.basis):
    print(f"  Site {i} ({site.label}): {site.frac}")

# Build FCC structure: 4×4×4 unit cells
region_fcc = BlockRegion(0, 4, 0, 4, 0, 4, coord_system="lattice")
builder_fcc = CrystalBuilder(lat_fcc)
structure_fcc = builder_fcc.build_block(region_fcc)

print(f"\nGenerated {len(list(structure_fcc.atoms))} atoms")
print(f"Expected: 4x4x4 (cells) x 4 (basis sites) = {4 * 4 * 4 * 4} atoms")
# Create FCC Ni lattice (lattice constant for Ni is 3.52 Å)
lat_fcc = Lattice.cubic_fcc(a=3.52, species="Ni")

print(f"Number of basis sites: {len(lat_fcc.basis)}")
print("Basis site positions (fractional coordinates):")
for i, site in enumerate(lat_fcc.basis):
    print(f"  Site {i} ({site.label}): {site.frac}")

# Build FCC structure: 4×4×4 unit cells
region_fcc = BlockRegion(0, 4, 0, 4, 0, 4, coord_system="lattice")
builder_fcc = CrystalBuilder(lat_fcc)
structure_fcc = builder_fcc.build_block(region_fcc)

print(f"\nGenerated {len(list(structure_fcc.atoms))} atoms")
print(f"Expected: 4x4x4 (cells) x 4 (basis sites) = {4 * 4 * 4 * 4} atoms")

Number of basis sites: 4
Basis site positions (fractional coordinates):
  Site 0 (A): (0.0, 0.0, 0.0)
  Site 1 (B): (0.5, 0.5, 0.0)
  Site 2 (C): (0.5, 0.0, 0.5)
  Site 3 (D): (0.0, 0.5, 0.5)

Generated 256 atoms
Expected: 4x4x4 (cells) x 4 (basis sites) = 256 atoms

3. Rocksalt (NaCl) Structure¶

Rocksalt structure consists of two interpenetrating FCC sublattices.

The rocksalt (NaCl) structure is composed of two interpenetrating FCC sublattices: - One sublattice for the cation (e.g., Na⁺) - One sublattice for the anion (e.g., Cl⁻)

Each unit cell contains 4 cations and 4 anions, giving a 1:1 stoichiometry.

In [5]:

Copied!





# Create NaCl structure (lattice constant for NaCl is 5.64 Å)
lat_nacl = Lattice.rocksalt(a=5.64, species_a="Na", species_b="Cl")

print(f"Total basis sites: {len(lat_nacl.basis)}")
species_counts = Counter([site.species for site in lat_nacl.basis])
print(f"Basis composition: {dict(species_counts)}")

# Build NaCl crystal: 3×3×3 unit cells
region_nacl = BlockRegion(0, 3, 0, 3, 0, 3, coord_system="lattice")
builder_nacl = CrystalBuilder(lat_nacl)
structure_nacl = builder_nacl.build_block(region_nacl)

print(f"\nGenerated {len(list(structure_nacl.atoms))} atoms")
symbols = structure_nacl.symbols
atom_counts = Counter(symbols)
print(f"Atomic composition: {dict(atom_counts)}")
print(f"Na:Cl ratio: {atom_counts['Na']}:{atom_counts['Cl']} = 1:1 ✓")
# Create NaCl structure (lattice constant for NaCl is 5.64 Å)
lat_nacl = Lattice.rocksalt(a=5.64, species_a="Na", species_b="Cl")

print(f"Total basis sites: {len(lat_nacl.basis)}")
species_counts = Counter([site.species for site in lat_nacl.basis])
print(f"Basis composition: {dict(species_counts)}")

# Build NaCl crystal: 3×3×3 unit cells
region_nacl = BlockRegion(0, 3, 0, 3, 0, 3, coord_system="lattice")
builder_nacl = CrystalBuilder(lat_nacl)
structure_nacl = builder_nacl.build_block(region_nacl)

print(f"\nGenerated {len(list(structure_nacl.atoms))} atoms")
symbols = structure_nacl.symbols
atom_counts = Counter(symbols)
print(f"Atomic composition: {dict(atom_counts)}")
print(f"Na:Cl ratio: {atom_counts['Na']}:{atom_counts['Cl']} = 1:1 ✓")

Total basis sites: 8
Basis composition: {'Na': 4, 'Cl': 4}

Generated 216 atoms
Atomic composition: {'Na': 108, 'Cl': 108}
Na:Cl ratio: 108:108 = 1:1 ✓

4. Custom Lattices¶

You can define custom lattices with arbitrary lattice vectors and basis sites.

"When to use custom lattices" Custom lattices are useful for: - Non-cubic crystal systems (orthorhombic, tetragonal, etc.) - Complex structures with multiple species - Research on novel materials - Educational purposes to understand crystal structures

In [6]:

Copied!





# Define custom lattice vectors (orthorhombic)
a1 = np.array([3.0, 0.0, 0.0])
a2 = np.array([0.0, 4.0, 0.0])
a3 = np.array([0.0, 0.0, 5.0])

# Define custom basis with multiple species
basis_custom = [
    Site(label="C1", species="C", frac=(0.0, 0.0, 0.0)),
    Site(label="N1", species="N", frac=(0.5, 0.5, 0.5)),
]

lat_custom = Lattice(a1=a1, a2=a2, a3=a3, basis=basis_custom)

# Build custom structure
region_custom = BlockRegion(0, 4, 0, 3, 0, 2, coord_system="lattice")
builder_custom = CrystalBuilder(lat_custom)
structure_custom = builder_custom.build_block(region_custom)

print(f"Generated {len(list(structure_custom.atoms))} atoms")
print(f"Composition: {dict(Counter(structure_custom.symbols))}")
# Define custom lattice vectors (orthorhombic)
a1 = np.array([3.0, 0.0, 0.0])
a2 = np.array([0.0, 4.0, 0.0])
a3 = np.array([0.0, 0.0, 5.0])

# Define custom basis with multiple species
basis_custom = [
    Site(label="C1", species="C", frac=(0.0, 0.0, 0.0)),
    Site(label="N1", species="N", frac=(0.5, 0.5, 0.5)),
]

lat_custom = Lattice(a1=a1, a2=a2, a3=a3, basis=basis_custom)

# Build custom structure
region_custom = BlockRegion(0, 4, 0, 3, 0, 2, coord_system="lattice")
builder_custom = CrystalBuilder(lat_custom)
structure_custom = builder_custom.build_block(region_custom)

print(f"Generated {len(list(structure_custom.atoms))} atoms")
print(f"Composition: {dict(Counter(structure_custom.symbols))}")

Generated 48 atoms
Composition: {'C': 24, 'N': 24}

5. Working with Generated Structures¶

The generated structures are Atomistic objects that integrate with MolPy's ecosystem.

Integration with MolPy Generated Atomistic structures support all standard MolPy operations: - Geometric transformations (move, rotate, scale) - File I/O (export to PDB, XYZ, LAMMPS, etc.) - Analysis tools - Simulation setup

In [7]:

Copied!





# Access coordinates and properties
structure = structure_fcc

print(f"Number of atoms: {len(list(structure.atoms))}")
positions = structure.xyz
print(f"Positions shape: {positions.shape}")
print("\nBox matrix (lattice vectors as rows):")
print(structure["box"].matrix)
print(f"\nBox volume: {structure['box'].volume:.2f} Å³")

# Geometric transformations
structure_copy = builder_fcc.build_block(region_fcc)
xyz = structure_copy.xyz
print(f"\nOriginal center: {xyz.mean(axis=0)}")
structure_copy.move([10.0, 5.0, 0.0])
xyz_moved = structure_copy.xyz
print(f"After translation: {xyz_moved.mean(axis=0)}")
# Access coordinates and properties
structure = structure_fcc

print(f"Number of atoms: {len(list(structure.atoms))}")
positions = structure.xyz
print(f"Positions shape: {positions.shape}")
print("\nBox matrix (lattice vectors as rows):")
print(structure["box"].matrix)
print(f"\nBox volume: {structure['box'].volume:.2f} Å³")

# Geometric transformations
structure_copy = builder_fcc.build_block(region_fcc)
xyz = structure_copy.xyz
print(f"\nOriginal center: {xyz.mean(axis=0)}")
structure_copy.move([10.0, 5.0, 0.0])
xyz_moved = structure_copy.xyz
print(f"After translation: {xyz_moved.mean(axis=0)}")

Number of atoms: 256
Positions shape: (256, 3)

Box matrix (lattice vectors as rows):
[[3.52 0.   0.  ]
 [0.   3.52 0.  ]
 [0.   0.   3.52]]

Box volume: 43.61 Å³

Original center: [6.16 6.16 6.16]
After translation: [16.16 11.16  6.16]

Summary¶

success "Key Features" The crystal builder module provides:

1. **Simple API**: `Lattice` → `Region` → `Builder.build_block()`
2. **Predefined structures**: SC, BCC, FCC, rocksalt
3. **Flexible coordinates**: Lattice or Cartesian units
4. **Integration**: Returns standard `Atomistic` objects
5. **Performance**: Fully vectorized NumPy operations

Next Steps¶

The generated structures can be:

Exported using molpy.io writers (PDB, XYZ, LAMMPS, etc.)
Transformed using molpy.op operations
Combined with other molecules using molpy.pack
Simulated using molpy.engine interfaces

Learn more For more details, see the API documentation.