Figure 1. The architecture of the HEOM modules (Python and C++) in Libra.
The C++ functions can be classified into 3 main groups: a) initialization (Figure 1, left, black); b) computation (Figure 1, left, blue); and c) auxiliary functions (Figure 1, left, green). The first group includes the functions: “gen_hierarchy” that constructs the list of multidimensional index vectors, Eq. 10a, “setup_bath” that computes the Matsubara frequencies, \(\left\{\gamma_{n}\right\}\), and the corresponding coefficients, \(\left\{c_{n}\right\}\), all according to Eqs. 7 and 8, respectively. “compute_matsubara_sum” that computes the sum of type Eq. 9, “initialize_el_phonon_couplings” that computes the matrices\(Q_{m}\) used in Eqs. 5 or 13 as\(Q_{m}=\left.\ |m\right\rangle\left\langle m|\right.\ \). The latter is useful when user does not desire to specify such matrices in the input and wants to use the assumption that each states is coupled to only one bath mode.
The generation of the ADMs’ indices, the indexing convention, and the calculation of the mappings by the function “gen_hierarchy” is worth additional discussion. In our implementation, the hierarchy of ADMs is constructed as shown in Figure 2. The function creates the tree structure representing the desired hierarchy (Figure 2a) via a recursive algorithm. Starting with a d-dimensional all-zero vector, one can increment the index of each of the d components by one to generate d new vectors. These vectors are the “proposed” children nodes for a given iteration (hierarchy level increment). Before each of the new proposed vectors is added to the final list of all the hierarchy indices\(\left\{\mathbf{n}\right\}\), it is compared to all of the previously added children of the same tier (level) to ensure it is unique. The duplication is not possible among the vectors of different tiers, so the comparison to all the previously added vectors is not needed and is not efficient computationally. Considering only the members of the same tier in the uniqueness determination is one of the time-consuming elements of the algorithm. Generating the hierarchy with all possible (non-unique) vectors first and discarding them later is very time and memory consuming algorithm, and it can become unpractical already for relatively low maximal tiers. Therefore, the uniqueness of the added vectors is verified on the fly during the construction. This greatly reduces both CPU and memory requirements of the algorithm.
The newly added unique vectors at a given tier level (children of the nodes of a lower tier) are used as the input for another iteration. Each of them would be considered a starting vector, to which the procedure described above would be applied. As a result, a tree structure is created and traversed as illustrated in Figure 2a. The traversal is organized in the width-first fashion, such that the lower-tier indices appear in the \(\left\{\mathbf{n}\right\}\) list before the d-dimensional indices of higher tiers (Figures 1a, 1b). According to Eq. 10, the dimensionality of each multi-index vectors (a node of the tree) is \(d=M(K+1)\).