Package for exactly solving low-dimensional fermionic Hamiltonians
Project description
FermPy
A Python package for exactly solving low-dimensional fermionic Hamiltonians. It was developed with the intention of constructing and diagonalizing matrix representations of low-dimensional, interacting fermionic Hamiltonians, in particular effective (surrogate) models that describe the low-energy (sub-gap) physics of the superconducting Anderson impurity model (see the arXiv preprint).
FermPy is based on NumPy and SciPy. Small systems use dense NumPy arrays as operator representations while large systems employ sparse SciPy arrays.
Functionality of this package includes:
- Construction of matrix representations of operators ($H, S_z$, etc.) from a basis (chosen from symmetry considerations)
- Diagonalization (numpy.linalg.eigh, scipy.sparse.linalg.eigh) of operators
- Calculation of expectation values and thermal averages of operators
Installation
FermPy is found on PyPI and easily installed by running
pip install fermpy
from the command-line.
Usage
As an example of the intended usage of FermPy, we calculate and plot the excitation spectrum and average spin on a quantum dot (QD), coupled to a superconducting (SC) lead which is described by an effective three-level discretization (see the arXiv preprint for more details).
Import the necessary packages.
import numpy as np
import matplotlib.pyplot as plt
from fermpy.operator_sum import OpSum
from fermpy.hamiltonian_constructors import build_ham_op_sum_S_D
from fermpy.operator_representation import OpRep, operator_exp_val
from fermpy.bases import BasisS2
Construct the operator sum for the Hamiltonian $H$ of an SC-QD junction with an effective three-level lead. The effective couplings gamma
and level positions xi
are collected in Table 1 of the Supplemental Material of arXiv preprint. Everything is in units of the SC gap Delta
.
Gamma = np.linspace(0, 10, 201)
h_os = build_ham_op_sum_S_D(epsilon_d=-7.5,
U=15,
Gamma=Gamma,
Delta=1,
gamma=[0.5080, 1.8454, 1.8454],
xi=[0, 2.7546, -2.7546])
Define the basis. Here, we choose the eigenstates of total spin $S^2$ (the Hamiltonian preserves spin rotational symmetry).
basis = BasisS2(1 + 3) # One quantum dot level + three lead levels
Build the matrix representation of the Hamiltonian from the operator sum h_os
and basis
. We only consider the singlet, doublet, and triplet states ($s = 0, 1/2, 1$) - the quantum numbers of $S^2$: $S^2 |s\rangle = s(s+1)|s\rangle$.
h_op = OpRep(h_os, basis, qns_to_consider=(0, 1, 2)) # quantum numbers: 2*s
Diagonalize $H$, giving the eigenenergies and eigenvectors. energies
and eigenvectors
are dictionaries with quantum numbers as keys ($2s = 0, 1, 2$) and NumPy arrays as values.
energies: dict[int, np.ndarray]
eigenvectors: dict[int, np.ndarray]
energies, eigenvectors = h_op.diagonalize_all_qns(return_eigenvectors=True)
Extract the ground state quantum number and energy.
s_gs = np.argmin([E[:, 0] for E in energies.values()], axis=0)
E_gs = np.array([energies[s][idx, 0] for idx, s in enumerate(s_gs)])
Calculate the z-projection of spin on the quantum dot in the ground state. The syntax of OpSum
is similar to ITensor.
sz_op = OpRep(opsum=OpSum([('sz', 0)]), basis=basis)
sz_exp = operator_exp_val(op_rep=sz_op, state_rep=eigenvectors, axis=1)
sz_exp_gs = [sz_exp[s][idx, 0] for idx, s in enumerate(s_gs)]
Plot the excitation spectrum and average dot spin in the ground state.
fig, (ax_top, ax_bot) = plt.subplots(2, 1,
figsize=(5.5, 5.5),
sharex=True,
constrained_layout=True
)
ax_top.set(ylabel=r'$(E - E_\mathrm{GS}) / \Delta$',
xlim=(0, 10),
ylim=(-0.05, 2)
)
ax_bot.set(xlabel=r'$\Gamma / \Delta$',
ylabel=r'$\langle S_{z, \mathrm{dot}} \rangle$',
xlim=(0, 10),
ylim=(-0.05, 0.55)
)
ax_top.hlines(1, 0, 10, linestyle='dashed', color='k', linewidth=0.5)
ax_bot.hlines(0, 0, 10, linestyle='dashed', color='k', linewidth=0.5)
colors = {0: 'tab:blue', 1: 'tab:orange', 2: 'tab:green'}
lines_to_label = []
for s, energies_s in energies.items():
lines_to_label.append(ax_top.plot(Gamma,
energies_s - E_gs[:, None],
color=colors[s])[0]
)
ax_top.legend(handles=lines_to_label,
labels=range(3),
loc='upper right',
title='s'
)
ax_bot.plot(Gamma, sz_exp_gs)
Here's the output figure.
Project details
Release history Release notifications | RSS feed
Download files
Download the file for your platform. If you're not sure which to choose, learn more about installing packages.