Efficient Calculation of the Maximal Rényi Divergence for a Matrix Product State via Generalized Eigenvalue Density Matrix Renormalization Group

TL;DR

The paper tackles the exponential cost of computing quantum mutual information $I(A;B)$ by focusing on the maximal Rényi divergence $D_{\infty}$, which can be formulated as a generalized eigenvalue problem. It develops a generalized DMRG framework with a Lanczos-based solver to efficiently compute $D_{\infty}$ for 1D states represented as matrix product states, including edge- and general-case strategies. In XXZ-chain benchmarks, $D_{\infty}$ shows distinct behavior from $I(A;B)$, signaling phase transitions while remaining correlated with entropy-based measures. The work provides a scalable, GP-aware tool for robust correlation diagnostics in 1D many-body systems and extends to other problems described by generalized eigenvalue problems.

Abstract

The study of quantum and classical correlations between subsystems is fundamental to understanding many-body physics. In quantum information theory, the quantum mutual information, $I(A;B)$, is a measure of correlation between the subsystems $A,B$ in a quantum state, and is defined by the means of the von Neumann entropy: $I\left(A;B\right)=S\left(ρ_{A}\right)+S\left(ρ_{B}\right)-S\left(ρ_{AB}\right)$. However, such a computation requires an exponential amount of resources. This is a defining feature of quantum systems, the infamous ``curse of dimensionality'' . Other measures, which are based on Rényi divergences instead of von Neumann entropy, were suggested as alternatives in a recent paper showing them to possess important theoretical features, and making them leading candidates as mutual information measures. In this work, we concentrate on the maximal Rényi divergence. This measure can be shown to be the solution of a generalized eigenvalue problem. To calculate it efficiently for a 1D state represented as a matrix product state, we develop a generalized eigenvalue version of the density matrix renormalization group algorithm. We benchmark our method for the paradigmatic XXZ chain, and show that the maximal Rényi divergence may exhibit different trends than the von Neumann mutual information.

Figures (6)

• Figure 1: Tensor network graphical notations. (a) Construction of an MPS from a high-dimensional quantum state vector. The state can be written as: $\ket{\Phi}=\sum_{\sigma_{1}\dots\sigma_{N}}\Phi_{\sigma_{1}\dots\sigma_{N}}\ket{\sigma_{1}\dots\sigma_{N}}=\sum_{j}A^{s_{1}\dots s_{i}}v^{i}A^{s_{i+1}\dots s_{N}}\ket{s_{1}\dots s_{N}}$, with $A^{s_{i}}$ site matrices of site $i$. Their contraction yields back the vector components $\Phi_i$, showing that the MPS does represent a state. (b) Left-orthogonality in an MPS. Each tensor to the left of the orthogonality center satisfies the condition of Eq. \ref{['eq:left-canonical']}, ensuring that contractions of these tensors with their conjugates yield the identity on the virtual bonds. This structure enables efficient and stable computation of observables and local optimizations. (c) Structure of an MPS with each set of Schmidt values explicitly depicted, representing Eq. \ref{['eq:MPS-definition-with-S']}. In this formulation, site tensors $\Gamma^{s_i}$ (circles) are not canonical on their own, but can be made left (right) canonical by contracting with the Schmidt values $\Lambda^i$ (diamonds) to their left (right). (d) Structure of an MPO. Each site tensor in the MPO carries two physical indices (input and output) and one or two virtual indices, allowing it to represent local operator actions and their correlations across sites. The MPO structure mirrors and complements the MPS ansatz, enabling efficient representation and application of many-body operators, such as Hamiltonians, within the tensor network framework.
• Figure 2: Two-site effective operator in the DMRG algorithm. The operator acts on two adjacent MPS tensors (sites $i$ and $i+1$) during a local update step, while the surrounding environment tensors appearing in (a) are contracted to form the left and right effective environments in (b). This construction enables the variational optimization of the two-site block within the full many-body context.
• Figure 3: Construction of the subsystem density matrix can be done by using the state MPS structure. Site matrices (circles) are saved in left and right canonical form. Schmidt values (diamonds) are saved as well. Using the Schmidt decomposition: $\ket{\Phi}=\sum_{i}S_{i}\ket{i_{A}}\bra{i_{E,B}}$ the reduced density matrix $\rho_{A}$ is simply $\rho_{A}=\sum_{i}S_{i}^2\ket{i_{A}}\bra{i_{A}}$. Therefore, the eigenvector decomposition of $\rho_{A}$ is given by the left canonical matrices of A, with eigenvalues being the Schmidt values at the subsystem edge (green diamonds). In terms of Eq. \ref{['eq:MPS-definition-with-S']}, the eigenvalues of $\rho_A$ are the diagonal elements of the matrix $\Lambda^i$ at its edge, with the corresponding eigenvectors $\sum_{s_i} \Gamma^{s_1}\Lambda^1\cdots\Gamma^i$. Consequently, manipulations of $\rho_{A}$ such as calculation of $\rho_{A}^{-1/2}$ can be done by raising the Schmidt values from the diagonal of $\Lambda^i$ to power $-1/2$; when contracting the bra and ket MPSs, we then square these values. (green diamonds), giving $\rho_{A}^{-1/2},\rho_{B}^{-1/2}$ (top). Chaining these gives $\sigma^{-1/2}=\rho_{A}^{-1/2}\otimes\rho_{B}^{-1/2}$, which can be used to construct $\sigma^{-1/2}\rho\sigma^{-1/2}$ (bottom).
• Figure 4: Two types of subsystems architectures were chosen to be studied. Top: An 'AEB' architecture where the subsystems are equally sized and are at the edges of the total system. This structure is discussed in Sec. \ref{['sec:System-Edge-Calculation']}, allowing an efficient calculation of the maximal divergence using a special case method. Bottom: An 'EAEBE' architecture where the subsystems are equal in size and lie at equal distances from each other and from the edges.
• Figure 5: Testing our method on the XXZ system for fixed $J=1,h=0$, for a system of $N=10$ sites, and subsystems with $N_S=3$ sites. Top: Calculation of $\braket{S_{i}^{z}S_{i+1}^{z}}$ in the ground state, displaying the XXZ phases: At low $\Delta$ the correlation is large and negative, meaning an anti-ferromagnet. At high $\Delta$ the correlation is high and positive, meaning a ferromagnet. In between there is a region of a paramagnetic phase. Center: Calculation in the 'AEB' subsystems architecture, shown in Fig. \ref{['fig:arch']}, displaying the von Neumann mutual information (blue), the maximal Rényi divergence (orange, green) calculated using the special case method discussed in Sec. \ref{['sec:System-Edge-Calculation']} (orange) and the general case method discussed in Sec. \ref{['sec:general-case-calculation']} (green). The mutual information and the Rényi divergence are not continuous at the paramagnetic-ferromagnetic phase transition and share asymptotic behavior, as explained in the main text, but otherwise differ. Bottom: Calculation in the 'EAEBE' subsystems architecture, displaying the von Neumann mutual information and the maximal Rényi divergence (blue, orange) calculated using brute force diagonalization of the density matrices, and the maximal Rényi divergence calculated using the general case method discussed in Sec. \ref{['sec:general-case-calculation']} (green). The maximal divergence shows a similar correlation to the mutual information as in the case of the 'AEB' architecture.