elphy

elphy, short for electron-phonon anharmonicity, is a C program to calculate free energies and forces for tight-binding models on supercells, complemented with linear electron-phonon coupling and a harmonic potential. It is a faster rewrite of the Python module elphmod.md and thus primarily a driver for the path-integral molecular-dynamics code i-PI. It is written in ANSI C (C89/C99).

Installation

The following command compiles the program using GCC without optimization:

make

Different compilers and optimization flags can be selected:

make CC=icx CFLAGS=-O2

LAPACK and BLAS are required.

Usage

The program accepts one, two, or three arguments:

elphy <data file>
elphy <data file> <socket>
elphy <data file> <init file> <radius>

With one argument, it alternately reads atomic positions from standard input and writes energies and forces to standard output, both using the XYZ format.

With two arguments, it exchanges these quantities with i-PI through its socket interface. <socket> is a host name optionally followed by a colon and a port number, e.g., localhost:31415. If the port number is omitted or zero, a UNIX socket is used for communication, otherwise an internet socket. The address must match the information in the i-PI input file input.xml.

With three arguments, it creates an <init file> with initial atomic positions in the XYZ format for i-PI, adding random displacements smaller than <radius>, and prints the corresponding energies and forces.

The <data file> is defined below:

<temperature kT>
<number of electrons per unit cell>
<number of orbitals per unit cell>
<maximum number of electrons per orbital>
A₀₀ A₀₁ A₀₂
A₁₀ A₁₁ A₁₂
A₂₀ A₂₁ A₂₂
a₀₀ a₀₁ a₀₂
a₁₀ a₁₁ a₁₂
a₂₀ a₂₁ a₂₂
<number of atoms per unit cell>
X₀ r₀₀ r₀₁ r₀₂ F₀₀ F₀₁ F₀₂
X₁ r₁₀ r₁₁ r₁₂ F₁₀ F₁₁ F₁₂
X₂ r₂₀ r₂₁ r₂₂ F₂₀ F₂₁ F₂₂
⋮
<number of lattice vectors>
R₀₀ R₀₁ R₀₂
R₁₀ R₁₁ R₁₂
R₂₀ R₂₁ R₂₂
⋮
<number of hopping parameters>
i₀ α₀ β₀ <0 α₀|H|Rᵢ₀ β₀>
i₁ α₁ β₁ <0 α₁|H|Rᵢ₁ β₁>
i₂ α₂ β₂ <0 α₂|H|Rᵢ₂ β₂>
⋮
<number of interatomic force constants>
j₀ x₀ y₀ ∂²E/[∂u(0, x₀) ∂u(Rⱼ₀, y₀)]
j₁ x₁ y₁ ∂²E/[∂u(0, x₁) ∂u(Rⱼ₁, y₁)]
j₂ x₂ y₂ ∂²E/[∂u(0, x₂) ∂u(Rⱼ₂, y₂)]
⋮
<number of electron-phonon matrix elements>
k₀ z₀ l₀ γ₀ δ₀ <0 γ₀|∂H/∂u(Rₖ₀, z₀)|Rₗ₀ δ₀>
k₁ z₁ l₁ γ₁ δ₁ <0 γ₁|∂H/∂u(Rₖ₁, z₁)|Rₗ₁ δ₁>
k₂ z₂ l₂ γ₂ δ₂ <0 γ₂|∂H/∂u(Rₖ₂, z₂)|Rₗ₂ δ₂>
⋮

The indices i, j, k, l run over lattice vectors, α, β, γ, δ over orbitals, and x, y, z over the three Cartesian displacement directions for all atoms. The primitive, position, and force vectors a, r, F are given in Cartesian, the supercell and lattice vectors A, R in integer crystal coordinates.

No unit conversions are performed, so any consistent energy and length units can be used. However, i-PI expects energies and forces in Hartree atomic units.

Note that the interatomic force constants are assumed to be partially screened: They shall exclude the harmonic term of the electronic potential-energy surface. Any forces that the model may generate at zero displacements can be compensated by adding a force correction specified next to the atomic positions.

Tests and examples

• make input.dat model=TaS2 (or graphene) creates an example input file.
• make test verifies that the computed free energy and forces are correct.
• make ipi lets elphy and i-PI perform a structural relaxation together. The files generated by i-PI can be cleaned up with make clean_ipi.
• make ipi_elphmod does the same job using elphmod instead of elphy.
• make show displays an animation of the relaxation process.

The Python packages elphmod, ipi, and matplotlib are required.

Acknowledgment

The method is described in the following paper:

• Arne Schobert, Jan Berges, Erik G. C. P. van Loon, Michael A. Sentef, Sergey Brener, Mariana Rossi, and Tim O. Wehling, SciPost Phys. 16, 046 (2024)