Cluster Expansion for Ni-Al with ICET#
This tutorial demonstrates building a cluster expansion for the Ni-Al binary system using ICET and VASP.
Note: Using reduced settings for demonstration. Production CE fitting requires more training structures.
import numpy as np
from ase import Atoms
from ase.build import bulk
# Check for ICET
try:
from icet import ClusterExpansion, ClusterSpace, StructureContainer
from icet.tools import enumerate_structures
from icet.tools.structure_generation import generate_sqs
HAS_ICET = True
except ImportError:
HAS_ICET = False
print("ICET not installed. Run: pip install icet")
from vasp import Vasp
Step 1: Define the Parent Lattice#
The parent lattice is FCC with sites that can be occupied by Ni or Al.
# FCC primitive cell as parent lattice
# Use experimental Ni lattice constant for consistency
a = 3.52 # Ni lattice constant in Angstrom
parent = bulk('Ni', 'fcc', a=a)
print(f"Parent lattice: {parent.get_chemical_formula()}")
print(f"Lattice constant: {a} A")
print(f"Atoms in primitive cell: {len(parent)}")
Parent lattice: Ni
Lattice constant: 3.52 A
Atoms in primitive cell: 1
Step 2: Create Cluster Space#
Define which clusters to include in the expansion.
if HAS_ICET:
# Define cluster space with smaller cutoffs for speed
cutoffs = [4.0, 3.0] # Reduced from [5.0, 4.0]
cs = ClusterSpace(
structure=parent,
cutoffs=cutoffs,
chemical_symbols=['Ni', 'Al']
)
print(cs)
else:
print("Skipping - ICET not installed")
====================================== Cluster Space ======================================
space group : Fm-3m (225)
chemical species : ['Al', 'Ni'] (sublattice A)
cutoffs : 4.0000 3.0000
total number of parameters : 5
number of parameters by order : 0= 1 1= 1 2= 2 3= 1
fractional_position_tolerance : 2e-06
position_tolerance : 1e-05
symprec : 1e-05
-------------------------------------------------------------------------------------------
index | order | radius | multiplicity | orbit_index | multicomponent_vector | sublattices
-------------------------------------------------------------------------------------------
0 | 0 | 0.0000 | 1 | -1 | . | .
1 | 1 | 0.0000 | 1 | 0 | [0] | A
2 | 2 | 1.2445 | 6 | 1 | [0, 0] | A-A
3 | 2 | 1.7600 | 3 | 2 | [0, 0] | A-A
4 | 3 | 1.4370 | 8 | 3 | [0, 0, 0] | A-A-A
===========================================================================================
Step 3: Create Training Structures#
We create a minimal training set with 3 structures for demo purposes.
def create_training_structures():
"""Create minimal training set for Ni-Al.
For cluster expansion, all structures must use the same parent lattice.
We use the FCC lattice with a=3.52 A for all structures.
"""
a = 3.52 # Same lattice constant as parent
structures = []
# 1. Pure Ni (FCC) - primitive cell
ni_pure = bulk('Ni', 'fcc', a=a)
structures.append(('Ni', ni_pure, -5.50))
# 2. Pure Al (FCC) - same lattice as parent (Al on Ni sites)
al_pure = bulk('Al', 'fcc', a=a)
structures.append(('Al', al_pure, -3.74))
# 3. L1_2 ordered structure (Ni3Al) - 4-atom conventional cell
# Use cubic FCC cell and replace 1 of 4 atoms with Al
l12 = bulk('Ni', 'fcc', a=a, cubic=True) # 4-atom cell
l12.symbols[0] = 'Al' # Replace corner atom with Al -> Ni3Al
structures.append(('Ni3Al_L12', l12, -5.10))
return structures
training_set = create_training_structures()
print("Training structures:")
print("-" * 50)
for name, atoms, energy in training_set:
n_ni = sum(1 for s in atoms.symbols if s == 'Ni')
n_al = sum(1 for s in atoms.symbols if s == 'Al')
x_al = n_al / len(atoms)
vol_per_atom = atoms.get_volume() / len(atoms)
print(f"{name:12s} x_Al={x_al:.2f} n_atoms={len(atoms)} vol/atom={vol_per_atom:.4f}")
Training structures:
--------------------------------------------------
Ni x_Al=0.00 n_atoms=1 vol/atom=10.9036
Al x_Al=1.00 n_atoms=1 vol/atom=10.9036
Ni3Al_L12 x_Al=0.25 n_atoms=4 vol/atom=10.9036
Step 4: Calculate DFT Energies with VASP#
Run VASP calculations for each training structure.
def calculate_energy(atoms, label):
"""Calculate energy with VASP."""
# Adjust k-points based on cell size
n_atoms = len(atoms)
if n_atoms <= 2:
kpts = (8, 8, 8)
elif n_atoms <= 4:
kpts = (6, 6, 6)
else:
kpts = (4, 4, 4)
calc = Vasp(
label=f'results/ce_train/{label}',
atoms=atoms,
xc='PBE',
encut=300, # Reduced for demo
kpts=kpts,
ismear=1,
sigma=0.1,
)
return calc.potential_energy / len(atoms) # eV/atom
# Calculate energies for all training structures
print("Calculating DFT energies...")
print("-" * 50)
dft_energies = []
for name, atoms, _ in training_set:
e = calculate_energy(atoms, name)
dft_energies.append((atoms, e))
print(f"{name:12s}: E = {e:.4f} eV/atom")
print()
print("DFT calculations complete.")
Calculating DFT energies...
--------------------------------------------------
Ni : E = -5.4164 eV/atom
Al : E = -2.7806 eV/atom
Ni3Al_L12 : E = -5.4663 eV/atom
DFT calculations complete.
Step 5: Fit the Cluster Expansion#
if HAS_ICET:
from trainstation import CrossValidationEstimator
# Build structure container
sc = StructureContainer(cluster_space=cs)
for (name, atoms, _), (_, energy) in zip(training_set, dft_energies):
sc.add_structure(
structure=atoms,
user_tag=name,
properties={'energy': energy}
)
print(f"Structure container: {len(sc)} structures")
print(sc)
else:
print("Skipping - ICET not installed")
Structure container: 3 structures
==================== Structure Container ====================
Total number of structures: 3
-------------------------------------------------------------
index | user_tag | n_atoms | chemical formula | energy
-------------------------------------------------------------
0 | Ni | 1 | Ni | -5.4164
1 | Al | 1 | Al | -2.7806
2 | Ni3Al_L12 | 4 | AlNi3 | -5.4663
=============================================================
if HAS_ICET:
# Simple fit (not enough data for proper CV)
from sklearn.linear_model import LinearRegression
A, y = sc.get_fit_data(key='energy')
reg = LinearRegression(fit_intercept=False)
reg.fit(A, y)
print(f"Fit R^2: {reg.score(A, y):.4f}")
# Create cluster expansion
ce = ClusterExpansion(
cluster_space=cs,
parameters=reg.coef_
)
print("\nEffective Cluster Interactions (ECIs):")
print(ce)
else:
print("Skipping - ICET not installed")
Fit R^2: 1.0000
Effective Cluster Interactions (ECIs):
================================================ Cluster Expansion ================================================
space group : Fm-3m (225)
chemical species : ['Al', 'Ni'] (sublattice A)
cutoffs : 4.0000 3.0000
total number of parameters : 5
number of parameters by order : 0= 1 1= 1 2= 2 3= 1
fractional_position_tolerance : 2e-06
position_tolerance : 1e-05
symprec : 1e-05
total number of nonzero parameters : 5
number of nonzero parameters by order : 0= 1 1= 1 2= 2 3= 1
-------------------------------------------------------------------------------------------------------------------
index | order | radius | multiplicity | orbit_index | multicomponent_vector | sublattices | parameter | ECI
-------------------------------------------------------------------------------------------------------------------
0 | 0 | 0.0000 | 1 | -1 | . | . | -2.15 | -2.15
1 | 1 | 0.0000 | 1 | 0 | [0] | A | -1.83 | -1.83
2 | 2 | 1.2445 | 6 | 1 | [0, 0] | A-A | 0.196 | 0.0327
3 | 2 | 1.7600 | 3 | 2 | [0, 0] | A-A | -2.15 | -0.716
4 | 3 | 1.4370 | 8 | 3 | [0, 0, 0] | A-A-A | 0.513 | 0.0641
===================================================================================================================
Step 6: Generate SQS Structure#
A Special Quasirandom Structure (SQS) mimics a random alloy in a small supercell.
if HAS_ICET:
from icet.tools.structure_generation import generate_sqs
# Generate small SQS for Ni0.5Al0.5
target_concentrations = {'Ni': 0.5, 'Al': 0.5}
sqs = generate_sqs(
cluster_space=cs,
max_size=8, # Small for demo
target_concentrations=target_concentrations,
n_steps=1000,
)
print(f"SQS structure: {sqs.get_chemical_formula()}")
print(f"Number of atoms: {len(sqs)}")
print(f"Composition: Ni={sum(s=='Ni' for s in sqs.symbols)/len(sqs):.2f}")
# Predict energy with cluster expansion
e_sqs = ce.predict(sqs)
print(f"CE predicted energy: {e_sqs:.4f} eV/atom")
else:
print("Skipping - ICET not installed")
print("\nManual SQS example (without ICET):")
# Create a simple disordered structure manually
sqs_manual = bulk('Ni', 'fcc', a=3.5, cubic=True).repeat((2, 2, 2))
# Replace half with Al randomly
symbols = list(sqs_manual.symbols)
n_replace = len(symbols) // 2
np.random.seed(42)
indices = np.random.choice(len(symbols), n_replace, replace=False)
for i in indices:
symbols[i] = 'Al'
sqs_manual.set_chemical_symbols(symbols)
print(f"Manual SQS: {sqs_manual.get_chemical_formula()}")
print(f"Atoms: {len(sqs_manual)}")
sqs = sqs_manual
SQS structure: Al4Ni4
Number of atoms: 8
Composition: Ni=0.50
CE predicted energy: -2.1473 eV/atom
Summary#
print("=" * 60)
print("Summary")
print("=" * 60)
print()
print("This demo showed:")
print(" 1. Cluster space setup - Define parent lattice and cutoffs")
print(" 2. Training structures - Pure elements + ordered compound")
print(" 3. DFT calculations - VASP energies for training")
print(" 4. CE fitting - Linear regression on cluster correlations")
print(" 5. SQS generation - Random alloy structure")
print()
print("Demo settings used:")
print(" - 3 training structures (minimal)")
print(" - ENCUT = 300 eV")
print(" - Small cutoffs [4.0, 3.0] A")
print()
print("For production use:")
print(" - 20-100 training structures")
print(" - Include relaxed structures")
print(" - LASSO with cross-validation")
print(" - Validate CE predictions with DFT")
print(" - Use Monte Carlo for thermodynamics")
============================================================
Summary
============================================================
This demo showed:
1. Cluster space setup - Define parent lattice and cutoffs
2. Training structures - Pure elements + ordered compound
3. DFT calculations - VASP energies for training
4. CE fitting - Linear regression on cluster correlations
5. SQS generation - Random alloy structure
Demo settings used:
- 3 training structures (minimal)
- ENCUT = 300 eV
- Small cutoffs [4.0, 3.0] A
For production use:
- 20-100 training structures
- Include relaxed structures
- LASSO with cross-validation
- Validate CE predictions with DFT
- Use Monte Carlo for thermodynamics