PARPHOM: PARallel PHOnon calculator for Moiré systems

This paper introduces PARPHOM, a parallel computational code designed to overcome the challenges of calculating phonon properties in large-scale twisted two-dimensional material systems by enabling force constant generation, band structure analysis, and finite-temperature dynamics investigations.

The Gist

Imagine you have two sheets of graphene (a material as thin as a single atom). If you stack them perfectly on top of each other, they behave like a standard sheet of metal. But, if you twist one sheet slightly relative to the other, something magical happens: the atoms create a giant, repeating pattern called a Moiré pattern. It looks like the interference pattern you see when you hold two window screens slightly out of alignment.

Scientists call this field "Twistronics." By changing the twist angle, they can turn these materials into superconductors, insulators, or semiconductors.

However, there's a big problem. When you twist these sheets, the resulting pattern is huge. Instead of a tiny unit cell with just a few atoms, a twisted system might contain thousands of atoms.

The Problem: The "Too Big to Fit" Dilemma

To understand how these materials work, scientists need to study their phonons. Think of phonons as the "vibrations" or "jiggles" of the atoms. Just like a guitar string vibrates to make a sound, atoms vibrate to conduct heat and electricity.

Calculating these vibrations for a normal crystal is easy. But for a twisted Moiré system with thousands of atoms, the math becomes so massive that it crashes standard computers. It's like trying to solve a Sudoku puzzle where the grid is the size of a football field, and you only have a calculator the size of a matchbox. The memory required is simply too high.

The Solution: PARPHOM (The Parallel Powerhouse)

The authors of this paper built a new software tool called PARPHOM. Think of PARPHOM as a giant construction crew instead of a single worker.

• The Old Way: A single worker tries to build a skyscraper brick by brick. It takes forever, and they get tired (run out of memory).
• The PARPHOM Way: PARPHOM hires thousands of workers (processors) who all work on different parts of the building simultaneously. It breaks the massive math problem into tiny, independent chunks that can be solved in parallel.

How It Works (The Analogy)

Here is the step-by-step process of what PARPHOM does, translated into everyday terms:

1. The Setup (The Dance Floor):
   First, the software takes the twisted structure and lets the atoms "relax." Imagine a crowded dance floor where everyone is bumping into each other. The atoms shuffle around until they find the most comfortable, stable spot. This is done using a simulation tool called LAMMPS.

2. The Push Test (Force Constants):
   To understand how the atoms vibrate, the software performs a "push test." It gently nudges one atom slightly to the left, then to the right, and measures how the other atoms react.

   • Analogy: Imagine a row of dominoes. If you push the first one, the others fall in a specific pattern. PARPHOM calculates exactly how hard you need to push one atom to make the whole system wiggle in a specific way. Because there are thousands of atoms, doing this one by one would take years. PARPHOM does it all at once using its army of processors.

3. The Symphony (Phonon Dispersion):
   Once it knows how the atoms react to pushes, it calculates the "music" of the material. It figures out the specific frequencies (notes) the atoms can vibrate at. This creates a "band structure," which is like a musical score showing which notes are allowed and which are forbidden.

4. The Heat Check (Temperature Dynamics):
   Real materials aren't at absolute zero; they are hot. Heat makes atoms jiggle more chaotically. PARPHOM can simulate this "heat dance." It tracks the atoms over time to see how their vibrations change when the temperature rises, revealing how the material conducts heat.

5. The Spin (Chirality):
   Some vibrations in these twisted materials have a "handedness" or chirality. They spin either clockwise or counter-clockwise, like a corkscrew. PARPHOM can detect this spin, which is crucial for future technologies like quantum computing.

Why This Matters

Before PARPHOM, scientists were stuck. They could study small, simple crystals, but the exciting, complex twisted materials were too big to analyze accurately.

PARPHOM is the key that unlocks the door. It allows researchers to:

• Study materials with 10,000+ atoms (which was previously impossible).
• Predict how these materials will behave at different temperatures.
• Design new materials with specific properties by simply changing the twist angle.

In short, PARPHOM is a high-speed, parallel computing engine that lets us listen to the "music" of the most complex, twisted atomic structures in the universe, helping us build the electronics of the future.

Technical Summary

1. Problem Statement

The emergence of twistronics (the study of 2D materials with a relative twist angle) has revealed unique electronic and phononic properties, such as significant phonon renormalization and enhanced electron-phonon interactions. However, calculating these properties poses a severe computational challenge:

• System Size: Moiré patterns formed at small twist angles require unit cells containing hundreds to thousands of atoms (often exceeding 10,000 atoms).
• Memory and Scalability: Standard Density Functional Theory (DFT) is limited to unit cells with a few hundred atoms. Even classical force-field-based approaches using existing codes (e.g., PHONOPY, PhonoLAMMPS) struggle because they typically generate force constants on a single core. The memory requirements for diagonalizing dynamical matrices for such large systems often exceed available resources, making these calculations infeasible with current tools.

2. Methodology

The authors present PARPHOM, a massively parallel software package designed to overcome these bottlenecks. The workflow is divided into three main stages:

A. Force Constant Generation

• Input: Uses structures relaxed via classical force fields in the LAMMPS package.
• Algorithm: Employs the finite displacement method with central finite differences.
  • Atoms are displaced by a default distance of $10^{-4}$ Å.
  • Forces are calculated for all atoms in response to these displacements.
• Parallelization: The generation of force constants is "embarrassingly parallel." The code distributes the $N$ atoms across $N_p$ processors. Each processor independently computes the force constant matrix elements resulting from displacing specific atoms.
• Symmetry Handling:
  • Instead of using symmetry operations to reduce the number of displacements (which adds overhead for massive parallelism), PARPHOM performs all $6N$ displacements.
  • It enforces translational symmetry (Acoustic Sum Rule) and index symmetry via a computationally inexpensive iterative method.
  • It explicitly notes that discrete rotational symmetries are often broken in relaxed low-angle moiré systems and are not explicitly enforced, relying instead on the force field's inherent symmetry.
• Data Storage: Force constants are stored in HDF5 format using parallel h5py to handle large datasets efficiently.

B. Phonon Dispersion Calculation

• Dynamical Matrix Construction: The dynamical matrix is constructed at specific $q$-points in reciprocal space using the generated force constants.
• Parallel Diagonalization: The code utilizes ScaLAPACK routines (specifically PZHEEVX) to diagonalize the Hermitian dynamical matrices. The matrix is distributed in a block-cyclic fashion across processors, allowing the calculation of eigenvalues (phonon frequencies) and eigenvectors for systems with thousands of atoms.
• Group Velocity: Velocities are calculated using the derivative of the dynamical matrix with respect to the wave vector ($q$), employing either analytic derivatives or central finite differences to avoid issues at band crossings.

C. Post-Processing Utilities (Python Package pyphutil)

The package includes a Python suite for advanced analysis:

• Density of States (DOS): Computed via Linear Triangulation (breaking the Brillouin zone into Delaunay triangles) or via Gaussian/Lorentzian smearing. This method is parallelized.
• Finite Temperature Dynamics: Calculates the Velocity Autocorrelation Function (VACF) from Molecular Dynamics (MD) trajectories to determine temperature-dependent DOS and mode lifetimes.
• Mode-Projected VACF: Analyzes the evolution of individual phonon modes with temperature, providing linewidths and frequency shifts (anharmonic effects).
• Chirality Analysis: Computes the chirality of phonon modes (circular polarization) in the $x, y, z$ directions, specifically targeting chiral modes at $K$ and $K'$ valleys in Transition Metal Dichalcogenides (TMDs).

3. Key Contributions

1. Massive Parallelism: Unlike existing tools, PARPHOM parallelizes both the force constant generation and the dynamical matrix diagonalization, enabling calculations on systems with >10,000 atoms.
2. Scalability: Benchmarks show the force constant generation scales efficiently with system size (near $O(n)$ to $O(n^{1.5})$ depending on the regime) using 64 MPI processes.
3. Comprehensive Analysis Suite: It is not just a calculator but a full analysis tool that handles:
   • Zero-temperature band structures and DOS.
   • Finite-temperature anharmonic effects (softening, linewidths).
   • Topological properties like phonon chirality.
4. Open Source: The code is released under GPLv3, written in Fortran (core) and Python (interface/post-processing), and integrates with standard libraries (LAMMPS, ScaLAPACK, HDF5, Spglib).

4. Results and Validation

The authors validated PARPHOM using twisted bilayer graphene (tBLG) and twisted bilayer WSe$_2$:

• Twisted Bilayer Graphene (21.8°): Successfully generated phonon band structures and compared DOS at 0K and 300K.
• Twisted Bilayer WSe$_2$ (2°): Calculated DOS for a system with 4,902 atoms in the moiré unit cell on an $18 \times 18$ $q$-grid.
• Anharmonicity: For 1.08° tBLG, the code analyzed a supercell with 401,904 atoms to compute the mode-projected power spectrum. It successfully detected a softening of the layer breathing mode by 2.6 cm$^{-1}$ at 300K compared to 0K, demonstrating its ability to capture anharmonic effects.
• Chirality: The code successfully identified chiral phonon modes near the $K$ and $K'$ valleys in 38.2° twisted bilayer MoS$_2$, confirming theoretical predictions.

5. Significance

PARPHOM fills a critical gap in the computational study of moiré materials. By enabling phonon calculations for systems with tens of thousands of atoms, it allows researchers to:

• Accurately model the renormalization of phonon frequencies due to atomic relaxation in twisted systems.
• Investigate finite-temperature dynamics and anharmonic effects in large supercells, which are crucial for understanding thermal transport and electron-phonon coupling in twistronics.
• Explore topological phonon properties (chirality) in complex moiré lattices.

The tool democratizes access to high-fidelity phonon calculations for large-scale 2D heterostructures, moving beyond the limitations of DFT and single-core classical codes.