Hierarchical Equations of Motion

The hierarchical equations of motion (HEOM) approach was originally developed by Tanimura and Kubo (1989) in the context of physical chemistry to "exactly" solve a quantum system (labeled as $\text{s}$) in contact with a bosonic environment, encapsulated in the following total Hamiltonian:

$${\hat{H}}_{\text{total}}={\hat{H}}_{\text{s}}+\sum _k{\omega }_k{\hat{b}}_k^{†}{\hat{b}}_k+{\hat{V}}_{\text{s}}\sum _k{g}_k({\hat{b}}_k+{\hat{b}}_k^{†}),$$

where ${\hat{b}}_k$ (${\hat{b}}_k^{†}$) is the bosonic annihilation (creation) operator associated to the $k$th mode (with frequency ${\omega }_k$), ${\hat{V}}_{\text{s}}$ refer to the coupling operator acting on the system's degree of freedom, and $g_k$ are the coupling strengths.

As in other solutions to this problem, the properties of the bath are encapsulated by its temperature and its spectral density,

$$J(\omega )=2\pi \sum _k{g}_k^2\delta (\omega -{\omega }_k).$$

In the HEOM approach, for bosonic baths, one typically chooses a Drude-Lorentz spectral density:

$$J_{\text{DL}}(\omega )=\frac{4\mathrm{\Delta }W\omega }{{\omega }^2+W^2},$$

or an under-damped Brownian motion spectral density,

$$J_{\text{U}}(\omega )=\frac{2{\mathrm{\Delta }}^{2}W\omega }{({\omega }^2-{\omega }_0^2{)}^2+{\omega }^2{W}^2}.$$

Here, $\mathrm{\Delta }$ represents the coupling strength between the system and the bosonic bath with band-width $W$ and resonance frequency ${\omega }_0$.

We introduce an efficient Julia framework for HEOM approach called HierarchicalEOM.jl. This package is built upon QuantumToolbox.jl and provides a user-friendly and efficient tool to simulate complex open quantum systems based on HEOM approach. For a detailed explanation of this package, we recommend to read its documentation and also the article Huang et al. (2023).

Given the spectral density, the HEOM approach requires a decomposition of the bath correlation functions in terms of exponentials. In the documentation of HierarchicalEOM.jl, we not only describe how this is done for both bosonic and fermionic environments with code examples, but also describe how to solve the time evolution (dynamics), steady-states, and spectra based on HEOM approach.