Pure bosonic extension¶
Given a quantum state $\rho_{AB}$, we know that bosonic extension can be used to compute the lower bounds of REE, as the following:
$$ \min_{\sigma_{AB}} S(\sigma_{AB} \| \rho_{AB})$$
$$ s.t.\; \begin{cases} \sigma_{AB^k} \succeq 0, ~\operatorname{Tr}[\sigma_{AB^k}]&=1,\\ \operatorname{Tr}_{B_2 \cdots B_k}\left[\sigma_{A B^k}\right] &=\sigma_{AB},\\ \left(\mathbb{I}_A \otimes P_{i j}\right) \sigma_{A B^k} & =\sigma_{A B^k},\\ \end{cases} $$
Pure bosonic extension is a special case of bosonic extension, where the optimized global state $\sigma_{AB^k}$ is pure, as the following:
$$ \min_{\sigma_{AB}} S( \sigma_{AB} \| \rho_{AB})$$
$$ s.t.\; \begin{cases} \langle\psi_{A B^k}|\psi_{A B^k}\rangle &=1,\\ \operatorname{Tr}_{B_2 \cdots B_k}\left[|\psi_{A B^k}\rangle \langle\psi_{A B^k}| \right] &=\sigma_{AB},\\ \left(\mathbb{I}_A \otimes P_{i j}\right)|\psi_{A B^k}\rangle & =|\psi_{A B^k}\rangle,\\ \end{cases} $$ where $\sigma_{AB^k}$ is denoted by the pure state $|\psi_{AB^k}\rangle$.
Pure bosonic extension reveals that when $k > k_0$, where $k_0$ is related to the local dimensions of $A$ and $B$, the above two optimization problems are equivalent generically. $k_0$ can be estimated by the following formula:
$$ \left(\begin{array}{c} k_0+d_B-1 \\ d_B-1 \end{array}\right)=\frac{1}{2} d_A d_B^2 $$
Since the pure bosonic state can be parametrized by the following form:
$$ \left|\psi_{A B_1 \ldots B_k}\right\rangle=\sum_{i, \vec{w}} p_{i, \vec{w}}|i\rangle \otimes\left|D_{\vec{w}}^k\right\rangle . $$ where $\left|D_{\vec{w}}^k\right\rangle$ is the generalized Dicke state doc-link. So the corresponding marginal state $\sigma_{AB}$ can be parametrized by the following form:
$$ \begin{aligned} \sigma_{A B} & =\operatorname{Tr}_{B_2 \ldots B_k}\left[\left|\psi_{A B_1 \ldots B_k}\right\rangle\left\langle\psi_{A B_1 \ldots B_k}\right|\right] \\ & =\sum_{i, j, \vec{w}, \vec{w}^{\prime}, r, s} p_{i, \vec{w}} p_{j, \vec{w}^{\prime}}^* \mathcal{B}_{\vec{w} \vec{w}^{\prime} r s}^k|i\rangle\langle j|\otimes| r\rangle\langle s|, \end{aligned} $$ where $\mathcal{B}^k$ is the marginal bosonic coefficient tensor satisfying
$$ \begin{aligned} \mathcal{B}_{\vec{w} \vec{w}^{\prime} r s}^k & =\left\langle r\left|\operatorname{Tr}_{B_2 \cdots B_k}\left[\left|D_{\vec{w}}^k\right\rangle\left\langle D_{\vec{w}^{\prime}}^k\right|\right]\right| s\right\rangle \\ & =\frac{1}{k} \sqrt{w_r w_s^{\prime}} \delta_{w_0, w_0^{\prime}} \ldots \delta_{w_r-1, w_r^{\prime}} \ldots \delta_{w_s, w_s^{\prime}-1} \ldots \delta_{w_{d-1}, w_{d-1}^{\prime}}. \end{aligned} $$
So the optimization problem can be reformulated as the following:
$$ \min_{\sigma_{AB}} S( \sigma_{AB} \| \rho_{AB})$$
$$ s.t.\; \sigma_{A B} = \sum_{i, j, \vec{w}, \vec{w}^{\prime}, r, s} p_{i, \vec{w}} p_{j, \vec{w}^{\prime}}^* \mathcal{B}_{\vec{w} \vec{w}^{\prime} r s}^k|i\rangle\langle j|\otimes| r\rangle\langle s| $$ which is a non-convex optimization over the parameters $p_{i, \vec{w}}$. This can be much more efficient than the SDP method since we only need to optimize over the pure states instead of the whole density matrices.
import numpy as np
try:
import numqi
except ImportError:
%pip install numqi
import numqi
dim = 2
alpha_list = np.linspace(0, 1, 100, endpoint=False) # alpha=1 is unstable for matrix logarithm
dm_target_list = [numqi.state.Werner(dim, x) for x in alpha_list]
# Analytical REE of Werner states
ree_analytical = np.array([numqi.state.get_Werner_ree(dim, x) for x in alpha_list])
# REE of Werner states obtained by PPT
ree_ppt = numqi.entangle.get_ppt_ree(dm_target_list, dim, dim)
kext = [6, 12, 18]
ree_pureb = []
for k in kext:
ree_kext = []
for dm_target in dm_target_list:
model = numqi.entangle.PureBosonicExt(dim, dim, kext=k)
model.set_dm_target(dm_target)
tmp0 = numqi.optimize.minimize(model, theta0='uniform', tol=1e-12, num_repeat=3, print_every_round=0).fun
ree_kext.append(tmp0)
ree_pureb.append(ree_kext)
ree_pureb = np.array(ree_pureb)
import matplotlib.pyplot as plt
fig, ax = plt.subplots()
ax.plot(alpha_list, ree_analytical, "x", label="analytical")
ax.plot(alpha_list, ree_ppt, label="PPT")
for k, ree_kext in zip(kext, ree_pureb):
ax.plot(alpha_list, ree_kext, label=f"pure {k}-bosonic extension")
ax.set_xlim(min(alpha_list), max(alpha_list))
ax.set_ylim(1e-13, 1)
ax.set_yscale('log')
ax.set_xlabel(r"$\alpha$")
ax.set_ylabel("REE")
ax.set_title(f'Werner({dim})')
ax.legend()
fig.tight_layout()
Werner state (d=2) boundary $\alpha^{\star}$
| kext | PureB-ext | qetlab-boson-ext | qetlab-symmetric-ext |
|---|---|---|---|
| 4 | 0.64585 | 0.6668 | 0.6668 |
| 5 | 0.63661 | 0.63650 | 0.63650 |
| 6 | 0.61206 | 0.61556 | 0.61554 |
| 7 | 0.60019 | 0.60016 | 0.60017 |
| 8 | 0.58846 | 0.58840 | NA |
| 9 | 0.57918 | 0.57919 | NA |
| 10 | 0.57165 | 0.57168 | NA |
| 11 | 0.56547 | 0.56548 | NA |
| 12 | 0.56024 | NA | NA |
| 16 | 0.54575 | NA | NA |
| 512 | 0.50175 | NA | NA |
| 8192 | 0.50040 | NA | NA |
| 65536 | 0.50034 | NA | NA |