GME by Seesaw algorithm¶
Geometric Measure of Entanglement (GME) is an important measure of entanglement in quantum information theory. GME of an entangled state is strictly positive, and it is zero for separable states. It is a measure of the distance of a state from the set of separable states (SEP, see doc-link for more details).
This tutorial will reproduce the seesaw algorithm for calculating GME, proposed in paper doi-link
Simple algorithm for computing the geometric measure of entanglement
import numpy as np
import scipy.special
from tqdm.notebook import tqdm
import matplotlib.pyplot as plt
try:
import numqi
except ImportError:
%pip install numqi
import numqi
Pure states¶
For a given pure state $|\psi\rangle$, the GME is given by
$$ E_{G}\left(|\psi\rangle\right)=\max_{\sigma\in\mathrm{SEP}}F\left(\sigma,|\psi\rangle\right)=\max_{|\phi\rangle\in\mathrm{ext}\left(\mathrm{SEP}\right)}\left|\langle\phi|\psi\rangle\right|^{2} $$
where $\mathrm{ext}(\mathrm{SEP})$ is the set of extremal states in the set of separable states, i.e., the pure product states.
Specifically, for bipartite pure states $|\psi\rangle_{AB}$, the GME is given by
$$ |\psi\rangle\in\mathcal{H}_{d_{1}}\otimes\mathcal{H}_{d_{2}}\;\Rightarrow\;E_{G}\left(|\psi\rangle\right)=1-\lambda_{\max}^{2} $$
where $\lambda_{\max}$ is the largest singular value (Schmidt coefficient ) of $|\psi\rangle$ (treated as a $d_1$-by-$d_2$ matrix). For example, the GME for the Bell state $|\Phi^{+}\rangle=\frac{1}{\sqrt{2}}\left(|00\rangle+|11\rangle\right)$ is $1-\frac{1}{2}=0.5$.
psi = numqi.state.Bell()
gme = numqi.entangle.get_GME_pure_seesaw(psi.reshape(2,2))[0]
print('Bell state:', psi)
print('GME:', gme)
Bell state: [0.70710678 0. 0. 0.70710678] GME: 0.4999999999999999
But for multipartite pure states, few analytical results are known. The seesaw algorithm is a simple numerical algorithm to calculate GME for multipartite pure states. The GME for Dicke states (basis for symmetric states) are known analytically doi-link, such as:
$$ \left|D_{4,2}\right\rangle =\frac{1}{\sqrt{6}}\left(\left|0011\right\rangle +\left|0101\right\rangle +\left|1001\right\rangle +\left|0110\right\rangle +\left|1010\right\rangle +\left|1100\right\rangle \right) $$
$$ E_G(\left|D_{4,2}\right\rangle)=\frac{5}{8}$$
psi = numqi.state.Dicke(2,2)
gme = numqi.entangle.get_GME_pure_seesaw(psi.reshape(2,2,2,2))[0]
print('Dicke(4,2):', psi)
print('GME:', gme)
Dicke(4,2): [0. 0. 0. 0.40824829 0. 0.40824829 0.40824829 0. 0. 0.40824829 0.40824829 0. 0.40824829 0. 0. 0. ] GME: 0.6250000010538573
Density matrices¶
For mixed states, the GME is given by
$$ E_{G}(\rho)=1-\max_{\sigma\in\mathrm{SEP}}F(\rho,\sigma) $$
where $F(\rho,\sigma)=\mathrm{Tr}\sqrt{\sqrt{\rho}\sigma\sqrt{\rho}}$ is the fidelity between states $\rho$ and $\sigma$. The seesaw algorithm can be extended to mixed states by considering the purification of the SEP states (a rough derivation is put at the end of this tutorial). Below we demonstrate the seesaw algorithm on some examples.
Isotropic state¶
$$ \rho_d(\alpha)=\frac{1}{d^2-d\alpha}I-\frac{\alpha}{d^2-d\alpha}\sum_{ij}\left|ij\right\rangle \left\langle ji\right|,\quad\alpha\in\left[-1,1\right] $$
Analytical GME is known for isotropic states: doi-link (eq-54)
Geometric measure of entanglement and applications to bipartite and multipartite quantum states
dim = 3
alpha_list = np.linspace(0, 1, 100)
ret_ = numqi.state.get_Isotropic_GME(dim, alpha_list)
kwargs = dict(dim_list=(dim,dim), converge_eps=1e-8)
ret0 = np.array([numqi.entangle.get_GME_seesaw(numqi.state.Isotropic(dim,x), **kwargs) for x in tqdm(alpha_list)])
fig,ax = plt.subplots()
ax.plot(alpha_list, ret_, 'x', label='analytical')
ax.plot(alpha_list, ret0, label='seesaw')
ax.set_xlabel('alpha')
ax.set_title(f'Isotropic(d={dim}), GME')
ax.legend()
ax.set_yscale('log')
fig.tight_layout()
GHZ state¶
$$ |\mathrm{GHZ}\rangle=\frac{1}{\sqrt{2}}\left(|000\rangle+|111\rangle\right) $$
$$ \rho(p)=(1-p) \rho_0 + p |\mathrm{GHZ}\rangle\langle\mathrm{GHZ}| $$
where $\rho_0=\frac{I}{8}$ is the maximally mixed state.
The boundary of the set of separable states for the GHZ state is known: $p=\frac{1}{5}$ doi-link (eq-38)
tmp0 = numqi.state.GHZ(3)
rho_ghz = tmp0.reshape(-1,1)*tmp0.conj()
kwargs = dict(dim_list=(2,2,2), converge_eps=1e-10)
plist = np.linspace(0, 1, 100)
ret_seesaw = [numqi.entangle.get_GME_seesaw(numqi.utils.hf_interpolate_dm(rho_ghz, alpha=x), **kwargs) for x in tqdm(plist)]
fig,ax = plt.subplots()
ax.plot(plist, ret_seesaw)
ax.set_xlabel('p')
ax.set_title('GHZ, GME')
ax.axvline(0.2, color='red')
ax.set_yscale('log')
fig.tight_layout()
XX model¶
XX model is also discussed in the paper, let's reproduce the results.
The Hamiltonian $H$ of three-qubits system is given by
$$ H=\frac{B}{2}\sum_{i=1}^{3}Z_{i}+J\sum_{i=1}^{3}\left(X_{i}X_{i+1}+Y_{i}Y_{i+1}\right) $$
with periodic boundary condition $X_{4}=X_{1},Y_{4}=Y_{1}$. The thermal state is given by
$$ \rho=\frac{1}{Z}e^{-H/kT},\;Z=\mathrm{Tr}\left[e^{-H/kT}\right]. $$
def get_XX_model_periodic_ham(num_qubit):
import functools
hf_kron = lambda *x: functools.reduce(np.kron, x)
ham0 = 0
for ind0 in range(num_qubit):
tmp0 = [numqi.gate.I for _ in range(num_qubit)]
tmp0[ind0] = numqi.gate.Z
ham0 = ham0 + hf_kron(*tmp0)
ham1 = 0
for ind0 in range(num_qubit):
tmp0 = [numqi.gate.I for _ in range(num_qubit)]
tmp0[ind0] = numqi.gate.X
tmp0[(ind0+1)%num_qubit] = numqi.gate.X
ham1 = ham1 + hf_kron(*tmp0)
tmp0 = [numqi.gate.I for _ in range(num_qubit)]
tmp0[ind0] = numqi.gate.Y
tmp0[(ind0+1)%num_qubit] = numqi.gate.Y
ham1 = ham1 + hf_kron(*tmp0)
return ham0, ham1
interactionJ = 0.5
kT_list = np.linspace(1e-4, 2, 100)
magneticB_list = [0, 0.5, 1, 1.2]
num_qubit = 3
kwargs = dict(dim_list=[2]*num_qubit, converge_eps=1e-7)
ham0, ham1 = get_XX_model_periodic_ham(num_qubit)
ret_list = []
for magneticB in magneticB_list:
tmp0 = (magneticB/2)*ham0 + interactionJ*ham1
EVL,EVC = np.linalg.eigh(tmp0)
for kT in tqdm(kT_list, desc=f'magneticB={magneticB}'):
rho = (EVC * scipy.special.softmax(-EVL/kT)) @ EVC.T.conj()
ret_list.append(numqi.entangle.get_GME_seesaw(rho, **kwargs))
ret_list = np.array(ret_list).reshape(-1, len(kT_list))
fig,ax = plt.subplots()
for ind0 in range(len(magneticB_list)):
ax.plot(kT_list, ret_list[ind0], label=f'B={magneticB_list[ind0]}')
ax.set_xlabel('kT')
ax.set_title('XX model, GME')
ax.legend()
fig.tight_layout()
interactionJ = 0.5
kT_list = [0.01, 0.1, 0.5]
magneticB_list = np.linspace(0, 2, 100)
num_qubit = 3
kwargs = dict(dim_list=[2]*num_qubit, converge_eps=1e-7)
ham0, ham1 = get_XX_model_periodic_ham(num_qubit)
ret_list = []
for magneticB in tqdm(magneticB_list):
tmp0 = (magneticB/2)*ham0 + interactionJ*ham1
EVL,EVC = np.linalg.eigh(tmp0)
for kT in kT_list:
rho = (EVC * scipy.special.softmax(-EVL/kT)) @ EVC.T.conj()
ret_list.append(numqi.entangle.get_GME_seesaw(rho, **kwargs))
ret_list = np.array(ret_list).reshape(-1, len(kT_list)).T
fig,ax = plt.subplots()
for ind0 in range(len(kT_list)):
ax.plot(magneticB_list, ret_list[ind0], label=f'kT={kT_list[ind0]}')
ax.set_xlabel('magnetic B')
ax.set_title('XX model, GME')
ax.legend()
fig.tight_layout()
Mathematical formulation¶
This part is mainly for a quick reference. The detailed derivation please refer to the original paper doi-link.
$$ \begin{aligned} \mathrm{Pure}\left(\mathrm{SEP}\right)&=\left\{ |\phi\rangle\in\mathcal{H}_{d_{1}}\otimes\mathcal{H}_{d_{2}}\otimes\cdots\otimes\mathcal{H}_{d_{a}},\mathrm{Tr}_{d_{a}}\left[|\phi\rangle\langle\phi\right]\in\mathrm{SEP}\right\} \\&=\left\{ \sum_{j\in[N]}\sqrt{q_{j}}\left(\otimes_{s}|\phi_{j}^{(s)}\rangle\right)\otimes U^{\dagger}|j\rangle:q\in\Delta_{+}^{N-1},|\phi_{i}\rangle\in\mathcal{H}_{d_{i}},U\in\mathrm{SU}(N)\right\} \end{aligned} $$
$$ \Delta_{+}^{N-1}=\left\{ x\in\mathbb{R}^{N}:x_{i}\geq0,x_{1}+\cdots+x_{N}=1\right\} $$
$$ \rho=\sum_{i\in[r]}p_{i}|\psi_{i}\rangle\langle\psi_{i}| $$
$$ |\psi\rangle=\sum_{i\in[r]}\sqrt{p_{i}}|\psi_{i}\rangle\otimes|i\rangle\;\rightarrow\;\rho=\mathrm{Tr}_{d_{a}}\left[|\psi\rangle\langle\psi|\right] $$
$$ E_{G}(\rho)=1-\max_{|\phi\rangle\in\mathrm{Pure}\left(\mathrm{SEP}\right)}\left|\langle\psi|\phi\rangle\right|^{2} $$