Self-testing quantum systems of arbitrary local dimension with minimal number of measurements

Bell nonlocality as a resource for device independent certification schemes has been studied extensively in recent years. The strongest form of device independent certification is referred to as self-testing, which given a device certifies the promised quantum state as well as quantum measurements performed on it without any knowledge of the internal workings of the device. In spite of various results on self-testing protocols, it remains a highly nontrivial problem to propose a certification scheme of qudit-qudit entangled states based on violation of a single $d$-outcome Bell inequality. Here we address this problem and propose a self-testing protocol for the maximally entangled state of any local dimension using the minimum number of measurements possible, i.e., two per subsystem. Our self-testing result can be used to establish unbounded randomness expansion, $\log_2d$ perfect random bits, while it requires only one random bit to encode the measurement choice.

Bell nonlocality as a resource for device independent certification schemes has been studied extensively in recent years. The strongest form of device independent certification is referred to as self-testing, which given a device certifies the promised quantum state as well as quantum measurements performed on it without any knowledge of the internal workings of the device. In spite of various results on self-testing protocols, it remains a highly nontrivial problem to propose a certification scheme of qudit-qudit entangled states based on violation of a single d-outcome Bell inequality. Here we address this problem and propose a self-testing protocol for the maximally entangled state of any local dimension using the minimum number of measurements possible, i.e., two per subsystem. Our self-testing result can be used to establish unbounded randomness expansion, log 2 d perfect random bits, while it requires only one random bit to encode the measurement choice.
Introduction.-The advent of quantum theory has not just changed the understanding of physics but has also given rise to new phenomena that would have never been possible in the classical world. Arguably one of the most interesting features of quantum theory is the existence of quantum correlations which cannot be explained by any local hidden variable model, a phenomenon commonly referred to as Bell nonlocality [1,2]. It has been understood that apart from its fundamental significance, nonlocality is a powerful resource for certain device-independent applications such as quantum cryptography [3], randomness generation [4,5], or, more recently, for device-independent certification methods [3][4][5][6].
The strength of the device-independent (DI) certification methods lies in the fact that given an unknown quantum system, one can make nontrivial statements on some of its key features based solely on the nonlocal correlations obtained from it. They are advantageous over standard certification methods such as those based on quantum tomography (cf. Ref. [7]) because they do not require making assumptions on the system under study, apart from that it is governed by quantum mechanics. An example of such a DI certification scheme would be verifying whether a quantum device produces entanglement [8] or certification of the dimension of a quantum system [9], both based on a violation of some Bell inequality by the corresponding quantum system.
The strongest form of DI certification is self-testing. First introduced in [6], it allows one to provide a full description, up to certain well-understood equivalences, of the considered quantum system and also the measurements performed on it based on observing the maximal violation of some Bell inequality. Such a form of certification is particularly interesting from the application point of view as it provides a way of verifying that a given quantum device functions properly without the need of knowing its internal working.
A lot of attention has thus been devoted to proposing self-testing schemes for entangled quantum states [10][11][12][13][14]. However, most of the obtained results focus on states that are locally qubits such as for instance the self-testing statement for any two-qubit entangled state [11,12] based on a violation of the tilted version of the famous Clauser-Horne-Shimony-Holt (CHSH) Bell inequality [15,16]. At the same time, quantum systems of higher local dimension remain mostly unexplored and, in fact, few results are devoted to them [11,[17][18][19]. In particular, in Refs. [17] the two-qubit results [11,12] were combined to design a self-testing protocol for any entangled state of arbitrary local dimension. Still, these results are based on a violation of many two-outcome Bell inequalities and the question whether one can design a self-testing statement for qudit quantum systems relying on violation of a single and truly d-outcome Bell inequality remains open. Moreover, the self-testing statement of Refs. [11,17] is not optimal in terms of the number of measurements that the observers need to perform (three and four, respectively) to certify the state and the corresponding measurements. Taking into account the possibility of experimental implementations of self-testing, it is a highly relevant question whether alternative protocols can be derived which rely on the minimal number of two measurements per observer.
The main aim of this work is to address the above questions. We provide the first self-testing statement for a maximally entangled state of local dimension d: from a violation of a truly d-outcome Bell inequality exploiting the minimal possible number of measurements per party, which is two. To this end, we use the Bell inequality introduced in Ref. [20] as a generalization of the well-known CHSH Bell inequality to d-outcome Bell scenarios. A straightforward implication of our selftesting statement is a novel and simpler, as compared to the previous approaches [21,22], scheme for parallel self-testing of N copies of the two-qubit maximally en-arXiv:1909.12722v2 [quant-ph] 2 Oct 2019 tangled state |φ + 2 . Another implication is that the outcomes of local measurements maximally violating the respective Bell inequality are perfectly random, which allows us to propose a quantum protocol for unbounded expansion of quantum randomness.
Preliminaries.-Before presenting our results let us set up the scenario and introduce the relevant notions.
Bell scenario. The simplest experimental set-up exhibiting quantum nonlocality, namely, the bipartite Bell set-up, comprises one preparation and two measurement devices. The latter are possessed by distant and noncommunicating parties, Alice and Bob, and can in general perform one of m measurements, denoted M x and N y with x, y = 1, . . . , m. In each run of the experiment, a bipartite quantum system ρ AB is prepared by the preparation device and subsequently each measurement device performs a measurement on one of its subsystems, yielding an outcome a and b, respectively, with Correlations obtained in this experiment are captured by a vector of probability distributions Here, p(a, b|x, y) is the probability of obtaining outcomes a and b by Alice and Bob after performing measurements M x and N y , respectively, and it is given by the well-known formula where M (a) x and N (b) y are the measurement operators defining the measurements M x and N y , respectively. In what follows we also refer to the vector p as correlations.
In a fixed Bell scenario with m d-outcome measurements such quantum correlations, that is, those obtained by performing local measurements on composite quantum systems, form a convex set Q. For further purposes it is important to recall that any point p ∈ Q can be achieved with a pure state |ψ AB , whose local dimensions might in general be higher than those of ρ AB , and projective measurements M x = {P (a) It turns out that the quantum set Q contains correlations that even if obtained from a quantum state, can be simulated by Alice and Bob in a purely classical way. Such correlations are said to admit a local-realistic description and for brevity we refer to them as local or classical. More formally, local correlations are those for which p(a, b|x, y) can be represented as a convex combination of product deterministic correlations where λ is a random variable distributed according to a distribution µ(λ) and p(a|x, λ), p(b|y, λ) ∈ {0, 1} for every x, y, a, b. For any m and d such local correlations form a convex polytope that we denote L.
Bell inequalities. As proven by Bell, the local set L is a proper subset of Q [1] and those quantum points that are outside L are termed nonlocal. The most natural way to show that a given point p ∈ Q is not an element of L, is to use Bell inequalities. Recall their generic form to be where t = {t abxy } is a collection of some real coefficients t abxy and β L = max p∈L I is the local bound of the inequality (5). Analogously, the quantum or the Tsirelson bound of (5) is defined as β Q = sup p∈Q I. Clearly, violation of a Bell inequality by correlations p implies that it is nonlocal. Moreover, any point p violating maximally some Bell inequality (or, equivalently, achieving its Tsirelson bound) necessarily belongs to the boundary of the quantum set Q.
Correlation picture. Let us finally notice that it is often beneficial to describe correlations obtained in the Bell experiment by expectation values instead of probability distributions. A convenient way to do so in Bell scenarios involving d-outcome measurements is to consider the two-dimensional Fourier transform of the conditional probabilities p(a, b|x, y): where ω is the d-th root of unity ω = exp(2πi/d) and k, l = 0, . . . , d − 1. Importantly, for any quantum point p, the expectation values (6) can be represented as y } being collections of unitary operators with eigenvalues ω i (i = 0, . . . , d − 1) defined as the Fourier transforms of the corresponding projective measurements It is not difficult to see that A (k) x is simply the kth power of A x (and similarly for Bob's operators) and thus in what follows we simply write A k x . In this correlator picture the unitary operators A x and B y are d-outcome observables measured in the Bell set-up.
Self-testing.-We are now ready to present our main result: the self-testing statement for the two-qudit maximally entangled state and certain d-outcome measurements usually referred to as Collins-Gisin-Linden-Massar-Popescu (CGLMP) measurements [23][24][25] (see Appendix B for their explicit form).
To recall the task of self-testing, or more generally, DI certification, let us consider a quantum device performing a Bell experiment on some state |ψ AB ∈ H A ⊗ H B with some quantum d-outcome observables A x and B y , where the dimensions of the H A and H B are unknown but finite. The only information accessible to the user about how this device functions are the observed correlations p. The aim of DI certification is to reveal the form of the state and observables from violation of some Bell inequality by p (given that the observed correlations violate some Bell inequality). However, there are two degrees of freedom which can never be detected from the observed statistics. One is the set of local unitaries U A , U B that act on H A , H B , that is, the state U A ⊗ U B |ψ AB together with {U A A x U † A }, {U B B y U † B } will generate the same p. Another one is the presence of an auxiliary system on which the measurements act trivially. A particular instance of such DI certification, termed selftesting, infers a unique state and two sets of unique measurements up to these equivalences.
Clearly, a necessary condition to derive an exact selftesting statement for a state and measurements is that the obtained correlations p lie on the boundary of Q, and thus violate some Bell inequality maximally. Consequently, in order to derive a self-testing statement for |φ + d for any d the first task is to identify a class of Bell inequalities in the bipartite scenario with the minimum number of measurements for which β Q is achieved by it for any d. The only known Bell inequality meeting these requirements is the Salavrakos-Augusiak-Tura-Wittek-Acín-Pironio (SATWAP) Bell inequality [20], which in the correlator picture reads where a k = (1 + i)ω k/4 /2 and the classical value is given by β L = (1/2)[3 cot(π/4d) − cot(3π/4d)] − 2.
As proven in Ref. [20], the maximal quantum value of this inequality is β Q = 2(d − 1) and it is achieved by the maximally entangled state of two qudits |φ + d and the aforementioned optimal CGLMP measurements (cf. Appendix B). In what follows we show that this is in fact the only quantum system, up to additional of freedom and local unitary operations, realizing the maximal quantum violation of the inequality (8).
Let us introduce Z d = diag[1, ω, . . . , ω d−1 ] to be the ddimensional generalization of the σ z -Pauli matrix in the standard basis and T d to be the following matrix where δ i,j is the Kronecker delta: δ i,j = 1 for i = j and δ i,j = 0 otherwise.  (1) and two pairs of d-outcome measurements (10)- (11). and where a 1 = (1 + i)ω 1/4 /2, Z d and T d are defined above, and 1 A , 1 B are identity matrices acting on H A and H B .
Finally, the state |ψ AB is equivalent to |φ + d in the following sense where |aux A B is some auxiliary state from H A ⊗ H B .
Before proving this theorem a few comments are in order. First, as proven in Appendix B, the observables Z d and T d are unitarily equivalent to the CGLMP measurements. Second, our self-testing statement exploits a single Bell inequality involving the minimal number of two measurements per observer for which one can observe nonlocality in quantum systems. The number of local observables one needs to measure in order to self-test a state is of key importance from the application point of view. Finally, in the case of d = 2 N , the theorem implies a parallel self-testing protocol for |φ + 2 ⊗N , i.e., N copies of maximally two-qubit entangled state (see Refs. [21,22] for the previous approaches).
Proof. The proof is quite technical and thus is deferred to Appendix A. Here we only sketch its main ingredients.
Let a state |ψ AB and observables A x and B y maximally violate the SATWAP Bell inequality. Without loss of generality, we can assume that the reduced states of |ψ AB are full-rank. Exploiting the sum-of-squares decompositions of the corresponding Bell operator we can establish that Alice's and Bob's observables satisfy for x, y = 1, 2 and any n < d which is a divisor of the number of outcomes d. These conditions imply that the multiplicities of the eigenvalues of all the observables A x and B y are equal, which has two consequences. First, their matrix dimensions are a multiple of d, meaning that they act on, respectively,  (10) and (11). Finally, to obtain an analogous statement for the state we again employ the sum-of-squares decomposition, which in virtue of the relations (10) and (11), imposes the following conditions on the state |ψ AB : where we have denoted |ψ AB = U A ⊗ U B |ψ AB . These conditions can easily be solved leading to Eq. (12), which completes the proof.
Applications to randomness expansion.-An interesting application of our self-testing statement is DI randomness certification [5]. Let us consider the distrustful scenario where one of the observers performing the Bell test, say Bob, wants to generate random bits from the outcomes of his measurements, and Eve who supplies the measurement devices wants to access that outcomes.
We show with the aid of the above self-testing result that the maximal violation of the SATWAP Bell inequality certifies log 2 d bits of perfect randomness in the outcomes of Alice's and Bob's measurements. Without loss of generality we can focus on Bob's measurements. One quantifies the randomness obtained in a Bell experiment as − log 2 G(x, p), where G(x, p) is the local guessing probability defined through the following formula (15) where |ψ ABE is a three-partite state shared by Alice, Bob and Eve, {Q y } is a projective measurement performed by Bob, and E = {E (c) } is a d-outcome measurement on H E whose result is Eve's best guess of Bob's outcome. Finally, S p is the set of all possible strategies that Eve implements to guess the outcome of Bob's measurements, consisting of the state |ψ ABE and the measurement E, which reproduce the correlations p observed by Alice and Bob, i.e., Assume now that the inequality (8) is maximally violated by p. This means that, up to local unitary op- are eigenprojectors of the observables in Eqs. (10) and (11), respectively. Taking all this into account along with the fact that Bob's observables are traceless, one finally finds that G(y, p) = 1/d. Thus, − log 2 G(x, p) = log 2 d bits of randomness can be certified using our self-testing statement. This is the maximum randomness that could be extracted from a system of local dimension d, whereas it requires one bit of randomness to encode the inputs. This gives rise to unbounded randomness expansion as d can be arbitrary large.
Conclusion and discussion.-We propose the first, to the best of our knowledge self-testing statement for quantum system of arbitrary local dimension that unlike the previous approaches [11,17] does not rely on self-testing results for qubit states and exploits a truly d-outcome Bell inequality. Moreover, our self-testing result exploits only two measurements per party, which is the minimal number allowing the parties to observe nonlocal correlations. This makes our results interesting from the experimental point of view. In fact, violation of the SATWAP Bell inequality has already been experimentally tested in Ref. [26].
Several follow-up questions arise from our work. First, it would be interesting to study robustness of our self-testing statement and compare it to that of the scheme of [11,17]. Another interesting direction is to explore whether the d-outcome SATWAP Bell inequality can be modified to be maximally violated by partially entangled state and whether the resulting modifications can be used to make self-testing statements for those states, again with the minimal number of measurements per observer. Finally, it is also interesting to investigate whether, analogously to Refs. [27,28], our self-testing statement can be used to establish certification of the optimal amount of 2 log 2 d bits of local or global randomness from the maximally entangled state of two qudits with the aid of a local POVM or two projective measurements, respectively. From a more general perspective, our self-testing result might provide non-trivial insights into the structure of the set of quantum correla-tions along the lines of [29], or might find applications in delegated quantum computation [30].
Note added.-While working on this project we became aware of the work [31] exploring parallel selftesting with minimal number of measurements.
Acknowledgement.-We would like to thank A. Acín, P. Nurowski and A. Sawicki for useful discussions. This work is supported by Foundation for Polish Science through the First Team project (No First TEAM/2017-4/31).
JK acknowledges support from the Homing programme of the Foundation for Polish Science (grant agreement no. POIR.04.04.00-00-5F4F/18-00) cofinanced by the European Union under the European Regional Development Fund.

Appendix A: Proof of self-testing
In this appendix we provide the full proof of theorem stated in the main text. However, before diving into more detailed considerations we recall certain facts that will be used throughout the proof. First, let us recall that the explicit form of the SATWAP Bell inequality with two measurements per observer is where and the classical bound β d C is given explicitly by Then, the maximal quantum value of I d turns out to be β d Q = 2(d − 1). Let us also recall that the sum-of-squares decomposition of the corresponding Bell operator B d is with Let us notice here that a d−k = a * k and therefore C i ] † for any k = 1, . . . , d − 1 and i = 1, 2. Introducing finally the following matrices we can now state and prove our main theorem. and where a 1 is defined in Eq. (A2), whereas Z d and T d are defined in Eq. (A7); 2. The state is |ψ equivalent to the maximally entangled state of two qudits |φ + d in the sense that where 1 A and 1 B are identity operators acting on H A and H B , respectively, and |ψ is some state from H A ⊗ H B .
Proof. The proof consists of three major steps to prove the theorem. First, we establish that the observables measured by Bob act on a Hilbert space of dimension multiple of d and simultaneously obtain their explicit form. Second, exploiting the symmetry of the Bell-operator we recover the desired form of A i and subsequently characterize the shared entangled state. Observables. We begin our proof by showing that with the aid of the sum-of-squares decomposition (A5) one can prove that the following sets of conditions hold true: To prove (A11) we first notice that for the state |ψ and observables A i and B i the sum-of-squares decomposition (A5) implies that are satisfied for i = 1, 2 and k = 1, . . . , d − 1. Exploiting the fact that A i are unitary we can rewrite the above as We then make use of (A13) for k = 1 to see that A † i |ψ , which allows us to rewrite the above relation as which is equivalent to the relation (A11) on the support of the second subsystem of |ψ .
To prove (A12) we again consider relation (A13) and apply to it. This due to the fact that A i is unitary gives which is equivalent to (A12) on the support of the second subsystem of |ψ .
In what follows we show that the conditions (A11) and (A12) are enough to fully characterize Bob's observables. First, in Lemma 1 we show that these conditions imply that for any positive integer n which is a proper divisor of d.
where H B is some auxiliary Hilbert space of unknown but finite dimension. Second, there exists a unitary matrix V B : H B → H B such that where Z d is defined in Eq. (A7) and 1 B is the identity acting on H B . Thirdly and finally, Lemma 3 implies that there exists another unitary operation V B that preserves B 1 = Z d ⊗ 1 B and Alice's observables. Let us now concentrate on Alice's observables. To determine them we will exploit the symmetries of the SATWAP Bell inequality. Precisely, let with the aid of the facts that a d−k = a * k and that A d i = B d i = 1, the SATWAP Bell expression I d can also be written as Moreover, one can construct another sum-of-squares decomposition for the SATWAP Bell inequality in terms of these combinations of Alice's observables. Precisely, one has where Thus, using the same reasoning as in the case of Bob's observables one can show that the C (k) i operators must satisfy the same relations as the C (k) i operators, that is, as well as with i = 1, 2 and k = 1, . . . , d − 1. The crucial step now is to realize that C with H A being some auxiliary Hilbert space corresponding to the party A of unknown but finite dimension. Moreover, there exists a unitary operation V A : H A → H A such that Thus, similarly to Bob's observables we can fully characterize Alice's measurements. In principle we could end our characterization of observables already here, however, for further benefits we will apply another local unitary operation W A : C d → C d to the qudit part of A i observables so that and The existence of W A is shown under i operators allows us to bring them to In particular, if Alice's observables are transformed by applying the unitary U A , then C (1) 1 = Z * d and C (1) 2 = T * d . As we will see next this choice of the A i observables as well as the C (1) i operators guarantees that the self-tested state (up to some additional degrees of freedom) is the maximally entangled state of two qudits written in the standard basis State. Let us finally prove Eq. (A10). From the fact that |ψ violates maximally the SATWAP Bell inequality we know, by virtue of the SOS decomposition (A23) that with i = 1, 2 and k = 1, . . . , d − 1, which after application of U A ⊗ U B can be rewritten as where by 1 A B we denote the identity acting on H A ⊗ H B .
Since the state |ψ belongs to where |ψ ij ∈ H A ⊗ H B are some in general nonorthogonal and not normalized states. With this decomposition the condition (A34) with k = 1 implies that This condition directly implies that for any pair i, j = 0, . . . , 0 such that i = j.
Let us now exploit the other condition (A35) with k = 1, which after projecting onto the ket ii| gives where we have used the fact i|T * d |k = k|T * d |i as T d is symmetric. Using the explicit form of T d given in (A7) we can rewrite the above as We now move the |ψ ii ket to the left-hand side of the equality (for d = 4 there is nothing to move) and obtain the following relation, from which it follows that all |ψ ii are equal. This together with (A38) allows us to conclude that We thus recover (A10), completing the proof.
Below we append three lemmas used in the proof of the above theorem. Remark that, series of Observations are employed inside the proof of Lemma 1 and Lemma 3.

Lemma 1.
Consider two unitary observables B i (i ∈ {1, 2}) acting on a finite-dimensional Hilbert space whose eigenvalues are ω l (l ∈ {0, . . . , d − 1}). If they satisfy the conditions (A11) and (A12), then for any proper divisor n of d, that is, a positive integer such that n < d and d/n ∈ N, Proof. Plugging the explicit form of C (k) 2 given in Eq. (A6) into Eq. (A12) and using the fact that Using then the fact that a * k = a d−k and a * k /a k = iω k 2 , a simple calculation leads us to the following condition Due to the fact that B d i = 1, the above relation (A45) extends to any integer k ∈ Z. We now establish some relations under observations 1.1 and 1.2 which will be used later in the proof.
Observation 1.1. First, we show that the following identities hold true for any non-negative integers s, x, y ∈ N ∪ {0} and We present only the proof of the first identity as that of the second one is analogous. It is pretty straightforward and consists of multiple application of the identity (A45). Let us thus focus on the right-hand side of Eq. (A46) and consider (A45) for k = 2sx, multiply it by B x 1 and trace both sides of the resulting equation. This gives We consider again (A45) for k = (2s − 1)x, multiply it by B x 2 and then trace both sides, which results in After plugging Eq. (A49) into Eq. (A48) we obtain We have thus lowered the power of B 1 from (2s + 1)x in Eq. (A48) to (2s − 1)x in Eq. (A50). We repeat this double substitution until the power of B 1 is x, each time acquiring the phase ω −x (notice that in this case the power of B 2 goes from −2sx in Eq. (A48) to zero). As a result we arrive at which is what we wanted to show.

Observation 1.2. Next, the following simple relation between traces of powers of B i operators holds true
In order to prove this observation let us consider (A11) for k = 2x, which we can equivalently state as This after taking into account (A6) implies We then multiply the above by B −x 1 and by taking trace on both sides we find By multiplying it by B x 2 and tracing both sides we obtain Further, substituting the expression of Tr(B −x 1 B 2x 2 ) from Eq. (A56) into Eq. (A55) leads us to (A52). Let us now make use of the above observations for finally prove the validity of (A43). We consider the cases of even and odd d separately.
Even d. Let n be a divisor of d, i.e., d/n ∈ N. There are two possibilities, d/n is even or odd, i.e., there exists integer s such that n = d/(2s) or n = d/(2s + 1), respectively. Whenever d/n is even, we consider Eq. (A46) for x = n = d/(2s), which reads (A57) Using then the facts that B d i = 1 and ω d/2 = −1, the above relation simplifies to, Thus, for any n such that d/n is even, Tr(B n 1 ) = 0. Similarly, by taking y = d/(2s) in Eq. (A47) one can conclude the same for B 2 , i.e., Tr(B n 2 ) = 0. Now, for any divisor n of d for which d/n is odd, we choose x = n = d/(2s + 1) in Eq. (A46), which leads us to Tr(B n 1 ) = (−1)ω −n/2 Tr (B n 2 ) .
The above relation together with Eq. (A52) imply Tr(B n i ) = 0 for any n such that d/n is odd and n ≤ d/2. Thus, we have shown that for any n which is a divisor of d and n = d, Tr(B n i ) = 0 for i = 1, 2. Odd d. Since d is odd, there exists an integer s such that 2s + 1 = d. Plugging s = (d − 1)/2 in Eq. (A46) we get which due to the facts that B d 2 = 1 and ω d/2 = −1 simplifies to Now, for any odd x such that x ≤ d/2 we obtain from Eqs. (A52) and (A61) that Tr(B x 1 ) = Tr(B x 2 ) = 0. Since all the divisors of an odd d (except d) are odd and less than d/2 , we conclude (A43).

Lemma 2. Consider a real polynomial
with rational coefficients λ i ∈ Q. Assume that ω n with ω = e 2πi/d is a root of W(x) for any n being a proper divisor of d, i.e., n = d such that d/n ∈ N. Then, λ 0 = λ 1 = . . . = λ d−1 .
To be able to prove this lemma we need to introduce the following two relevant concepts.

Euclidean division algorithm. Let F(x), G(x) ∈ F[x]
be two nonzero polynomials over the field F such that deg G(x) ≤ deg F(x). Then there exist two polynomials Q(x), Cyclotomic polynomials. Let n be any positive integer. The nth cyclotomic polynomial φ n (x) is defined as the unique irreducible polynomial over the field of rational numbers whose root is the primitive root of unity exp(2πi/n), and is given explicitly by the following formula Let us then consider a positive integer d and denote by n 1 < n 2 < . . . < n k its proper divisors with n 1 = 1 and n k < d. Then, it is known that the product of the cyclotomic polynomials φ d/n i over all these proper divisors can be represented by the following simple formula Let us finally notice for further purposes that ω n i for some divisor n i is a root of φ d/n i (x) but not of φ d/n j (x) for any other divisor n j = n i of d.
Proof. Let n 1 < n 2 < n 3 < . . . < n k < d with n 1 = 1 and n k < d be the proper divisors of d. From the assumption we know that ω n i is a root of W(x), that is We now use this last condition along with the Euclidean division algorithm to prove that the polynomial W(x) must be a product of cyclotomic polynomials φ d/n i (x). We do it in a recursive way Let us begin with n 1 = 1. The Euclidean division algorithm implies that . From Eq. (A65) we know that W(ω) = 0, implying that R 1 (ω) = 0. As φ d (x) is the unique polynomial of minimal degree over the rational field whose root is ω, the latter implies that R 1 (x) ≡ 0. Thus, Let us then move on to the next divisor of d, n 2 > n 1 . Due to the fact that ω n 2 is not a root of φ d (x), the condition (A65) implies that Q 1 (ω n 2 ) = 0. We therefore apply the Euclidean division algorithm to Q 1 (x), that is, . Due to the fact that ω n 2 is a root of Q 1 we see that R 2 (ω n 2 ) = 0. Again, due to the fact that φ d/n 2 (x) is the unique polynomial of minimal degree over the rational field whose root is ω n 2 , the latter implies R 2 (x) ≡ 0. Consequently, Applying the above procedure iteratively for every divisor n i of d it is now not difficult to see that our initial polynomial W(x) is of the following form with Q(x) is a rational polynomial. Owing to the fundamental relation satisfied by the cyclotomic polynomials given in Eq. (A64) one sees that the degree of the above product of cyclotomic polynomials is d − 1 and equals the degree of W(x), which directly implies that Q(x) is simply a constant factor, denoted Q. Moreover, by comparing Eqs. (A62) and (A70) and taking into account Eq. (A64), we immediately see that λ i = Q for i = 0, . . . , d − 1, meaning that λ 0 = λ 1 = . . . = λ d−1 , which completes the proof. Proof. We begin by proving the following relation for B 1 and B 2 matrices: We prove this relation by induction. For k = 1 it is not difficult to see that both its sides equal B 2 . Assuming then that (A71) holds true we will prove it for k + 1. To this end, let us rewrite (A11) as Eq. (A6)] can be rewritten as Having established (A71), let us now write B 2 as where F ij are some matrices acting on H B . Our aim now is to determine F ij using the conditions (A71). We first focus on F ii and then move on to F ij with i = j. Finding F ii . Taking the relation (A71) for k = d − 1, we obtain which after taking into account that B 1 = Z d ⊗ 1 B and that B 2 is given by Eq. (A75), can be rewritten as By projecting the first subsystem onto |i i|, we arrive at the following condition for F ii matrices Conjugating the above equation on both sides we obtain After plugging into the above relation F † ii from Eq. (A78) we obtain an equation for F ii whose solution is Finding F ij . Let us now move on to determining the F ij matrices for i = j. We formulate our derivation as a sequence of observations. First, by comparing the matrices appearing on |i j| position of both sides of (A77) with i = j, we obtain the following equation which after taking into account the fact that for i = j, reduces to Thus, (A76) only provides a relation between the symmetric elements of B 2 in the form (A83). To find the explicit form of F ij , we have to look for equations involving higher order terms in F ij . To this aim, let us prove the following observation.
Observation 3.1. The following conditions hold true Again using (A71) for B k 2 , the left-hand side of the above equation can be written as Let us now focus on the sum appearing in Eq. (A86). After substituting the explicit forms of B 1 and B 2 into it, we can rewrite it as After splitting the above summations into the cases of i = l and i = l, we obtain Let us then consider the sum appearing in Eq. (A85). After plugging the explicit forms of B 1 and B 2 into it we arrive at We finally substitute B k+1 2 from Eq. (A85) into Eq. (A86), and then plug the explicit forms of B 1 and B 2 as well as the expressions (A88) and (A89) into the resulting equations, which leads us to (A84). This completes the proof of observation 3.1. Now, having the formula (A84), we proceed as before. That is, we first look at the diagonal elements of it and solve the resulting equations. Then, we consider the off-diagonal elements of Eq. (A84) to have full characterization of F ij . As for the diagonal part of Eq. (A84) we establish the following observation.
To prove the above relation let us consider the |i i| elements of (A84) which, after some algebra, gives us the following relation which after substituting F ii from Eq. (A80) into it simplifies to This after some manipulations can further be rewritten as and after taking into account (A83) and changing the index l to j, we obtain the desired relation (A90).
The subsequent observation provides the solution of equations (A90).

Observation 3.3.
The solution of (A90) is given by We multiply Eq. (A90) by ω kn with k = 0, . . . , d − 1 and n = 1, . . . , d − 1 and then sum the resulting identity over all k s, which leads us to where we have exploited the fact that for any n = 1, . . . , d − 1, Applying then Eq. (A114) to the right-hand side of Eq. (A95) and the fact that to its left-hand side we can bring it to To complete the proof it suffices to realize that for any i = 0, . . . , d − 1 there exist n = 1, . . . , d − 1 such that i − n mod d is any number from {0, . . . , d − 1} different than i.
The relation (A94) is unfortunately not sufficient to fully characterize F ij . However, we can find a unitary operation U that preserves B 1 and allows us to obtain an explicit form of a few matrices F ij from Eq. (A94).
To be more precise, let us consider a unitary matrix U acting on C d ⊗ H B of the form where We can readily check that UB 1 U † = B 1 . Let us then denote Note that, all the algebraic relations for F ij obtained so far hold also forF ij , andF ii = F ii . Now, we see thatF where to obtain the last equality we employed Eq. (A94). Using Eq. (A83) we then obtainF j0 . Thus, the remaining part of the proof of Lemma 3 is to obtain the elements F ij such that i, j = 0 and i = j. To do so, we use the off-diagonal elements of (A90) to obtain a set of relations as follows. Taking the inner product with i| . |j (where i = j) on the both side of (A84) we obtain After rearranging some terms and using the formula for F ii given above, we express the above equation as Next, we set i = 0 and obtain and then substituteF 0j from (A102) which gives Let us then take the term corresponding to l = j out of the sum and use Eq. (A80) to express F jj . This after some simplifications allows us to rewrite Eq. (A107) as (A108) Finally, changing the index l to i leads us to (A103). Note that the equations derived for F ij are also valid forF ij , and F ii = F ii .
The following observation provides the solution of (A103).
Observation 3.5. The solution of the equation (A103) is We multiply Eq. (A103) by ω −kn with n ∈ {1, . . . , d − 1} such that n = j and then sum both sides of the resulting formula over k = 0, . . . , d − 1, obtaining We now notice that the first sum on the left-hand side of the above and the last two sums on the right-hand side simply vanish due to Eq. (A96) and the fact that n = j. Exploiting then Eq. (A115) as well as the fact that  Proof. To prove the first relation (A114) we first notice that its left-hand side can be rewritten as Observing then that the expression appearing on the right-hand side under the sum is actually a sum of a geometric series, one can rewrite the above equation as Let us finally notice that for any n = 1, . . . , d − 1 one has the following chain of equalities After plugging this last formula into Eq. (A118) we arrive at Eq. (A114).
To show the second identity we use the following relation Taking the derivative and multiplying by x on both sides, we get Substituting x = ω n and using the fact that ω d = 1 we obtain (A115).