User documentation for the HQS Spin Mapper

This is the user documentation for the HQS Spin Mapper package. The intended purpose of HQS Spin Mapper is the investigation of first-principles Hamiltonian descriptions of an input material (this includes both solids and molecules) for orbitals of spin-like character and the subsequent derivation of an effective model Hamiltonian description for this material as a spin-bath model, i.e. a system of spins embedded in a material-specific fermionic bath.

Installation of the HQS Spin Mapper

Installation of the binary package (Recommended)

Several packages need to be installed in a given environment as dependencies of the HQS Spin Mapper. If the HQS Spin Mapper is installed inside a conda environment, this preceding installation of dependencies is facilitated with the command:

conda install -y cxx-compiler boost cmake make hdf5 unzip mkl

The HQS Spin Mapper package is available on the HQS pypi server. To install it, use:

pip install --extra-index-url http://pypiserver.office.quantumsimulations.de/ --trusted-host pypiserver.office.quantumsimulations.de hqs_spin_mapper

Building from source code (for development purposes only!)

To build the package from source code inside a conda environment, the dependencies can be installed via:

conda install -y cxx-compiler boost cmake make hdf5 unzip mkl

If you are not installing the package inside a conda environment, check whether the corresponding packages have been installed on your system, including a sufficiently recent C++ compiler, i.e. g++ >= 10.

Building lattice_solver and lattice_functions (not required)

With the external dependencies installed, one needs to build the lattice_solver and lattice_functions packages. The lattice_solver and lattice_functions packages are dependencies required for the Schrieffer-Wolff transformation / Spin-bath model derivation functionality of Spin Mapper.

The lattice_solver contains the Krylov SVD solver needed for the Schrieffer-Wolff transformation functionality. The repository can be obtained from gitlab and installed via:

git clone https://gitlab.office.quantumsimulations.de/hqs-software/solver/lattice_solver
cd lattice_solver/
make pip_install

The lattice_functions contains the core functionality for handling the fermionic operator algebras used in the Schrieffer-Wolff transformation. The lattice_functions repository is available from gitlab via:

git clone git@gitlab.office.quantumsimulations.de:pschmitteckert/lattice_functions.git

Then enter the lattice_functions directory and verify that the CmakeLists.txt file contains the command set(HDF5_PREFER_PARALLEL TRUE) and add it otherwise.

Then perform:

pip install --extra-index-url http://pypiserver.office.quantumsimulations.de/ --trusted-host pypiserver.office.quantumsimulations.de .

Alternatively, there is a lattice_functions binary on the HQS pypi-server:

pip install --extra-index-url http://pypiserver.office.quantumsimulations.de/ --trusted-host pypiserver.office.quantumsimulations.de lattice_functions

Clone the lattice_2_spin repository and perform:

git clone https://gitlab.office.quantumsimulations.de/hqs-software/lattice/hqs_spin_mapper
cd hqs_spin_mapper/
make pip_install

Theory behind HQS Spin Mapper

In Spin-like orbitals we define the notion of spin-like orbitals and in Local parity we propose the local parity as a metric for the spin-like character of orbitals. In Local parity optimization of the orbital basis, we present our parity optimization procedure as a way to determine the spin-like orbitals of a system. We discuss the basic concept of the Schrieffer-Wolff transformation in Perturbative similarity transformation. We detail our proposed extended Schrieffer-Wolff transformation method in sections Symmetry specification of block-offdiagonal operators and Schrieffer-Wolff transformation as a system of linear equations for unique block-offdiagonal operators . In Full workflow we show how the described methods are combined into our workflow for deriving effective spin-bath model Hamiltonians for materials with relevant low-energy spin physics.

Spin-like orbitals

We consider an orbital ${\varphi }_i$ as being spin-like only if the electron density $n_i$ contained in it is strictly equal to one. This requires that the average electron density in the orbital ${\varphi }_i$ satisfies $⟨n_i⟩=1$. It furthermore requires that the fluctuations around this average electron density satisfy $\delta n_i\to 0$. Negligible fluctuations around the average electron density imply that electron density of the orbital does not affect the low-energy dynamics of the system and vice versa. Just an average electron density $⟨n_i⟩=1$ places no restrictions on the orientation and dynamics of the electron spin in the orbital ${\varphi }_i$. When both requirements are simultaneously met by the states in the low-energy Hilbert space, the dynamics of the electron density in the orbital ${\varphi }_i$, often referred to as the charge degree of freedom, become superfluous to the description of the dynamics of the system. The electron density in the orbital ${\varphi }_i$ consequently couples to the remainder of the system exclusively via its spin degree of freedom. It is then sufficient to represent the degrees of freedom of the electron density contained in the orbital ${\varphi }_i$ as a pure spin degree of freedom and to employ the associated spin operator algebra. The local Hilbert space ${\mathcal{H}}_i$ of a spin degree of freedom is half the size of the local Hilbert space of a fermionic orbital and is furthermore naturally represented on a qubit. We therefore aim to determine each spin-like basis orbital or linear combinations thereof meeting both stated requirements, so that they can be represented as spins.

Local parity

We propose the ground state local parity $P_i$ as a measure for the spin-like character of an orbital ${\varphi }_i$. The operator representation of the local parity reads $$P_i=(-1{)}^{n_{i\uparrow }+n_{i\downarrow }}$$ where $n_{i\uparrow }=c_{i\uparrow }^{†}{c}_{i\uparrow }$ denotes the electron density in ${\varphi }_i$ with electron spin quantum number $s_i^z=+1/2$ and $n_{i\downarrow }=c_{i\downarrow }^{†}{c}_{i\downarrow }$ the electron density with quantum number $s_i^z=-1/2$ respectively. The local Hilbert space ${\mathcal{H}}_i$ of the orbital ${\varphi }_i$ is spanned by the states,

$$\left\{|0⟩,|\uparrow ⟩,|\downarrow ⟩,|\uparrow \downarrow ⟩\right\}$$

and the action of the local parity parity operator on these states reads $$\begin{array}{lcl}{P}_i|0⟩& =& +1|0⟩\\ P_i|\downarrow ⟩& =& -1|\downarrow ⟩\\ P_i|\uparrow ⟩& =& -1|\uparrow ⟩\\ P_i|\uparrow \downarrow ⟩& =& +1|\uparrow \downarrow ⟩.\end{array}$$ For states where the orbital contains a single electron, the local parity operator $P_i$ returns the eigenvalue $p_i=-1$. For the remaining two basis states, $P_i$ returns the eigenvalue $p_i=+1$. Any state $|\psi ⟩$ in the local Hilbert space ${\mathcal{H}}_i$ that contains contributions from non singly occupied basis states hence satisfies $$P_i|\psi ⟩=\left(-1+\alpha \right)|\psi ⟩,$$ with $2\ge \alpha >0$, because the resulting fluctuations in the electron density $\delta n_i\ne 0$ manifest themselves in strictly positive contributions to the local parity. An alternative and more useful operator representation of the local parity reads $$\begin{gathered} P_i=\left(1-2c_{i\uparrow }^{†}{c}_{i\uparrow }\right)\left(1-2c_{i\downarrow }^{†}{c}_{i\downarrow }\right) \\ =1-2\left(c_{i\uparrow }^{†}{c}_{i\uparrow }+c_{i\downarrow }^{†}{c}_{i\downarrow }\right)+4c_{i\uparrow }^{†}{c}_{i\downarrow }^{†}{c}_{i\downarrow }{c}_{i\uparrow }, \end{gathered}$$ and its corresponding expectation value with respect to the many-body state $|n⟩$ can be expressed as $$⟨n|P_i|n⟩=1-2{\rho }_{ii}^{(1)}+4{\rho }_{iiii}^{(2)},$$ where $${\rho }_{qp}^{(1)}=\sum _{\sigma }⟨n|c_{q\sigma }^{†}{c}_{p\sigma }|n⟩,$$ denotes the one-electron reduced density matrix (1-RDM) and $${\rho }_{qprs}^{(2)}=⟨n|c_{q\uparrow }^{†}{c}_{p\downarrow }^{†}{c}_{r\downarrow }{c}_{s\uparrow }|n⟩,$$ denotes the two-electron reduced density matrix (2-RDM) respectively. If we choose $|n⟩=|{\psi }_0⟩$, with $|{\psi }_0⟩$ a good approximation of the ground state of the system, we can identify orbitals ${\varphi }_i$ for which $⟨P_i{⟩}_0=⟨{\psi }_0|P_i|{\psi }_0⟩=-1+\epsilon$, with $\epsilon \to 0$, as spin-like orbitals of the system. In general, the spin-like orbitals of the system do not coincide with the basis orbitals. We therefore require a method to determine the set $\{{\varphi }_i\}$ of orthonormal linear combinations of basis orbitals, for which the local parities most closely approach $\{p_i\}=-1$.

Local parity optimization of the orbital basis

We propose an iterative procedure to determine the particular orbital basis in which the local parities are extremal. For this we attempt a sequence of unitary pairwise rotations of the orbital fermionic operators given by $$\begin{gathered} c_{q\sigma }=\cos \theta c_{i\sigma }+\sin \theta c_{j\sigma }, \\ {c}_{p\sigma }=-\sin \theta c_{i\sigma }+\cos \theta c_{j\sigma }, \end{gathered}$$ with the same rotation being performed for the hermitian conjugates of the operators. From the reduced density matrices ${\rho }^{(1)}$ and ${\rho }^{(2)}$, we can compute the local parity of an orbital ${\varphi }_q$, which results from the linear combination of orbitals ${\varphi }_i$ and ${\varphi }_j$, as $⟨P_q{⟩}_0(\theta )$. The local parity $⟨P_q{⟩}_0(\theta )$ is an analytic, $2\pi$-periodic function of the rotation angle $\theta$. We find extremal points ${\theta }_n$ of the function $⟨P_q{⟩}_0(\theta )$ in the domain $\theta \in [0,2\pi )$ from $${\frac{d⟨P_q{⟩}_0}{d\theta }|}_{{\theta }_n}=0,$$ and select solutions ${\theta }_n$ that satisfy $${\frac{d^2⟨P_q{⟩}_0}{d{\theta }^2}|}_{{\theta }_n}\ne 0.$$ If a solution ${\theta }_m$ exists which satisfies $d^2⟨P_q{⟩}_0/d{\theta }^2{|}_{{\theta }_m}>0$ and $⟨P_q{⟩}_0({\theta }_m)\le ⟨P_q{⟩}_0({\theta }_l)\forall {\theta }_l\in {\theta }_n$, we accept the rotation attempt $i\to q$ and $j\to p$ with rotation angle ${\theta }_m$ and we reject the rotation attempt otherwise. By repeating the procedure for each pair of basis orbitals ${\varphi }_i$ and ${\varphi }_j$, we arrive at at an orthonormal basis in which the local parities have taken up extremal values. We then identify the orbitals ${\varphi }_q$ of the resulting basis for which $⟨P_q{⟩}_0=-1+\epsilon$, with $\epsilon$ an arbitrary small, positive value, as the spin-like orbitals of the system. Basis orbitals ${\varphi }_p$ with a local parity $⟨P_p{⟩}_0\simeq +1$ can be regarded as beneficial for the purpose of separating the system's spin degrees of freedom and their respective environment since they experience exclusively even number particle transfer, such that the spin degree of freedom of electrons occupying the orbitals ${\varphi }_p$ becomes insignificant.

Perturbative similarity transformation

In principle, any Hamiltonian can be completely diagonalized by means of a particular unitary transformation $U$. In practice, finding the particular unitary transformation $U$ often requires a complete diagonalization of the Hamiltonian to begin with. Here we recap how a perturbative similarity transformation, namely the Schrieffer-Wolff transformation, can be used to determine an approximate transformation operator $U$, or the generator $S$ thereof, which does not diagonalize the Hamiltonian fully, but yields a block-diagonal Hamiltonian instead. These blocks consist of the orthonormal states in the Hilbert space $\mathcal{H}$ that share the a given choice of characteristics, e.g. the local particle quantum number $n_i$. We will denote the set of terms in the Hamiltonian that are already block-diagonal in the initial basis as $H_0$. The remaining terms connect different blocks, i.e are block-offdiagonal, and are denoted $V$. The complete Hamiltonian thus reads $$H=H_0+V.$$ A unitary similarity transformation of the Hamiltonian is given by $$H=U^{†}HU={\text{e}}^{-S}\left(H_0+V\right){\text{e}}^S={\text{e}}^{-S}{H}_0{\text{e}}^S+{\text{e}}^{-S}V{\text{e}}^S,$$ where $S$ is an anti-hermitian operator. One refers to it as the generator of the Schrieffer-Wolff transformation. The key problem of the Schrieffer-Wolff transformation is finding the generator $S$ such that the transformed Hamiltonian $H$ becomes entirely block-diagonal. In order to arrive at an equation for $S$ one makes use of the Campbell-Baker-Hausdorff formula to expand $$H={\mathrm{e}}^{S}H{\mathrm{e}}^{-S}=\sum _{m=0}^{\infty }\frac{1}{m!}[S,H{]}_m=H+[S,H]+\frac{1}{2}\left[S,[S,H]\right]+\dots ,$$ where for generators $S$ satisfying $\Vert S\Vert \ll \Vert H_0\Vert$, with $\Vert \cdot \Vert$ a suitable norm, one can approximate the expression as $$H=H_0+V+\left[S,H_0\right]+\left[S,V\right]+\mathcal{O}(S^2).$$ Considering that commutators of pairs of block-diagonal operators or pairs of block-offdiagonal operators respectively generally become block-diagonal, while the commutators of block-diagonal operators with block-offdiagonal operators become block-offdiagonal, one chooses the equation $$[S,H_0]=-V\iff [H_0,S]=V,$$ by which one can determine the generator $S$ which removes the block-offdiagonal terms $V$ of the Hamiltonian. If a solution $S$ exists, one can use $$\left[S,\left[S,H_0\right]\right]=\left[S,-V\right]=-\left[S,V\right]\notin \mathcal{O}(S^2),$$ to simplify the expression for the transformed Hamiltonian $$H={\mathrm{e}}^{-S}\left(H_0+V\right){\mathrm{e}}^S\simeq H_0+\frac{1}{2}\left[S,V\right]+\mathcal{O}(S^2),$$ where the terms originating from $[S,V]/2$ contain the perturbative corrections arising from the consecutive application of two block-offdiagonal operators. A subsequent projection to the subspaces $${\mathcal{P}}_n=\sum _{p\in {\mathcal{P}}_n}|p⟩⟨p|,$$ yields the block-diagonal Hamiltonian $$H_{\mathrm{block}-\mathrm{diagonal}}=\sum _n{\mathcal{P}}_n\tilde{H}{\mathcal{P}}_n,$$ where $n$ denotes the distinct blocks of the Hilbert space.

Symmetry specification of block-offdiagonal operators

The block-offdiagonal part of the Hamiltonian $V$ consists of a sum of block-offdiagonal terms. These in turn comprise products of individually block-offdiagonal operators $x$. In the following we detail a method to decompose generic block-offdiagonal operators $x$ into distinct components. Each of the components exclusively connects two distinct blocks of the Hilbert space $\mathcal{H}$, often associated with distinct quantum numbers of a symmetry of the system. A given block-offdiagonal operator $x:\mathcal{H}\to \mathcal{H}$ satisfies $$\left[H_0,x\right]=\epsilon z\ne 0,$$ where $\epsilon$ denotes an arbitrary scalar and $z:\mathcal{H}\to \mathcal{H}$ an arbitrary operator. Let $A$ be a diagonal operator in the initial basis. It can be identical to the symmetry operator differentiating the blocks of the Hilbert space, but it is not required to be. One can use the spectrum of $A$ to expand the operator $x$ as $$x=\sum _q{x}_q=\sum _q{\beta }_q\prod _{i\ne q}\left(A-a_i\right)x,$$ where the different $x_q$ couple the target subspace associated with the eigenvalue $a_q$ of $A$ to other subspaces of the Hilbert space. If the operator $A$ satisfies $$[x,A]=0,$$ then both subspaces, initial and final, coupled via $x_q$ are specified by the eigenvalue $a_q$. This is possible for the fermionic creation and annihilation operators $c_{i\sigma }^{†}$ and $c_{i\sigma }$. The coefficients ${\beta }_q$ are solutions to the equation $${\beta }_q\prod _{i\ne q}\left(a_q-a_i\right)=1.$$ The symmetry-specified block-offdiagonal operators $x_q$ satisfy $$\left[x_q,x_{q^{\prime }}^{†}\right]\propto {\delta }_{q{q}^{\prime }},$$ where the operator $x_q^{†}$ denotes the hermitian conjugate of the operator $x_q$.

Schrieffer-Wolff transformation as a system of linear equations for unique block-offdiagonal operators

Following the procedure outlined in the section Perturbative similarity transformation we separate the Hamiltonian of the system into its block-diagonal and block-offdiagonal contributions as $$H=H_0+V.$$ The block-offdiagonal contribution $V$ comprises each block-offdiagonal term $v$ $$\begin{gathered} V=\sum _{\left\{v\right\}}{\alpha }_v\left[\left(\prod _j{o}^j\prod _i{x}^i\right)+{\left(\prod _j{o}^j\prod _i{x}^i\right)}^{†}\right] \\ =\sum _{\left\{v\right\}}{\alpha }_v\left[\left(\prod _{i,j}{o}^j\sum _{q(i)}{x}_{q(i)}^i\right)+{\left(\prod _{i,j}{o}^j\sum _{q(i)}{x}_{q(i)}^i\right)}^{†}\right], \end{gathered}$$ where ${\prod }_i{x}^i$denotes sequences of individually block-offdiagonal operators, ${\prod }_j{o}^j$ denotes sequences of individually block-diagonal operators, and $x_{q(i)}^i$ denotes the symmetry-specified components of the operator $x^i$. We introduce a vector spaces ${\mathcal{V}}_0^h$ and ${\mathcal{V}}_0^a$, for which each unique pair of hermitian, or anti-hermitian respectively, symmetry-specified operator sequences in $V$ corresponds to a unique orthonormal basis vector $${\hat{e}}_v=\{\begin{array}{ll}\left({\prod }_j{o}^j{\prod }_i{x}_{q(i)}^i\right)+{\left({\prod }_j{o}^j{\prod }_i{x}_{q(i)}^i\right)}^{†}& \in {\mathcal{V}}_0^h\\ \left({\prod }_j{o}^j{\prod }_i{x}_{q(i)}^i\right)-{\left({\prod }_j{o}^j{\prod }_i{x}_{q(i)}^i\right)}^{†}& \in {\mathcal{V}}_0^a\end{array},$$ where ${\mathcal{V}}^h$ denotes the hermitian vector space and ${\mathcal{V}}^a$ the anti-hermitian vector space. We define the linear map $$L:\begin{array}{l}{\mathcal{V}}_0^h\to {\mathcal{V}}_1^a\\ {\mathcal{V}}_0^a\to {\mathcal{V}}_1^h\end{array},$$ where the action of the linear map on a vector ${\hat{e}}_v$ is given by $$L{\hat{e}}_v=\left[H_0,{\hat{e}}_v\right],$$ where $$\begin{array}{lcl}(L{\hat{e}}_v)\in {\mathcal{V}}_1^a& \mathrm{if}& {\hat{e}}_v\in V_0^h\\ (L{\hat{e}}_v)\in {\mathcal{V}}_1^h& \mathrm{if}& {\hat{e}}_v\in V_0^a\end{array},$$ with $$\dim \left({\mathcal{V}}_1\right)\ge \dim \left({\mathcal{V}}_0\right),$$ since ${\mathcal{V}}_1$ includes additional unique operator sequences generated by $[H_0,{\hat{e}}_v]$. In the vector spaces $\mathcal{V}$, we can translate the equation for the generator $S$ of the Schrieffer-Wolff transformation to $$L\vec{S}=\vec{V}\iff \left[H_0,S\right]=V,$$ with $\vec{S}\in {\mathcal{V}}_0^a$ and $\vec{V}\in {\mathcal{V}}_1^h$. Determining the generator $S$ of the Schrieffer-Wolff transformation becomes equivalent to solving the set of linear equations. In general one finds $$\mathrm{rank}(L)\le \dim \left({\mathcal{V}}_0\right)\le \dim \left({\mathcal{V}}_1\right),$$ and there is consequently no unique solution $\vec{S}$ to the set of equations. It is intuitive that the terms of $V$ approximately bring about transitions between distinct eigenstates of $H_0$ belonging to different blocks of the Hilbert space. This is reflected by $$\left[H_0,{\hat{e}}_v^h\right]=\left(\Delta E_{0,v}\prod _l{o}^l\right){\hat{e}}_v^a,$$ where $\Delta E_{0,v}$ denotes the eigenvalue difference between the two eigenstates of $H_0$ and ${\prod }_l{o}^l$ denotes an arbitrary sequence of individually block-diagonal operators. From the previous equations we can approximate $$\Vert \vec{S}\Vert \approx \sqrt{\sum _{\{v\}}{\left(\frac{{\alpha }_v}{\Delta E_{0,v}}\right)}^2},$$ which highlights the necessity for a significant energy gap $\Delta E_{0,v}$ between the separate subspaces of the Hilbert space coupled by $V$ in order for the series expansion of $H$ to $\mathcal{O}(S^2)$ to be considered a good approximation. We consider terms $v$ for which $$\frac{{\alpha }_v^2}{|\Delta E_{0,v}|}\ge 1,$$ to be the resonant terms of $V$ which should be retained in the transformed Hamiltonian $H$. To identify these resonant terms of $V$, we employ a singular value decomposition (SVD) of the linear map $L$ and arrive at $$L=U\Sigma W^{†}=U\left({\Sigma }_{>}+{\Sigma }_{<}\right)W^{†}=L_{\mathrm{gapped}}+L_{\mathrm{resonant}},$$ where ${\Sigma }_{>}$ comprises the singular values $${\sigma }_i>\frac{1}{N_v}\sum _{\{v\}}\sqrt{{\alpha }_v^2},$$ where $N_v$ denotes the number of terms $v$. We define the Moore-Penrose pseudoinverse of $L$ that acts exclusively on the gapped, i.e. non-resonant, terms in $V$ as $$L_{\mathrm{gapped}}^{+}=\left(W{\Sigma }_{>}^{+}{U}^{†}\right),$$ where ${\Sigma }_{>}^{+}$ denotes the pseudoinverse of ${\Sigma }_{>}$. In the absence of a unique solution to $\left[H_0,S\right]=V$, the closest approximate solution is given by $$\vec{S}=\left(L_{\mathrm{gapped}}^{+}\vec{V}\right)\in {\mathcal{V}}_0^a,$$ and the resulting transformed Hamiltonian reads $$\tilde{H}=H_0+V_{\mathrm{resonant}}+\left[S,V_{\mathrm{resonant}}\right]+\frac{1}{2}\left[S,V_{\mathrm{gapped}}\right],$$ where $V_{\mathrm{gapped}}=[H_0,S]$ and $V_{\mathrm{resonant}}=(V-V_{\mathrm{gapped}})$.

Full workflow

In the following, we outline the steps of the workflow that we use to derive an effective spin-bath model Hamiltonian from a first principles description of a material.

Computation of the required system information

We start with an ab-initio electronic structure calculation of the material to determine its ground state. The electronic structure method needs to be a post-Hartree-Fock or related method - this excludes density functional theory - to capture the effect of correlations in the two-particle reduced density matrix ${\rho }^{(2)}$. From the electronic structure calculation we obtain the basis orbitals in which the Hamiltonian of the system is formulated, and the one-electron and two-electron integrals which specify the Hamiltonian description. For the ground state of the calculation we compute the one-particle and two-particle reduced density matrices ${\rho }^{(1)}$ and ${\rho }^{(2)}$.

Determination of spin-like basis orbitals

We utilize the reduced density matrices to assign a local parity $⟨P_i{⟩}_0$ to the basis orbitals i${\varphi }_i$. We then perform pairwise rotations of the basis orbitals to determine the basis in which the local parities of the basis orbitals are extremized. If there exist optimized basis orbitals ${\varphi }_q$ with $⟨P_q{⟩}_0+1<\epsilon$, where we typically choose $\epsilon \le 10^{-1}$, we proceed with the subsequent steps of the workflow. If no spin-like orbitals are found, we terminate the workflow.

Schrieffer-Wolff transformation of the Hamiltonian

We use our Schrieffer-Wolff transformation approach to integrate out the valid terms of the Hamiltonian that modify the electron density in the spin-like orbitals, which leads to renormalized couplings of the electron spins in the spin-like orbitals to the environment. The valid terms are the ones that couple subspaces of the Hilbert space between which there exists a significant energy gap.

Construction of the effective spin-bath Hamiltonian

The transformed Hamiltonian $H$ is projected to the particular subspace of the Hilbert space where the electron density of spin-like orbitals ${\varphi }_q$ is fixed to $n_q\equiv 1$. Utilizing the identity $n_q=n_{q\uparrow }+n_{q\downarrow }\equiv 1$ fermionic operators acting on the spin-like operators are substituted with the corresponding spin operators. The resulting Hamiltonian $\mathcal{P}H\mathcal{P}$ is the effective spin-bath representation of the material.

Representation on a device (optional)

The effective spin-bath model Hamiltonian $\mathcal{P}\tilde{H}\mathcal{P}$ is re-expressed in terms of the spin operators that are realized on the specific device.

Applications of the HQS Spin Mapper package

Outline

• Accepted input formats
• Functionality: Spin-like orbital basis optimization
• Preparation of input data for the Schrieffer-Wolff transformation
• Functionality: Effective spin-bath model derivation
• Extraction of tensor descriptions
• Serialization with struqture
• Best practices for HQS Spin Mapper

Applications of the HQS Spin Mapper software package

The HQS Spin Mapper software package contains a selection of modules which provide two major functionalities. The first functionality is the identification of orbital bases containing spin-like orbitals, i.e. orbitals with an occupancy of strictly a single electron. The second one derives effective spin-bath model Hamiltonians for the determined spin-like orbitals. The remaining modules handle input and output tasks.

Accepted input formats (hqs_spin_mapper.transformables)

HQS Spin Mapper currently accepts several different input formats for the Hamiltonian descriptions of a system (i.e. crystals and molecules). Some of the input formats store additional information about the system, which is, for example, needed for the determination of spin-like orbital bases or to reduce the computational of the spin-bath model derivation. There is also the possibility for the user to provide their own input format through the use of the supplied protocols.

Currently supported by default are:

• LatticeBuilder configuration object instances
• LatticeBuilder configuration data files as *.yml
• Struqture FermionHamiltonianSystem object instances
• Data container (Transformable_QuantumChemistry) including the one- and two-particle integrals specifying the Hamiltonian stored as numpy binary *.npy files
• Data container (Transformable_Matrices) including the single particle picture matrices of the Hamiltonian as numpy array instances

For each default input format, we provide a specific class derived from the Transformable class. The Transformable object serves as the necessary input argument of the essential methods of the HQS Spin Mapper package and are manipulated by these methods. The Transformable classes satisfy the Supports_SW_Transformation protocol. Any object satisfying the Supports_SW_Transformation protocol can be used with the HQS Spin Mapper functions.

LatticeBuilder object instances and data files (Transformable_LatticeBuilder)

Given an instance of a LatticeBuilder object, it is easy to instantiate a Transformable_LatticeBuilder object from, which can then be used as the input for the subsequent functions. To instantiate the Transformable_LatticeBuilder object, use:

from hqs_spin_mapper.transformables import Transformable_LatticeBuilder

system = Transformable_LatticeBuilder(input_builder)

If the LatticeBuilder description of the system is stored in a valid *.yml file, the Transformable_LatticeBuilder object is instantiated in the same way.

from hqs_spin_mapper.transformables import Transformable_LatticeBuilder

system = Transformable_LatticeBuilder(input_lattice_builder_yaml)

Struqture FermionHamiltonianSystem instances (Transformable_Struqture)

A Transformable_Struqture object can be directly instantiated from a struqture_py.fermions.FermionHamiltonianSystem instance.

from hqs_spin_mapper.transformables import Transformable_Struqture

system = Transformable_Struqture(input_struqture)

Hamiltonians as numpy array instances (Transformable_Matrices)

The description of the Hamiltonian of a given system can be provided as a collection of matrices or tensors respectively, in the single particle picture. The matrices and tensors first need to be stored in a data container class called Hamiltonian_Matrices.

from hqs_spin_mapper.spin_mapper_protocols import Hamiltonian_Matrices

hamiltonian_matrices = Hamiltonian_Matrices(H0_uu=H0_uu,
                                            H0_dd=H0_dd,
                                            H0_ud=H0_ud,
                                            HU=HU,
                                            Jz=None,
                                            Jp=None,
                                            Jc=None)

The spin indices of the quadratic $c_{i\sigma }^{†}{c}_{j{\sigma }^{\prime }}$ fermionic terms H0_uu, H0_dd, H0_ud indicate the spin component of the electrons. If provided as a matrix, HU denotes density-density interaction terms. If HU is provided as a tensor, more complex four index interaction terms can be specified. The matrix Jz specifies $S_i^z{S}_j^z$ interactions, Jp specifies $S_i^{+}{S}_j^{-}$ interactions, and Jc specifies ${\vec{S}}_i\times {\vec{S}}_j$ interactions.

The quadratic Hamiltonian matrices H0_uu, H0_dd, H0_ud have to follow the convention $${\left({\left[H_0\right]}_{\sigma {\sigma }^{\prime }}\right)}_{ij}={\left({\left[H_0\right]}_{\sigma {\sigma }^{\prime }}\right)}_{ij}{c}_{i\sigma }^{†}{c}_{j{\sigma }^{\prime }}$$ The matrices H0_uu, H0_dd, H0_ud have to be provided for three different spin direction combinations $\sigma {\sigma }^{\prime }$, i.e. $(\uparrow ,\uparrow )$, $(\downarrow ,\downarrow )$ and $(\uparrow ,\downarrow )$.

If provided as a matrix, a HU has to conform to $${\left(HU\right)}_{ij}=(HU{)}_{ij}\sum _{\sigma {\sigma }^{\prime }}{n}_{i\sigma }^{†}{n}_{j{\sigma }^{\prime }}^{†}$$ and if provided as a tensor as $${\left(HU\right)}_{ijkl}=(HU{)}_{ijkl}{c}_{i\sigma }^{†}{c}_{j{\sigma }^{\prime }}^{†}{c}_{k{\sigma }^{\prime }}{c}_{l{\sigma }^{\prime }}$$

Using the Hamiltonian_Matrices instance, one can instantiate a Transformable_Matrices object that satisfies the Supports_SW_Transformation protocol via

from hqs_spin_mapper.transformables import Transformable_Matrices

system = Transformable_Matrices(hamiltonian_matrices)

Electron integrals from electronic structure codes (Transformable_QuantumChemistry)

If the output of a preceding electronic structure calculations contains the one-electron and two-electron integrals specifying the Hamiltonian stored as numpy binary *.npy files, one can construct a Transformable from the corresponding matrices via instantiation of a Transformable_Matrices object as shown in the previous section.

from hqs_spin_mapper.spin_mapper_protocols import Hamiltonian_Matrices
from hqs_spin_mapper.transformables import Transformable_Matrices

H0_uu = np.load('h0_uu.npy')
H0_dd = np.load('h0_dd.npy')
H0_ud = np.load('h0_ud.npy')
HU = np.load('hu.npy')

hamiltonian_matrices = Hamiltonian_Matrices(H0_uu=H0_uu,
                                            H0_dd=H0_dd,
                                            H0_ud=H0_ud,
                                            HU=HU,
                                            Jz=None,
                                            Jp=None,
                                            Jc=None)

system = Transformable_Matrices(hamiltonian_matrices)

If the results from the preceding electronic structure calculations also include the one-electron density matrix (1-RDM) and the two-electron density matrix (2-RDM), one can also instantiate a Transformable_QuantumChemistry object. The Transformable_QuantumChemistry object also satisfies the Supports_Spin_Optimization protocol. This means that the object can be used as the input argument for the spin basis optimization functionality of Spin Mapper.

from hqs_spin_mapper.transformables import Transformable_QuantumChemistry
from hqs_spin_mapper.spin_mapper_protocols import Hamiltonian_Matrices

hamiltonian_matrices = Hamiltonian_Matrices(H0_uu=H0_uu,
                                            H0_dd=H0_dd,
                                            H0_ud=H0_ud,
                                            HU=HU,
                                            Jz=None,
                                            Jp=None,
                                            Jc=None)

rdm1_uu = np.load('rdm1_uu.npy')
rdm1_dd = np.load('rdm1_dd.npy')
rdm1_ud = np.load('rdm1_ud.npy')
rdm2 = np.load('rdm2.npy')

system = Transformable_QuantumChemistry(hamiltonian_matrices,
                                        rdm1=(rdm1_uu, rdm1_dd, rdm1_ud),
                                        rdm2=rdm2,
                                        tolerance=1e-1)

The input reduced density matrices have to adhere to specific conventions. The spin orbital optimization will otherwise yield non-physical results. We expect the following conventions for the reduced density matrices:

• single 1-RDM: $${\rho }_{ij}^{(1)}=\sum _{\sigma }⟨c_{i\sigma }^{†}{c}_{j\sigma }⟩$$

• spin-resolved 1-RDMs: $${{\rho }_{ij}^{(1)}|}_{\sigma {\sigma }^{\prime }}=⟨c_{i\sigma }^{†}{c}_{j{\sigma }^{\prime }}⟩$$

• single 2-RDM: $${\rho }_{ijkl}^{(2)}=⟨c_{i\uparrow }^{†}{c}_{j\downarrow }^{†}{c}_{k\downarrow }{c}_{l\uparrow }⟩$$

Input exclusively for spin-like orbital optimization (SpinFinder)

If the spin-bath model derivation functionality is not required, the user can instantiate a SpinFinder object, which is the valid minimum input for the functions pertaining to the determination of the spin-like orbital basis.

from hqs_spin_mapper.transformables import SpinFinder

rdm1_uu = np.load('rdm1_uu.npy')
rdm1_dd = np.load('rdm1_dd.npy')
rdm1_ud = np.load('rdm1_ud.npy')
rdm2 = np.load('rdm2.npy')

system = SpinFinder(rdm1=(rdm1_uu, rdm1_dd, rdm1_ud),
                    rdm2=rdm2)

Keyworded constructor arguments for Transformable objects

Fields described by the Supports_SW_Transformation protocol:

• system_size: Optional[int] = None: Number of orbitals in the system. This is automatically derived from the model description. (available for Transformable_Matrices, Transformable_QuantumChemistry)
• rotation_matrix: Optional[numpy.ndarray] = None: Basis transformation matrix $U_{iq}$ between different orbital bases. The index $i$ denotes the orbital in the old basis and $q$ the orbital in the new basis. Default is the identity. (available for Transformable_LatticeBuilder, Transformable_Matrices, Transformable_QuantumChemistry)
• prefactor_cutoff: float = 1e-6: Terms of the Hamiltonian with coupling constant of absolute value smaller than the chosen value are discarded in the derivation of the spin-bath model (Schrieffer-Wolff transformation). (available for Transformable_LatticeBuilder, Transformable_Struqture, Transformable_Matrices, Transformable_QuantumChemistry)
• site_type: str = "spinful_fermions": Type of particles occupying the orbitals. (available for Transformable_Matrices, Transformable_QuantumChemistry)
• spin_indices: Optional[Set[int]] = None: User choice of orbitals to be treated as spin-like. The spin-like orbital optimization will add indices to this set. (available for Transformable_LatticeBuilder, Transformable_Struqture, Transformable_Matrices, Transformable_QuantumChemistry)

Fields described by the Supports_Spin_Optimization protocol:

• tolerance: float = 1e-1: Chosen discrepancy $\Delta$ used in $P_i\le -1+\Delta$ for a basis orbital $i$ to be considered spin-like, where $P_i$ denotes the local parity. (available for Transformable_QuantumChemistry, SpinFinder)
• optimization_steps: int = 4: Number of optimization loops to be performed. Parities are typically converged after 4 loops. (available for Transformable_QuantumChemistry, SpinFinder)

Other mutable fields of Transformable objects

Fields described by the Supports_SW_Transformation protocol:

• generator: Union[ExpressionSpinful, ExpressionSpinful_complex]: Result for the generator of the Schrieffer-Wolff transformation.
• transformed_hamiltonian: Union[ExpressionSpinful, ExpressionSpinful_complex]: Result for the transformed fermionic Hamiltonian.

Fields described by the Supports_Spin_Optimization protocol:

• immutable_indices: Set[int]: Indices of orbitals that are exempt from parity optimization.

Custom Transformable objects (hqs_spin_mapper.spin_mapper_protocols)

The user can create their own custom input format by implementing a class that satisfies the protocols Supports_Spin_Optimization and/or Supports_SW_Transformation.

Functionality: Spin-like orbital basis optimization (hqs_spin_mapper.orbital_optimization)

HQS Spin Mapper provides the functionality to determine the optimal single-particle orbital basis, in which a subset of the basis orbitals exhibits the maximum possible spin-like character. We quantify the spin-like character of basis orbitals on the basis of their local parity $P_i=(-1{)}^{n_{i\uparrow }+n_{i\downarrow }}$. For more details about the use of the local parity as a measure for spin-like character, please consult the Theory section. An orbital $q$ in the optimized basis is stored as spin-like if its local parity satisfies $P_q\le -1+\Delta$, where $\Delta$ corresponds to the value stored in the tolerance field of the Transformable.

The spin_orbital_optimization function in the orbital_optimization module determines a basis transformation which would optimize the local parities and stores the corresponding (proposed) basis transformation in the rotation_matrix field of the Transformable. The respective indices of the spin-like orbitals in this optimized basis are stored in the spin_indices field of the Transformable. The spin_orbital_optimization function requires Transformable objects that satisfy the Supports_Spin_Optimization protocol.

from hqs_spin_mapper.orbital_optimization import spin_orbital_optimization

spin_orbital_optimization(system)

Keyworded arguments

• space_reduction: bool = True: If flag is set, the optimization is restricted to the indices not contained in the immutable_indices field of the Transformable.
• fast_mode: bool = False: If flag is set, the rotations of the 2-RDM are performed approximately to reduce computation time.
• target: str = "minimum": Type of extremum to be searched for. Options are: minimum, maximum, extremum.
• _delta_sufficiently_spin_like: float = 1e-3: Exclude orbitals with initial parities $P_i<-1+\delta$ and $P_i>+1-\delta$ from the optimization procedure.

The function spin_orbital_optimization rotates the rdm1 and rdm2 to the optimized basis (in order to save memory), but not the matrix descriptions of the Hamiltonians.

The results of the optimized system can subsequently be accessed via:

spin_indices = system.spin_indices
rotation_matrix = system.rotation_matrix

rdm1 = system.rdm1
rdm2 = system.rdm2

Preparation of input data for the Schrieffer-Wolff transformation (hqs_spin_mapper.preconditioning)

For the derivation of the spin-bath model, the spin_indices: Set[int] field of the Transformable instance first requires a finite amount of indices. If the spin orbital optimization feature was used, spin_indices of the Transformable have been set to feature the obtained spin orbitals automatically. Otherwise spin_indices needs to be set manually by the user as in:

spin_like_orbital_indices: Set[int] = {4, 9}

system.spin_indices = spin_like_orbital_indices

An intuitive choice for spin-like orbitals are the orbitals with particularly strong repulsive on-site interactions ${\left[HU\right]}_{iiii}\gg 0$.

The necessary preconditioning of the input data is done automatically by calling the preconditioning function on the Transformable once.

from hqs_spin_mapper.preconditioning import preconditioning

preconditioning(system)

The following tasks are performed by the preconditioning function:

• Generalization of the input interaction terms HU, Jz, Jp, Jc to a four index tensor description.

• Rotation of the matrix and tensor descriptions to the optimized basis by rotation with the rotation_matrix $U_{iq}$. $${\left[H_0\right]}_{qp}=\sum _{ij}{\left[H_0\right]}_{ij}{U}_{iq}{U}_{jp}$$ $${\left[HU\right]}_{qprs}=\sum _{ijkl}{\left[HU\right]}_{ijkl}{U}_{iq}{U}_{jp}{U}_{kr}{U}_{ls}$$

• Mean-field decoupling (if Supports_Spin_Optimization) / Removal of bath exclusive interaction terms (if apply_pre_diagonalization = True). $${\tilde{H}}_0^{\mathrm{bath}}=\sum _{qp\in \mathrm{bath}}\sum _{\sigma {\sigma }^{\prime }}{t}_{qp}^{\sigma {\sigma }^{\prime }}{c}_{q\sigma }^{†}{c}_{p{\sigma }^{\prime }}-\sum _{qprs\in \mathrm{bath}}{V}_{qprs}\left(⟨c_{q\uparrow }^{†}{c}_{s\uparrow }⟩c_{p\downarrow }^{†}{c}_{r\downarrow }-⟨c_{p\downarrow }^{†}{c}_{r\downarrow }⟩c_{q\uparrow }^{†}{c}_{s\uparrow }-⟨c_{q\uparrow }^{†}{c}_{r\downarrow }⟩c_{p\downarrow }^{†}{c}_{s\uparrow }-⟨c_{p\downarrow }^{†}{c}_{s\uparrow }⟩c_{q\uparrow }^{†}{c}_{r\downarrow }\right)$$

• Diagonalization of the quadratic bath Hamiltonian ${\tilde{H}}_0$ and rotation of the remaining Hamiltonians to the diagonal basis (if apply_pre_diagonalization = True).

• Placing the spin indices in the middle of the system at $x_0+q$ (if apply_cross_geometry = True), with $x_0$ defined as

x0 = (system.system_size // 2) - (len(system.spin_indices) // 2)

Keyworded arguments

• apply_pre_diagonalization: bool = True: Perform simplification of bath exclusive interactions and diagonalization of the resulting quadratic bath Hamiltonian.
• apply_cross_geometry: bool = False: Move the spin indices to the center of the system (speeds up subsequent DMRG calculations).

Functionality: Effective spin-bath model derivation (hqs_spin_mapper.schrieffer_wolff)

With the input data prepared, the effective spin-bath model derivation becomes a simple call of schrieffer_wolff on the Transformable object.

from hqs_spin_mapper.schrieffer_wolff import schrieffer_wolff

schrieffer_wolff(system)

The results of the derivation, namely the transformed Hamiltonian $H$ and the generator of the transformation $S$ in $$H=H_0+V+\left[S,H_0\right]+\left[S,V\right]+\mathcal{O}(S^2)$$ are stored in the Transformable object instances as:

transformed_hamiltonian = system.transformed_hamiltonian
generator = system.generator

Keyworded arguments

• _max_length: int = 5: Size limit for the terms included in the generator.
• _vector_space_cap: int = 1500000: Size limit of the vector space that the generator is determined from.
• _max_krylov_space_dimension: int = 200: Size restriction of the Krylov space used in the SVD.
• rdm2: Optional[numpy.ndarray] = None: 2RDM for use in decoupling (currently not activated).

Extract tensor descriptions for the spin-bath model (hqs_spin_mapper.extract_couplings)

Using the functions contained in the extract_couplings module, one can retrieve the matrix/tensor description of the spin-bath model Hamiltonian.

from hqs_spin_mapper.extract_couplings import extract_quadratic, extract_quartic

(H0_uu, H0_dd, H0_ud) = extract_quadratic(system)
HU = extract_quartic(system)

The matrices and tensor only contain the quadratic/quartic terms of the Hamiltonian. The terms containing more fermionic operators $(n\ge 5)$ are not extracted!

Serialize the spin-bath model Hamiltonian with struqture (hqs_spin_mapper.struqture_interface)

The transformed Hamiltonian can be exported to an instance of struqture_py.mixed_systems.MixedHamiltonianSystem for use in other HQS software or for storage in a file.

from hqs_spin_mapper.struqture_interface import get_struqture_description

struqture_description = get_struqture_description(system)

Best practices for HQS Spin Mapper

Avoiding local extrema in the orbital optimization

Depending on the initial basis, it is possible that a single call to spin_orbital_optimization fails to determine the true parity optimized basis. The pairwise optimization algorithm can get stuck in local parity extrema of the parameter manifold. This behavior can be discouraged by using the following procedure:

• Call the spin_orbital_optimization function on the input data with target="minimum".
• Add the orbital indices that became spin-like to the immutable_indices of the Transformable instance.
• Call the spin_orbital_optimization function a second time on the updated input data with either target="extremum" or target="maximum".
• Call the spin_orbital_optimization function on the input data a final time with target="minimum".

Choosing a sensible value for prefactor_cutoff (Important to avoid memory overflow!)

To limit the system memory required for the Schrieffer-Wolff transformation, we recommend adjusting the value of prefactor_cutoff. The terms of the input Hamiltonian in the optimized basis with an absolute coupling constant smaller than prefactor_cutoff will be dropped from the Hamiltonian prior to the Schrieffer-Wolff transformation step. This is done because a large number of small coupling constant terms slows down the SVD step of the transformation significantly, while not significantly changing the result. In the evaluation of the commutators that yield the transformed Hamiltonian, a larger number of terms can cause memory usage to exceed the system memory and resulting in a fatal error. We recommend inspecting the size of the coupling constants of the input Hamiltonian and choosing a value of prefactor_cutoff larger than the coupling constants that can be considered insignificant. When in doubt, choose a larger value for prefactor_cutoff first and then rerun the calculation with a slightly smaller value of prefactor_cutoff to check whether the result is sensitive to the change in prefactor_cutoff. If the result remains sensitive, we recommend gradually decreasing prefactor_cutoff and observing the memory usage.

We recommend to run the scripts or notebooks using the HQS Spin Mapper in an environment with a set ulimit. The ulimit can be set in the terminal before running the python interpreter via:

ulimit -v {maximum_desired_memory_in_kilobytes}

This way the process gets safely terminated without potentially crashing the system.

Estimation of the required memory

We provide a method for a rough estimate of the required memory to perform the Schrieffer-Wolff transformation at the currently set prefactor_cutoff. The function can be called on the Transformable object using:

from hqs_spin_mapper.preconditioning import (
    estimate_required_memory
)

estimate_required_memory(data)

The function estimate_required_memory prints an estimate of the necessary memory in GB. If the estimated value exceeds the realistically available system memory, we recommended to raise the value of the field prefactor_cutoff and to check the estimated memory requirement again.

Post-processing of the transformed Hamiltonian

The result of the model Hamiltonian derivation is stored as an ExpressionSpinful object. These objects offer a set of methods that can be used for further manipulation of the transformed Hamiltonian. Here we list some useful methods of the ExpressionSpinful class.

• normalize(eps = 0.0): Orders the operators in the terms of the Hamiltonian in a standardized fashion and removes terms with an absolute prefactor smaller than eps.
• normal_order: Rewrites operators as creators and annihilators and normal orders them (creators to the left of annihilators).
• sort_descending: Sorts the terms in the expression by the absolute value of the coupling constants in descending order.
• sort_ascending: Sorts the terms in the expression by the absolute value of the coupling constants in ascending order.
• find(fermion): Returns the prefactor of the term containing the input operator.
• project_on_local_hilbert_subspace(x, condition): Projects the Hamiltonian on a subspace of the Hilbert space where the condition is satisfied for site x. condition is the integer corresponding to a bit string, where the 1-bit indicates the x-local $|0⟩$ state, the 2-bit indicates $|\uparrow ⟩$, the 4-bit indicates $|\downarrow ⟩$, and the 8-bit indicates $|\uparrow \downarrow ⟩$. The argument condition=6=4+2 thus corresponds to $\text{span}(|\downarrow ⟩,|\uparrow ⟩)$, namely the singly occupied site x.
• convert_to_spin_model(positions): Replace fermionic operators acting on the positions with spin operators assuming the identity $n_{i\uparrow }+n_{i\downarrow }\equiv 1$.

A simple script to convert the transformed fermionic Hamiltonian to a true spin-bath Hamiltonian description is given by:

spin_bath_system = system.transformed_hamiltonian

for n in cast(List[int], system.spin_indices):
    spin_bath_system.project_on_local_hilbert_subspace(n, 6)

spin_bath_system = spin_bath_system.convert_to_spin_model(system.spin_indices)

Examples

Some example python scripts of a spin-bath model derivation for the well-known Single Impurity Anderson Model (SIAM) and the Hubbard Model are provided in the examples folder.

For the real valued SIAM:

cd ./examples
python3 siam_example.py

For the complex valued SIAM:

cd ./examples
python3 complex_siam_example.py

For the real valued Hubbard model:

cd ./examples
python3 hubbard_example.py

For a Kondo lattice chain model

cd ./examples
python3 kondo_lattice_example.py