Tapping into Permutation Symmetry for Improved Detection of k-Symmetric Extensions

Abstract

Symmetric extensions are essential in quantum mechanics, providing a lens through which to investigate the correlations of entangled quantum systems and to address challenges like the quantum marginal problem. Though semi-definite programming (SDP) is a recognized method for handling symmetric extensions, it struggles with computational constraints, especially due to the large real parameters in generalized qudit systems. In this study, we introduce an approach that adeptly leverages permutation symmetry. By fine-tuning the SDP problem for detecting k-symmetric extensions, our method markedly diminishes the searching space dimensionality and trims the number of parameters essential for positive-definiteness tests. This leads to an algorithmic enhancement, reducing the complexity from $O(d^{2k})$ to $O(k^{d^2})$ in the qudit k-symmetric extension scenario. Additionally, our approach streamlines the process of verifying the positive definiteness of the results. These advancements pave the way for deeper insights into quantum correlations, highlighting potential avenues for refined research and innovations in quantum information theory.

1. Introduction

In the intricate domain of quantum mechanics, symmetric extensions stand out as a cornerstone, providing a structured mathematical lens through which to explore the nature and behavior of quantum states. A bipartite state ${\rho }_{AB}$ is deemed symmetrically extendible if there exists a multi-partite density matrix ${\rho }_{A_1{A}_2\dots A_m{B}_1{B}_2\dots B_n}$ such that each of its reduced density matrices, when traced over its complements, matches ${\rho }_{AB}$:

$${\mathrm{Tr}}_{{\left(A_j{B}_k\right)}^c}\left({\rho }_{A_1{A}_2\dots A_m{B}_1{B}_2\dots B_n}\right)={\rho }_{AB},\forall j,k.$$

(1)

Delving into the importance of symmetric extensions, they serve as a tangible framework to probe the nature of quantum entanglement, offering a means to understand the profound correlations present in entangled quantum systems [1,2,3]. Furthermore, they pave the way for addressing the quantum marginal problem, which investigates the necessary conditions under which a set of density matrices can correspond to a global state [4,5]. This problem’s universality is showcased by its resonance with the N-representability problem in quantum chemistry [6,7].

A common approach for identifying a k-extension is to cast the problem as a semi-definite programming (SDP) problem [8,9,10]. SDP, a form of convex optimization, involves minimizing a linear function subject to the constraints that the solution lies in the intersection of the cone of positive semi-definite matrices and an affine space. Given that density matrices are inherently semi-definite, SDP has found extensive application in quantum information problems [11,12]. By leveraging the properties of SDP, we can devise efficient algorithms for detecting k-extensions. For example, the algorithm used in QETLAB that determines if a bipartite quantum state ${\rho }_{AB}$ is k-symmetric extendible for the d-dimensional subsystem B has the form:

$$\begin{array}{cc}\hfill & \mathrm{find}\tilde{\rho }\hfill \\ \hfill s.t.& \left\{\begin{array}{c}\tilde{\rho }⪰0\hfill \\ {\mathrm{Tr}}_{{\left(A{B}_1\right)}^c}\left(\tilde{\rho }\right)={\rho }_{AB}\hfill \\ \left({\mathbb{1}}_A\otimes P_{ij}\right)\tilde{\rho }\left({\mathbb{1}}_A\otimes P_{ij}\right)=\tilde{\rho },\forall i,j\hfill \end{array}\right.\hfill \end{array}$$

(2)

where $\tilde{\rho }={\rho }_{A{B}_1{B}_2\dots B_k}$ is the k-symmetric extension of ${\rho }_{AB}$ and the operator $P_{ij}$ is an element in the permutation group $S_k$. However, the substantial number of real parameters, notably in the general qudit scenario, can present formidable computational obstacles. This is largely due to the requirement of the entire extended Hilbert space ${\mathcal{H}}_{A{B}_1{B}_2\dots B_k}$, which scales as $O({(d^2)}^k)$. Such exponential scaling can make calculations intractable for larger systems or higher dimensions.

In this work we are going to present a new optimization scheme, which not only considers the permutation symmetry to reduce the total parameters, but also optimizes the subroutine to determine positive definiteness, where the parameters for single-time optimization grow no faster than

$$\prod _{m=1}^{d-1}\frac{1}{m!}{\left(1+\frac{2k}{d(d+1)}\right)}^{d(d-1)/2},$$

for large k and d.

A testament to our methodology’s effectiveness is its application to the renowned bipartite Werner state, where it exhibits a pronounced acceleration in comparison to the established QETLAB software. This enhancement equips us to approach larger k-extension challenges with unparalleled efficiency.

Furthermore, our calculations have explicitly determined the dimensions of the searching space and the number of parameters required for positive-definiteness tests. This efficiency stems from our algorithm’s ability to undergo multiple distinct positive-definiteness tests, each correlating to a unique Young diagram. Each individual test, though involving a significantly smaller matrix, culminates in a comprehensive and efficient analysis.

Our findings contribute to a clearer understanding of quantum systems, potentially aiding in the design of more proficient quantum algorithms and enhancing our grasp of quantum information theory.

The structure of our paper is as follows: Section 2 delves into the intricacies of the 3-extension of the qutrit case as an illustrative example. Section 3 elucidates our methodology to compute the reduced density matrix of global states for a general k-extendible state, and underscores our rationale for dimensionality reduction. Section 4 shows a comparison of our new algorithm and the traditional one. Concluding insights and discussions are furnished in Section 5.

2. Qutrit Example

Before starting to solve the general problem, we first take a look at two simple examples, 2- and 3-extensions of qutrits. In fact, these two cases clearly demonstrate why our new algorithm can greatly reduce the size of the searching space. We are going to investigate how many real parameters are needed to fully describe the global symmetric extended matrix ${\rho }_{A\overrightarrow{B}}$, which lies in Hilbert space $V=V_A\otimes \mathcal{T}$ constituted by parts A and $\overrightarrow{B}$, with the constraints that ${\mathrm{Tr}}_{{\left(A{B}_1\right)}^c}\left({\rho }_{A\overrightarrow{B}}\right)={\rho }_{AB}$.

2.1. 2-Qutrit Case

In this case, $\mathcal{T}\equiv V^{\left(1\right)}\otimes V^{\left(2\right)}$ is constituted by the two qutrits $B_1$ and $B_2$, where $V^{\left(1\right)}$ and $V^{\left(2\right)}$ represent $B_1$ and $B_2$, respectively. $\mathcal{T}$ is spanned by nine vectors $\{|00\rangle ,|01\rangle ,\cdots ,|22\rangle \}$. According to different permutation symmetry, $\mathcal{T}$ can be decomposed as a 2-invariant orthogonal subspace, a 6-dimensional bosonic one $V^B$, and a 3-dimensional fermionic one $V^F$. It is clear that there does not exist any cross term of bosonic subspace and fermionic subspace, therefore, we only have to consider the density matrix ${\overline{\rho }}_{A\overrightarrow{B}}$ supporting on $\mathrm{End}\left(V_A\right)\otimes \mathrm{End}\left(V^B\right)$ and ${\tilde{\rho }}_{A\overrightarrow{B}}$ supporting on $\mathrm{End}\left(V_A\right)\otimes \mathrm{End}\left(V^F\right)$.

The general form of ${\overline{\rho }}_{A\overrightarrow{B}}$ reads

$${\overline{\rho }}_{A\overrightarrow{B}}=\sum _{\alpha ,\beta =1}^6{\rho }_A^{(\alpha ,\beta )}\otimes {\overline{p}}_{\alpha ,\beta }|{\varphi }_{\overrightarrow{B}}^{\alpha }\rangle \langle {\varphi }_{\overrightarrow{B}}^{\beta }|,$$

(3)

where ${\rho }_A^{(\alpha ,\beta )}\in \mathrm{End}\left(V_A\right)$, ${\overline{p}}_{\alpha ,\beta }$ is a complex number and

$$\begin{array}{ccc}\hfill \{|{\varphi }_{\overrightarrow{B}}^{\alpha }\rangle \}& =& \left\{|00\rangle ,|11\rangle ,|22\rangle ,\frac{1}{\sqrt{2}}(|01\rangle +|10\rangle ),\frac{1}{\sqrt{2}}(|02\rangle +|20\rangle ),\frac{1}{\sqrt{2}}(|12\rangle +|21\rangle ),\right\}.\hfill \end{array}$$

The reduced density matrix can be obtained by performing a partial trace over $B_2$

$$\begin{array}{c}\hfill 2\times {\mathrm{Tr}}_{V^{\left(2\right)}}\left({\overline{\rho }}_{A\overrightarrow{B}}\right)=\sum _{a,b=0}^2{\overline{M}}_{ab}\otimes |a\rangle \langle b|,\end{array}$$

(4)

where ${\overline{M}}_{ab}$ is given by

$$\begin{array}{ccc}\hfill {\overline{M}}_{00}& =& 2{\overline{p}}_{1,1}{\rho }_A^{(1,1)}+{\overline{p}}_{4,4}{\rho }_A^{(4,4)}+{\overline{p}}_{5,5}{\rho }_A^{(5,5)},\hfill \\ \hfill {\overline{M}}_{11}& =& 2{\overline{p}}_{2,2}{\rho }_A^{(2,2)}+{\overline{p}}_{4,4}{\rho }_A^{(4,4)}+{\overline{p}}_{6,6}{\rho }_A^{(6,6)},\hfill \\ \hfill {\overline{M}}_{22}& =& 2{\overline{p}}_{3,3}{\rho }_A^{(3,3)}+{\overline{p}}_{5,5}{\rho }_A^{(5,5)}+{\overline{p}}_{6,6}{\rho }_A^{(6,6)},\hfill \\ \hfill {\overline{M}}_{01}& =& {\overline{M}}_{10}^{†}=\sqrt{2}{\overline{p}}_{4,1}{\rho }_A^{(4,1)}+\sqrt{2}{\overline{p}}_{2,4}{\rho }_A^{(2,4)}+{\overline{p}}_{6,5}{\rho }_A^{(6,5)},\hfill \\ \hfill {\overline{M}}_{02}& =& {\overline{M}}_{20}^{†}=\sqrt{2}{\overline{p}}_{5,1}{\rho }_A^{(5,1)}+\sqrt{2}{\overline{p}}_{3,5}{\rho }_A^{(3,5)}+{\overline{p}}_{6,4}{\rho }_A^{(6,4)},\hfill \\ \hfill {\overline{M}}_{12}& =& {\overline{M}}_{21}^{†}=\sqrt{2}{\overline{p}}_{6,2}{\rho }_A^{(6,2)}+\sqrt{2}{\overline{p}}_{3,6}{\rho }_A^{(3,6)}+{\overline{p}}_{5,4}{\rho }_A^{(5,4)}.\hfill \end{array}$$

(5)

It is noticed that each term of the right-hand side in Equation (4) does not contain every ${\overline{p}}_{\alpha ,\beta }$. In fact, the nonzero coefficients before ${\overline{p}}_{\alpha ,\beta }$ in the term $|a\rangle \langle b|$ are exactly the nonzero entries of the representation matrix $T_{\left(ab\right)}$ over this bosonic invariant subspace (You may find the representation of $su\left(3\right)$ for this case and the following case in standard textbooks on group theory, such as [13,14]).

Similarly, one can write down the general form of ${\tilde{\rho }}_{A\overrightarrow{B}}$ supporting on $\mathrm{End}\left(V_A\right)\otimes \mathrm{End}\left(V^F\right)$

$${\tilde{\rho }}_{A\overrightarrow{B}}=\sum _{\alpha ,\beta =1}^3{\rho }_A^{(\alpha ,\beta )}\otimes {\tilde{p}}_{\alpha ,\beta }|{\psi }_{\overrightarrow{B}}^{\alpha }\rangle \langle {\psi }_{\overrightarrow{B}}^{\beta }|,$$

(6)

where ${\rho }_A^{(\alpha ,\beta )}\in \mathrm{End}\left(V_A\right)$, ${\tilde{p}}_{\alpha ,\beta }$ is a complex number and $\{|{\psi }_{\overrightarrow{B}}^{\alpha }\rangle \}$ is

$$\begin{array}{c}\hfill \left\{\frac{1}{\sqrt{2}}(|01\rangle -|10\rangle ),\frac{1}{\sqrt{2}}(|02\rangle -|20\rangle ),\frac{1}{\sqrt{2}}(|12\rangle -|21\rangle ),\right\}.\end{array}$$

The reduced density matrix can be obtained by performing a partial trace over $B_2$

$$\begin{array}{c}\hfill 2\times {\mathrm{Tr}}_{V^{\left(2\right)}}\left({\tilde{\rho }}_{A\overrightarrow{B}}\right)=\sum _{a,b=0}^2{\tilde{M}}_{ab}\otimes |a\rangle \langle b|,\end{array}$$

(7)

where ${\tilde{M}}_{ab}$ is given by

$$\begin{array}{ccc}\hfill {\tilde{M}}_{00}& =& 2{\tilde{p}}_{1,1}{\rho }_A^{(1,1)}+{\tilde{p}}_{2,2}{\rho }_A^{(2,2)},\hfill \\ \hfill {\tilde{M}}_{11}& =& {\tilde{p}}_{1,1}{\rho }_A^{(1,1)}+{\tilde{p}}_{3,3}{\rho }_A^{(3,3)},\hfill \\ \hfill {\tilde{M}}_{22}& =& {\tilde{p}}_{2,2}{\rho }_A^{(2,2)}+{\tilde{p}}_{3,3}{\rho }_A^{(3,3)},\hfill \\ \hfill {\tilde{M}}_{01}& =& {\tilde{M}}_{10}^{†}={\tilde{p}}_{3,2}{\rho }_A^{(3,2)},\hfill \\ \hfill {\tilde{M}}_{02}& =& {\tilde{M}}_{20}^{†}=-{\tilde{p}}_{3,1}{\rho }_A^{(3,1)},\hfill \\ \hfill {\tilde{M}}_{12}& =& {\tilde{M}}_{21}^{†}={\tilde{p}}_{2,1}{\rho }_A^{(2,1)}.\hfill \end{array}$$

(8)

Taking both ${\overline{\rho }}_{A\overrightarrow{B}}$ and ${\tilde{\rho }}_{A\overrightarrow{B}}$ into account, far fewer real parameters are needed: the original algorithm searches the entire Hilbert space and, thus, the number of parameters is $9^2{d}_A^2$, while usage of permutation symmetry can reduce this number to $(6^2+3^2)d_A^2$. It should be stressed that such a simplification comes from the fact that the cross terms between subspaces corresponding to different permutation symmetries are forbidden.

However, the simple usage of symmetry, such as simply symmetrizing the Gell–Mann matrices over $B_1$ and $B_2$, has to determine the positive definiteness of one matrix with dimensions $(6^2+3^2)d_A^2$. As a comparison, our method involves determining the positive definiteness of two matrices, whose dimensions are $6^2{d}_A^2$ and $3^2{d}_A^2$, respectively.

2.2. 3-Qutrit Case

In this case, $\mathcal{T}\equiv V^{\left(1\right)}\otimes V^{\left(2\right)}\otimes V^{\left(3\right)}$ is constituted by the three qutrits $B_1$, $B_2$, and $B_3$, where $V^{\left(i\right)}$ represents $B_i$. Similar to the procedure in the previous subsection, one can decompose $\mathcal{T}$ as the direct sum of the subspace according to different permutation symmetries, and further more, there does not exist a cross term between subspaces corresponding to different permutation symmetries. The permutation symmetry of the 3-qutrit case is much more complicated than the 2-qutrit case. It is easy to verify that there exist a 10-dimensional bosonic subspace and a 1-dimensional fermionic subspace. Therefore, one can solve this problem by imitating the previous subsection and obtaining the constraint equations. In this situation, the dimensions of the searching space can be reduced from ${27}^2{d}_A^2$ to $({10}^2+{16}^2+1^2)d_A^2$.

However, more room is left for simplification. According to Weyl duality, the 16-dimensional subspace can be further decomposed as two orthogonal 8-dimensional invariant subspaces ${\mathcal{T}}_1^{[2,1]}$ and ${\mathcal{T}}_2^{[2,1]}$, and both are loaded with an equivalent $su\left(3\right)$-irreducible representation, described by a two-row Young diagram $[2,1]$. Here, $\left[\lambda \right]\equiv \{{\lambda }_1,{\lambda }_2,\cdots ,{\lambda }_n\}$ is a partition of integer k, where all ${\lambda }_i$ are integers satisfying ${\lambda }_1\ge {\lambda }_2\ge \cdots \ge {\lambda }_n\ge 0,{\sum }_{i=1}^n{\lambda }_i=k$. Such a partition is denoted by an n-row Young diagram.

${\mathcal{T}}_1^{[2,1]}$ and ${\mathcal{T}}_2^{[2,1]}$ are spanned by the vectors $\left\{{\phi }_{\overrightarrow{B}}^{\alpha ,\left(1\right)}\right\}$ and $\left\{{\phi }_{\overrightarrow{B}}^{\alpha ,\left(2\right)}\right\}$:

$$\begin{array}{ccc}\hfill & & \left\{\frac{1}{\sqrt{6}}(2|001\rangle -|010\rangle -|100\rangle ),\frac{1}{\sqrt{6}}(2|002\rangle -|020\rangle -|200\rangle ),\frac{1}{\sqrt{6}}(|011\rangle -|101\rangle -2|110\rangle ),\right.\hfill \\ \hfill & & \frac{1}{\sqrt{12}}(2|012\rangle -|021\rangle +2|102\rangle -|120\rangle -|201\rangle -|210\rangle ),\frac{1}{\sqrt{4}}(|021\rangle -|120\rangle +|201\rangle -|210\rangle ),\hfill \\ \hfill & & \left.\frac{1}{\sqrt{6}}(|022\rangle +|202\rangle -2|220\rangle ),\frac{1}{\sqrt{6}}(2|112\rangle -|121\rangle -|211\rangle ),\frac{1}{\sqrt{6}}(|122\rangle +|212\rangle -2|221\rangle )\right\},\hfill \end{array}$$

$$\begin{array}{ccc}\hfill & & \left\{\frac{1}{\sqrt{2}}(|010\rangle -|100\rangle ),\frac{1}{\sqrt{2}}(|020\rangle -|200\rangle ),\frac{1}{\sqrt{2}}(|011\rangle -|101\rangle ),\right.\hfill \\ \hfill & & \frac{1}{\sqrt{4}}(|021\rangle +|120\rangle -|201\rangle -|210\rangle ),\frac{1}{\sqrt{12}}(2|012\rangle +|021\rangle -2|102\rangle -|120\rangle -|201\rangle +|210\rangle ),\hfill \\ \hfill & & \left.\frac{1}{\sqrt{2}}(|022\rangle -|202\rangle ),\frac{1}{\sqrt{2}}(|121\rangle -|211\rangle ),\frac{1}{\sqrt{2}}(|122\rangle -|212\rangle ),\right\}.\hfill \end{array}$$

Under the constraint condition imposed in Equation (1), the general form of global state ${\widehat{\rho }}_{A\overrightarrow{B}}$ supporting on $\mathrm{End}\left(V_A\right)\otimes \mathrm{End}({\mathcal{T}}_1^{[2,1]}\oplus {\mathcal{T}}_2^{[2,1]})$ must be of the following form:

$${\widehat{\rho }}_{A\overrightarrow{B}}=\sum _{\alpha ,\beta =1}^8{\rho }_A^{(\alpha ,\beta )}\otimes {\widehat{p}}_{\alpha ,\beta }\left(|{\phi }_{\overrightarrow{B}}^{\alpha ,\left(1\right)}\rangle \langle {\phi }_{\overrightarrow{B}}^{\beta ,\left(1\right)}|+|{\phi }_{\overrightarrow{B}}^{\alpha ,\left(2\right)}\rangle \langle {\phi }_{\overrightarrow{B}}^{\beta ,\left(2\right)}|\right),$$

(9)

The reduced density matrix can be obtained by performing a partial trace over $B_2$

$$\begin{array}{ccc}\hfill & & \frac{3}{2}\times {\mathrm{Tr}}_{V^{\left(2\right)}}\left({\widehat{\rho }}_{A\overrightarrow{B}}\right)=\sum _{a,b=0}^2{\widehat{M}}_{ab}\otimes |a\rangle \langle b|,\hfill \end{array}$$

(10)

where ${\tilde{M}}_{ab}$ is given by

$$\begin{array}{ccc}\hfill {\widehat{M}}_{00}& =& 2{\widehat{p}}_{1,1}{\rho }_A^{(1,1)}+2{\widehat{p}}_{2,2}{\rho }_A^{(2,2)}+{\widehat{p}}_{1,1}{\rho }_A^{(1,1)}+{\widehat{p}}_{4,4}{\rho }_A^{(4,4)}+p_{5,5}{\rho }_A^{(5,5)}+{\widehat{p}}_{6,6}{\rho }_A^{(6,6)},\hfill \\ \hfill {\widehat{M}}_{11}& =& {\widehat{p}}_{1,1}{\rho }_A^{(1,1)}+2{\widehat{p}}_{3,3}{\rho }_A^{(3,3)}+{\widehat{p}}_{4,4}{\rho }_A^{(4,4)}{\widehat{p}}_{5,5}{\rho }_A^{(5,5)}+2{\widehat{p}}_{7,7}{\rho }_A^{(7,7)}+{\widehat{p}}_{8,8}{\rho }_A^{(8,8)},\hfill \\ \hfill {\widehat{M}}_{22}& =& {\widehat{p}}_{2,2}{\rho }_A^{(2,2)}+{\widehat{p}}_{4,4}{\rho }_A^{(4,4)}+{\widehat{p}}_{5,5}{\rho }_A^{(5,5)}+2{\widehat{p}}_{6,6}{\rho }_A^{(6,6)}+{\widehat{p}}_{7,7}{\rho }_A^{(7,7)}+2{\widehat{p}}_{8,8}{\rho }_A^{(8,8)},\hfill \\ \hfill {\widehat{M}}_{01}={\widehat{M}}_{10}^{†}& =& -\frac{1}{\sqrt{2}}{\widehat{p}}_{4,1}{\rho }_A^{(4,1)}+\sqrt{\frac{3}{2}}{\widehat{p}}_{5,1}{\rho }_A^{(5,1)}+{\widehat{p}}_{6,2}{\rho }_A^{(6,2)}+\frac{1}{\sqrt{2}}{\widehat{p}}_{8,4}{\rho }_A^{(8,4)}-\sqrt{\frac{3}{2}}{\widehat{p}}_{8,5}{\rho }_A^{(8,5)}-{\widehat{p}}_{7,3}{\rho }_A^{(7,3)},\hfill \\ \hfill {\widehat{M}}_{02}={\widehat{M}}_{20}^{†}& =& {\widehat{p}}_{2,1}{\rho }_A^{(2,1)}+\frac{1}{\sqrt{2}}{\widehat{p}}_{4,3}{\rho }_A^{(4,3)}+\sqrt{\frac{3}{2}}{\widehat{p}}_{5,3}{\rho }_A^{(5,3)}+\frac{1}{\sqrt{2}}{\widehat{p}}_{6,4}{\rho }_A^{(6,4)}+\sqrt{\frac{3}{2}}{\widehat{p}}_{6,5}{\rho }_A^{(6,5)}+{\widehat{p}}_{2,1}{\rho }_A^{(2,1)},\hfill \\ \hfill {\widehat{M}}_{12}={\widehat{M}}_{21}^{†}& =& {\widehat{p}}_{3,1}{\rho }_A^{(3,1)}+\sqrt{2}{\widehat{p}}_{4,2}{\rho }_A^{(4,2)}+{\widehat{p}}_{7,4}{\rho }_A^{(7,4)}+{\widehat{p}}_{8,7}{\rho }_A^{(8,7)}.\hfill \end{array}$$

(11)

Due to the permutation requirement, the number of real parameters is less than $({10}^2+8^2+1^2)d_A^2$, since some pairs of $\alpha$ and $\beta$ may contribute nothing when computing the one-body reduced density matrix (but these numbers cannot be set to 0 directly, since they may affect the positive definiteness). It should be stressed that the simplification comes not only from the fact that the cross terms between subspaces corresponding to different permutation symmetries are forbidden, but also arises from the fact that the majority of the cross terms within subspaces corresponding to identical permutation symmetries are also forbidden. It should also be noticed that, via our method, one can check the positive definiteness of the global state by successively checking the positive definiteness of the density matrix corresponding to different permutation symmetries.

3. Complexity of Improved SDP

In this section, we are going to give the general form of a global state that corresponds to the given quantum marginals ${\rho }_{AB}$.

Consider the symmetric extension problem described in Equation (2). It is required that the global state ${\rho }_{A{B}_1\cdots B_k}$ is invariant under any exchange of $B_i$ and $B_j$, but it does not require that ${\rho }_{A{B}_1\cdots B_k}$ must support on a subspace with specific permutation symmetry. E.g., for a 2-symmetric extendible state, its extension can be bosonic, which supports on the symmetric subspace only, or fermionic, whose support only resides on the antisymmetric subspace, or more generally, can be a mixture of both.

Consider a Hilbert space $\mathcal{T}={⨂}_{i=1}^k{V}^{\left(i\right)}$ constituted by the k qudits $B_1,B_2,\cdots ,B_k$, whose computational basis is $\left\{{\mathsf{\Phi }}_{i_1,i_2,\cdots ,i_k}\equiv \left|i_1,i_2,\cdots ,i_k\right.〉\right\}$, where $i_1,i_2,\cdots ,i_k=0,1,\cdots ,d-1$.

Each subsystem $V^{\left(i\right)}$ is invariant under $SU\left(d\right)$ ‘rotation’, and transforms according to the d-dimensional fundamental irreducible representation $D^{\left[1\right]}$, which corresponds to the Young diagram $\left[1\right]$ (this is the one block Young diagram □).

Therefore, the Lie algebra $su\left(d\right)$, which is constituted by three series of zero-trace Hermitian matrices and describes the infinitely small rotation of $SU\left(d\right)$, has the following matrix forms on each identical $V^{\left(j\right)}$ if we set $\left|i\right.〉$ to be the natural basis,

$$\begin{array}{ccc}\hfill {\left(T_{mn}^{\left(1\right)}\right)}_{st}& =& \frac{1}{2}({\delta }_{ms}{\delta }_{nt}+{\delta }_{ns}{\delta }_{mt}),\hfill \\ \hfill {\left(T_{mn}^{\left(2\right)}\right)}_{st}& =& \frac{-i}{2}({\delta }_{ms}{\delta }_{nt}+{\delta }_{ns}{\delta }_{mt}),\hfill \\ \hfill {\left(T_p^{\left(3\right)}\right)}_{st}& =& \left\{\begin{array}{cc}\hfill ll{\delta }_{st}{\left[2(p+1)p\right]}^{-\frac{1}{2}},& s<p,\\ \hfill -{\delta }_{st}{[p/(2p+2)]}^{\frac{1}{2}},& s=p,\\ \hfill 0,& s>p,\end{array}\right.,\hfill \end{array}$$

where $m<n,$ and $1\le p\le d-1$. Taking the global phase into account, one should also include the identity matrix. Therefore, one can obtain a new basis for Lie algebra $u\left(d\right)$ by

$$\{T_{ab}|{\left(T_{ab}\right)}_{st}={\delta }_{as}{\delta }_{bt},0\le a,b\le d-1\}.$$

$\mathcal{T}$ is also invariant under the global $U\left(d\right)$ transformation, whose corresponding Lie algebra is given by $\{{\mathbf{T}}_{ab}|{\mathbf{T}}_{ab}={\sum }_i{T}_{ab}^{\left(i\right)}\}$ (here, $T_{ab}^{\left(i\right)}$ denotes that the i-th subsystem transforms according to $T_{ab}$ while others transform according to the identity operator). $\mathcal{T}$ transforms under representation ${\otimes }^k\left[1\right]$, which is not irreducible, but can be decomposed as the direct sum of a series of irreducible representations,

$$\begin{array}{c}\hfill {⨂}^k{D}^{\left[1\right]}=\underset{\left[\lambda \right]}{⨁}{m}_{\left[\lambda \right]}{D}^{\left[\lambda \right]},\end{array}$$

(12)

where $m_{\left[\lambda \right]}$ is the multiplicity of the irreducible representation $D^{\left[\lambda \right]}$. This is equivalent to saying that $\mathcal{T}$ can be partitioned as the direct sum of subspaces. (Please note that subspaces corresponding to different Young diagrams are orthogonal to each other, while those corresponding to the same Young diagram are not. However, it is guaranteed that the intersection of two such different subspaces is zero.)

$$\mathcal{T}=\underset{\left[\lambda \right]}{⨁}{m}_{\left[\lambda \right]}{\mathcal{T}}^{\left[\lambda \right]}.$$

(13)

It can be easily shown that such a ${\mathcal{T}}^{\left[\lambda \right]}$ has a particular permutation symmetry described by Young diagram $\left[\lambda \right]$, and the multiplicity $m_{\left[\lambda \right]}$ equals the dimension of the irreducible representation of $S_k$ corresponding to the identical Young diagram $\left[\lambda \right]$, which gives the equation

$$d^k=\sum _{\left[\lambda \right]}{m}_{\left[\lambda \right]}{D}^{\left[\lambda \right]}.$$

(14)

Two irreducible representation spaces ${\mathcal{T}}_{\mu }^{\left[\lambda \right]}$ and ${\mathcal{T}}_{\nu }^{\left[\lambda \right]}$ corresponding to the same Young diagram but different Young tableaus are orthogonal to each other. Although there might probably be multiplicity in some weight subspace for a general irreducible subspace, one can uniquely label a vector within an arbitrary given irreducible subspace by its weight $\overrightarrow{\omega }$ in $su\left(d\right)$ and the subgroup chain $su\left(d\right)\supset su(d-1)\supset \cdots \supset su\left(2\right)$ [15]. Thus, one can safely use the weight $\overrightarrow{\omega }$ to label different states inside an irreducible subspace ${\mathcal{T}}_{\mu }^{\left[\lambda \right]}$. Therefore, $\left\{\left|\left[\lambda \right],\mu ,\overrightarrow{\omega }\right.〉\right\}$ labels a complete basis of $\mathcal{T}$ one by one, where $\left[\lambda \right]$ denotes non-equivalent $su\left(d\right)$ representations, and $\mu$ differentiates equivalent ones. Together they determine an orthogonal irreducible subspace, and $\overrightarrow{\omega }$ labels every different vector inside.

On the other hand, $\left\{\left|\left[\lambda \right],\mu ,\omega \right.〉\right\}$ can be interpreted in another way: $\overrightarrow{\omega }$ describes the weight, and $\left[\lambda \right]$ denotes non-equivalent $S_k$ representations, thus, these two parameters differentiate orthogonal invariant subspaces, while $\mu$ labels vectors inside. From now on, we shall use $\left|{\overrightarrow{\omega }}_{\mu }^{\left[\lambda \right]}\right.〉$ as short for $\left|\overrightarrow{\omega },\left[\lambda \right],\mu \right.〉$.

Any matrix ${\rho }_{A{B}_1{B}_2\cdots B_n}\in \mathrm{End}\left(V_A\right)\otimes \mathrm{End}\left(\mathcal{T}\right)$ can be expressed as

$${\rho }_{A{B}_1{B}_2\cdots B_k}=\sum _{\alpha ,{\alpha }^{\prime }}\sum _{\left[\lambda \right],\left[{\lambda }^{\prime }\right]}\sum _{\mu ,{\mu }^{\prime }}\sum _{\overrightarrow{\omega },\overrightarrow{{\omega }^{\prime }}}\left|{\psi }_{\overrightarrow{\omega },\left[\lambda \right],\mu }^{\alpha }\right.〉\left.〈{\psi }_{\overrightarrow{{\omega }^{\prime }},\left[{\lambda }^{\prime }\right],{\mu }^{\prime }}^{{\alpha }^{\prime }}\right|\otimes \left|{\overrightarrow{\omega }}_{\mu }^{\left[\lambda \right]}\right.〉\left.〈{\overrightarrow{{\omega }^{\prime }}}_{{\mu }^{\prime }}^{\left[{\lambda }^{\prime }\right]}\right|,$$

(15)

where $\left|{\psi }_{\overrightarrow{\omega },\left[\lambda \right],\mu }^{\alpha }\right.〉$ are non-normalized states in $V_A$ and $\alpha$ labels different states in $V_A$.

Inserting Equation (15) into Equation (2), $\forall \pi \in S_k$, we obtain a series of constraints for ${\rho }_{A{B}_1{B}_2\cdots B_k}$:

$$\begin{array}{ccc}& & \forall \left[\lambda \right],\left[{\lambda }^{\prime }\right]\overrightarrow{\omega },\overrightarrow{{\omega }^{\prime }}\mathrm{and}\mu ,{\mu }^{\prime },\hfill \\ & & \sum _{\alpha ,{\alpha }^{\prime }}\left|{\psi }_{\overrightarrow{\omega },\left[\lambda \right],\mu }^{\alpha }\right.〉\left.〈{\psi }_{\overrightarrow{{\omega }^{\prime }},\left[{\lambda }^{\prime }\right],{\mu }^{\prime }}^{{\alpha }^{\prime }}\right|\sum _{\nu ,{\nu }^{\prime }}\mathcal{A}{\left(\pi \right)}_{\mu ,\nu }^{\left[\lambda \right]}\mathcal{A}{\left(\pi \right)}_{{\nu }^{\prime },{\mu }^{\prime }}^{\left[{\lambda }^{\prime }\right]*}\left|{\overrightarrow{\omega }}_{\nu }^{\left[\lambda \right]}\right.〉\left.〈{\overrightarrow{{\omega }^{\prime }}}_{{\nu }^{\prime }}^{\left[{\lambda }^{\prime }\right]}\right|\hfill \\ & & =\sum _{\alpha ,{\alpha }^{\prime }}\left|{\psi }_{\overrightarrow{\omega },\left[\lambda \right],\mu }^{\alpha }\right.〉\left.〈{\psi }_{\overrightarrow{{\omega }^{\prime }},\left[{\lambda }^{\prime }\right],{\mu }^{\prime }}^{{\alpha }^{\prime }}\right|\left|{\overrightarrow{\omega }}_{\mu }^{\left[\lambda \right]}\right.〉\left.〈{\overrightarrow{{\omega }^{\prime }}}_{{\mu }^{\prime }}^{\left[{\lambda }^{\prime }\right]}\right|,\hfill \end{array}$$

(16)

where ${\mathcal{A}}^{\left[\lambda \right]}$ and ${\mathcal{A}}^{\left[{\lambda }^{\prime }\right]}$ are irreducible representations of the permutation group $S_k$.

Define matrix

$$M(\overrightarrow{\omega },\overrightarrow{{\omega }^{\prime }},\left[\lambda \right],\left[{\lambda }^{\prime }\right])\equiv \sum _{\mu ,{\mu }^{\prime }}M{(\overrightarrow{\omega },\overrightarrow{{\omega }^{\prime }},\left[\lambda \right],\left[{\lambda }^{\prime }\right])}_{\mu {\mu }^{\prime }}\left|{\overrightarrow{\omega }}_{\mu }^{\left[\lambda \right]}\right.〉\left.〈{\overrightarrow{{\omega }^{\prime }}}_{{\mu }^{\prime }}^{\left[{\lambda }^{\prime }\right]}\right|,$$

(17)

where

$$M{(\overrightarrow{\omega },\overrightarrow{{\omega }^{\prime }},\left[\lambda \right],\left[{\lambda }^{\prime }\right])}_{\mu {\mu }^{\prime }}\equiv \sum _{\alpha ,{\alpha }^{\prime }}\left|{\psi }_{\overrightarrow{\omega },\left[\lambda \right],\mu }^{\alpha }\right.〉\left.〈{\psi }_{\overrightarrow{{\omega }^{\prime }},\left[{\lambda }^{\prime }\right],{\mu }^{\prime }}^{{\alpha }^{\prime }}\right|,$$

(18)

thus, $\forall \pi \in S_k$

$${\mathcal{A}}^{\left[\lambda \right]}\left(\pi \right)M(\overrightarrow{\omega },\overrightarrow{{\omega }^{\prime }},\left[\lambda \right],\left[{\lambda }^{\prime }\right]){\mathcal{A}}^{\left[{\lambda }^{\prime }\right]}{\left(\pi \right)}^{†}=M(\overrightarrow{\omega },\overrightarrow{{\omega }^{\prime }},\left[\lambda \right],\left[{\lambda }^{\prime }\right]).$$

(19)

Schur’s lemma guarantees that,

a.

when $\left[\lambda \right]\ne \left[{\lambda }^{\prime }\right]$, $M=0$;

b.

when $\left[\lambda \right]=\left[{\lambda }^{\prime }\right]$, M is invertible.

$\left|{\overrightarrow{\omega }}_{\mu }^{\left[\lambda \right]}\right.〉$ should be chosen carefully such that the representations ${\mathcal{A}}^{\left[\lambda \right]}$ are identical, not just an isomorphic matrix, for different weights $\omega$. Then, all $M(\overrightarrow{\omega },\overrightarrow{{\omega }^{\prime }},\left[\lambda \right],\left[\lambda \right])$ can be proportional to the corresponding identity matrix. Therefore, one can eliminate the majority of cross terms and restrict ${\rho }_{A{B}_1{B}_2\cdots B_k}$ to

$${\rho }_{A{B}_1{B}_2\cdots B_k}=\sum _{\left[\lambda \right]}\sum _{\overrightarrow{\omega },\overrightarrow{{\omega }^{\prime }}}f(\left[\lambda \right],\overrightarrow{\omega },\overrightarrow{{\omega }^{\prime }})\sigma (\left[\lambda \right],\overrightarrow{\omega },\overrightarrow{{\omega }^{\prime }})\otimes \sum _{\mu }\left|{\overrightarrow{\omega }}_{\mu }^{\left[\lambda \right]}\right.〉\left.〈{\overrightarrow{\omega }}_{\mu }^{\prime \left[\lambda \right]}\right|,$$

(20)

where $f(\left[\lambda \right],\overrightarrow{\omega },\overrightarrow{{\omega }^{\prime }})$ is the coefficient and $\sigma (\left[\lambda \right],\overrightarrow{\omega },\overrightarrow{{\omega }^{\prime }})\in \mathrm{End}\left(V_A\right)$ (this does not have to be a density matrix!), both of which correspond to the $S_k$-irreducible representation described by the Young diagram $\left[\lambda \right]$ and different weights $\overrightarrow{\omega }$ and $\overrightarrow{{\omega }^{\prime }}$.

Our next task is to determine the RDM of the global state given by Equation (20). For every given $\left[\lambda \right],\overrightarrow{\omega }$, and $\overrightarrow{{\omega }^{\prime }}$, one can temporally ignore system A and concentrate on group $\{B_1,B_2,\cdots ,B_k\}$,

$$\begin{array}{ccc}\hfill \sum _{i,j=0}^{d-1}{B}_{ij}\left|i\right.〉\left.〈j\right|& =& {\mathrm{Tr}}_{B_1^c}\left(\sum _{\mu }\left|{\overrightarrow{\omega }}_{\mu }^{\left[\lambda \right]}\right.〉\left.〈{\overrightarrow{\omega }}_{\mu }^{\prime \left[\lambda \right]}\right|\right)\hfill \\ & =& \sum _{i,j=0}^{d-1}\left|i\right.〉\left.〈j\right|\mathrm{Tr}\left(T_{ji}^{\left(1\right)}\sum _{\mu }\left|{\overrightarrow{\omega }}_{\mu }^{\left[\lambda \right]}\right.〉\left.〈{\overrightarrow{\omega }}_{\mu }^{\prime \left[\lambda \right]}\right|\right)\hfill \\ & =& \sum _{i,j=0}^{d-1}\left|i\right.〉\left.〈j\right|\mathrm{Tr}\left(\frac{1}{k}{\mathbf{T}}_{ji}\sum _{\mu }\left|{\overrightarrow{\omega }}_{\mu }^{\left[\lambda \right]}\right.〉\left.〈{\overrightarrow{\omega }}_{\mu }^{\prime \left[\lambda \right]}\right|\right)\hfill \\ & =& \sum _{i,j=0}^{d-1}\left|i\right.〉\left.〈j\right|\left(\frac{m_{\left[\lambda \right]}}{k}\left.〈{\overrightarrow{\omega }}^{\prime \left[\lambda \right]}\right|{\mathbf{T}}_{ji}\left|{\overrightarrow{\omega }}^{\left[\lambda \right]}\right.〉\right),\hfill \end{array}$$

(21)

where $\left.〈{\overrightarrow{\omega }}^{\prime \left[\lambda \right]}\right|{\mathbf{T}}_{ji}\left|{\overrightarrow{\omega }}^{\left[\lambda \right]}\right.〉$ is exactly the matrix element of irreducible representation corresponding to Young diagram $\left[\lambda \right]$ for generator $T_{ji}$ in Lie algebra $u\left(d\right)$ (do not worry about this part, the general matrix form of ${\mathbf{T}}_{ji}$ in irreducible representation $D^{\left[\lambda \right]}$ of $u\left(d\right)$ has been calculated by mathematicians and you can refer to [15]). Taking subsystem A into account, one can obtain

$$\begin{array}{ccc}\hfill & & {\mathrm{Tr}}_{{\left(A{B}_1\right)}^c}\left({\rho }_{A{B}_1{B}_2\cdots B_k}\right)\hfill \\ & =& \sum _{\left[\lambda \right]}\sum _{m,n=0}^{d_A-1}\sum _{i,j=0}^{d-1}\left|m\right.〉\left.〈n\right|\otimes \left|i\right.〉\left.〈j\right|\hfill \\ & & \times \sum _{\overrightarrow{\omega },{\overrightarrow{\omega }}^{\prime }}\sigma {(\left[\lambda \right],\overrightarrow{\omega },\overrightarrow{{\omega }^{\prime }})}_{mn}\frac{m_{\left[\lambda \right]}}{k}\left.〈{\overrightarrow{\omega }}^{\prime \left[\lambda \right]}\right|{\mathbf{T}}_{ji}\left|{\overrightarrow{\omega }}^{\left[\lambda \right]}\right.〉f(\left[\lambda \right],\overrightarrow{\omega },\overrightarrow{{\omega }^{\prime }}).\hfill \end{array}$$

(22)

For every given $\left[\lambda \right]$, the number of different values for $\overrightarrow{\omega }$ and $\overrightarrow{{\omega }^{\prime }}$ is just the dimension of the $u\left(d\right)$-irreducible representation $D^{\left[\lambda \right]}$. Therefore, ignoring the size of subsystem A, the size of the searching space in dealing with symmetric extension is given by

$$\sum {\left(D^{\left[\lambda \right]}\right)}^2=\left(\genfrac{}{}{0pt}{}{d^2-1+k}{k}\right)<d^{2k},$$

(23)

where the summation runs over all possible proper Young diagrams. One may conclude that the dimensions of the entire searching space grow no faster than $O\left(k^{d^2}\right)$, which is significantly smaller than the original $O\left(d^{2k}\right)$, therefore, the efficiency of SDP can be greatly improved when dealing with symmetric extension problems.

To guarantee the positive definiteness of the solved global matrix, one can test whether the density matrix corresponding to a different permutation is positive definite, and hence each testing needs far fewer resources.

To investigate the amount of resources needed for each single testing, one should focus on the growth rates of $D^{\left[\lambda \right]}$ and $m_{\left[\lambda \right]}$. The asymptotic behavior of the upper bound of $D^{\left[\lambda \right]}$ is given by

$$\prod _{m=1}^{d-1}\frac{1}{\left(m\right)!}{\left(1+\frac{2k}{d(d+1)}\right)}^{\frac{d(d-1)}{2}},$$

(24)

which corresponds to the irreducible representation that satisfies ${\lambda }_i-{\lambda }_{i+1}\approx 2k/d(d-1)$ [16]. For a given k, the number of different valid Young diagrams whose number of rows is less than or equal to d is hard to compute analytically, but for a sufficiently large k, the asymptotic value is $\frac{1}{d!}\left(\genfrac{}{}{0pt}{}{k+d-1}{k}\right)$. To find a partition (not necessarily a partition that corresponds to a valid Young diagram!) that satisfies ${\sum }_i^k{\lambda }_i=k$ is equivalent to inserting $d-1$ separators between a line of k balls, which reads $\left(\genfrac{}{}{0pt}{}{k+d-1}{k}\right)$. Since some ${\lambda }_i$’s might be identical, the number of valid different Young diagrams should be less than $\frac{1}{d!}\left(\genfrac{}{}{0pt}{}{k+d-1}{k}\right)$, but when k is sufficiently large, the odds of any ${\lambda }_i$’s being identical approaches 0, so the asymptotic number of valid Young diagrams is given by $\frac{1}{d!}\left(\genfrac{}{}{0pt}{}{k+d-1}{k}\right)$.

4. Numerical Results

First, we apply our algorithm to the famous bipartite Werner state ${\rho }_{W,d}\left(\alpha \right)\in {\mathcal{H}}_d\otimes {\mathcal{H}}_d$

$${\rho }_{W,d}\left(\alpha \right)=\frac{1}{d^2-d\alpha }I-\frac{\alpha }{d^2-d\alpha }\sum _{ij}\left|ij\right.〉\left.〈ji\right|,\alpha \in \left[-1,1\right].$$

Previous work [17] proved that the Werner state is $(1,k)$-extendible for $\alpha \in [-1,\frac{k+d^2-d}{kd+d-1}]$. As k goes to infinity, it gives the separable Werner state $\alpha \in \left[-1,\frac{1}{d}\right]$. To obtain such a $(1,k)$-extendible boundary ${\alpha }_k^{*}$, one can solve the following semi-definite programming

$$\begin{array}{cc}\hfill & maxc,\hfill \\ \hfill s.t.& \left\{\begin{array}{c}{\rho }_{A{B}_1\cdots B_k}⪰0,\hfill \\ \left({\mathbb{1}}^A\otimes P_{ij}\right){\rho }_{A{B}_1\cdots B_k}\left({\mathbb{1}}^A\otimes P_{ij}\right)={\rho }_{A{B}_1\cdots B_k},\hfill \\ {\mathrm{Tr}}_{B_1^c}\left({\rho }_{A{B}_1\cdots B_k}\right)=\left(1-c\right){\rho }_o+c{\rho }_{W,d}\left(1\right),\hfill \end{array}\right.\hfill \end{array}$$

where ${\rho }_0$ denotes the maximally mixed state (as semi-definite programming requires linear or affine equation constraints, we convert the non-linear expression $\alpha$ in the Werner state into a linear interpolation $(1-c){\rho }_0+c{\rho }_{W,d}\left(1\right)$ used in the optimization), and the boundary can be calculated from the optimal value ${\alpha }_k^{*}=\frac{c^{*}d}{c^{*}+d-1}$.

The results are shown in

Table 1. Time usage for calculating the Werner $(1,k)$-extendible boundary. The dashed line “-” indicates the optimization failed due to memory limitation or intolerable time usage.

($\mathit{d},\mathit{k}$)	QETLAB (s)	irrep (s)	${\mathit{\alpha }}_{\mathit{k}}^{*}$
$(2,8)$	0.19	0.16	0.588
$(2,10)$	12.60	0.16	0.571
$(2,16)$	-	0.32	0.545
$(2,32)$	-	3.18	0.523
$(2,32)$	-	51.96	0.512
$(3,3)$	0.62	0.51	0.818
$(3,4)$	7.96	2.38	0.714
$(3,5)$	-	11.56	0.647
$(3,6)$	-	55.60	0.6

. We compare the time required with the software QETLAB [18], a widely used MATLAB package in the quantum information community. The benchmark is performed on a standard laptop, AMD R7-5800H, 16 CPU cores (hyperthread enabled), 16 GB memory, and our algorithm is implemented in the CVXPY package [19] with the MOSEK solver [20]. The solved boundary ${\alpha }_k^{*}$ is within an absolute error of ${10}^{-8}$ compared with the analytical results. From the results, a significant speedup can be observed and a much larger k-extension problem can be handled by our algorithm.

We explicitly calculate the dimensions of the searching space and the number of parameters required to be tested for positive definiteness, which demonstrates the efficiency of our algorithm, as shown in

Table 2. Dimensions and number of parameters needed in positive definiteness.

($\mathit{d},\mathit{k}$)	#Searching Space (QETLAB, irrep)	#Positive Definiteness (QETLAB, irrep)
$(3,3)$	$(729,165)$	$(729, {10}^2+8^2+1^2)$
$(3,4)$	$(6561,495)$	$(6561, {15}^2+{15}^2+6^2+3^2)$
$(3,5)$	(59,049, 1287)	(59,049, ${21}^2+{24}^2+{15}^2+6^2+3^2)$
$(3,6)$	(531,441, 3003)	(531,441, ${28}^2+{35}^2+{27}^2+{10}^2+{10}^2+8^2+1^2)$
$(4,3)$	$(4096, 816)$	$(4096, {20}^2+{20}^2+4^2)$
$(4,4)$	(65,536, 3876)	(65,536, ${35}^2+{45}^2+{20}^2+{15}^2+1^2)$
$(4,5)$	(1,048,576, 15,504)	(1,048,576, ${56}^2+{84}^2+{60}^2+{36}^2+{20}^2+4^2)$

(our algorithm needs to undergo multiple different positive-definiteness tests, where each different Young diagram corresponds to its own test, but each individual test involves a significantly smaller matrix, and hence the efficiency is improved).

5. Discussion

The complexity of our new algorithm for dealing with k-symmetric extensions of quantum states is $O\left(k^{d^2}\right)$, which is an improvement over the original algorithm with $O\left(d^{2k}\right)$ complexity. However, it is important to note that the complexity of detecting entanglement is a QMA problem, which means that it is generally considered to be computationally hard. Although our new algorithm reduces the computational complexity of the problem, it does not change the fundamental difficulty of detecting entanglement. This is due to the fact that the size of the input of this problem is given by $O(logk,d)$, and hence the resources needed in our algorithm still grow exponentially relative to the input. Therefore, while our algorithm presents some advance, it does not contradict the known fact that detecting entanglement is a QMA problem. The challenge of detecting entanglement remains an important area of research, with many open questions and opportunities for new breakthroughs.