Introduction to ferrodispcalc

What is ferrodispcalc?

ferrodispcalc is a Python library for computing order parameters in ferroelectric materials, including ionic displacements, local polarization, and octahedral tilting. It is designed for perovskite-type structures but can be extended to other systems as long as the coordination environment is well-defined.

In addition, ferrodispcalc provides high-performance I/O utilities for LAMMPS, including fast trajectory reading (read_lammps_dump) and log file parsing (read_lammps_log). See the LAMMPS output tutorial and the ferrodispcalc.io API reference for details.

The Core Idea: Coordination and Neighbor Lists

Note

This tutorial assumes basic familiarity with ASE’s Atoms object. See Introduction to ASE if you are new to ASE.

All calculations in ferrodispcalc boil down to one thing: finding the right neighbors.

Take a perovskite ABO₃ as an example. The B-site cation sits at the center of an oxygen octahedron, and the A-site cation occupies the corner of the unit cell. If you can correctly identify which O atoms coordinate a given B atom (and vice versa), you can compute:

Ionic displacement: the offset of the center atom from the geometric center of its neighbors.
Local polarization: displacement weighted by Born effective charges, converted to C/m².
Octahedral tilt: the rotation angles of the BO₆ octahedron about the Cartesian axes.

ferrodispcalc uses build_neighbor_list to construct these coordination lists. You specify the center element, the neighbor element, a cutoff distance, and the expected coordination number. Everything else follows from there.

from ferrodispcalc.neighborlist import build_neighbor_list

# B-O coordination (octahedron): 6 neighbors
nl_bo = build_neighbor_list(atoms, ["Ti"], ["O"], cutoff=4, neighbor_num=6)

# B-A coordination: 8 neighbors
nl_ba = build_neighbor_list(atoms, ["Ti"], ["Pb"], cutoff=4, neighbor_num=8)

# A-O coordination: 12 neighbors
nl_ao = build_neighbor_list(atoms, ["Pb"], ["O"], cutoff=4, neighbor_num=12)

Example: PbTiO₃ Displacement and Polarization

This example computes ionic displacements and local polarization for a PbTiO₃ (PTO) structure.

from ase.io import read
from ferrodispcalc.neighborlist import build_neighbor_list
from ferrodispcalc.compute import calculate_displacement, calculate_polarization
from ferrodispcalc.config import BEC

# 1. Read structure
atoms = read("PTO.vasp")

# 2. Build neighbor lists
nl_bo = build_neighbor_list(atoms, ["Ti"], ["O"], cutoff=4, neighbor_num=6)
nl_ba = build_neighbor_list(atoms, ["Ti"], ["Pb"], cutoff=4, neighbor_num=8)
nl_ao = build_neighbor_list(atoms, ["Pb"], ["O"], cutoff=4, neighbor_num=12)

# 3. Calculate displacements
disp_bo = calculate_displacement(atoms, nl_bo)   # Ti off-center in O octahedron
disp_ao = calculate_displacement(atoms, nl_ao)   # Pb off-center in O cage

# 4. Calculate polarization
bec = BEC.get("PTO")   # {'Pb': 3.44, 'Ti': 5.18, 'O': -2.87}
P = calculate_polarization(atoms, nl_ba, nl_bo, bec)

calculate_displacement returns an array of shape (n_centers, 3), where each row is the displacement vector of a center atom relative to the geometric center of its neighbors. For trajectories (list[Atoms]), the output shape is (n_frames, n_centers, 3).

calculate_polarization converts these displacements into local polarization (in C/m²) using Born effective charges. It requires two neighbor lists: B-A and B-X.

Beyond Perovskites

ferrodispcalc is not limited to simple ABO₃ perovskites. As long as you can define the coordination environment, you can use the same workflow:

Identify center and neighbor elements.
Choose an appropriate cutoff and coordination number.
Build the neighbor list and compute.

This makes it applicable to:

Solid solutions (e.g. PbTiO₃-SrTiO₃): use mixed element lists like ["Ti"] for B-site and ["Pb", "Sr"] for A-site, and provide averaged Born effective charges. See Example: PTO-STO Superlattice Polarization.
Molecular ferroelectrics: use bond vectors as dipole proxies. See Example: Molecular Ferroelectric (TMCM-CdCl3) for a TMCM-CdCl₃ example.
Other coordination geometries: any system where the order parameter can be defined as a displacement from the coordination center. See Example: Wurtzite Ferroelectric (ScAlN) for a wurtzite example with a custom neighbor list.

# Solid solution example: PTO-STO superlattice
nl_ba = build_neighbor_list(atoms, ["Ti"], ["Pb", "Sr"], cutoff=4, neighbor_num=8)
nl_bo = build_neighbor_list(atoms, ["Ti"], ["O"], cutoff=4, neighbor_num=6)

bec = {"Pb": 3.44, "Sr": 2.56, "Ti": 5.18, "O": -2.87}
P = calculate_polarization(atoms, nl_ba, nl_bo, bec)