Entanglement Trajectory and its Boundary

In this article, we present a novel approach to investigating entanglement in the context of quantum computing. Our methodology involves analyzing reduced density matrices at different stages of a quantum algorithm's execution and representing the dominant eigenvalue and von Neumann entropy on a graph, creating an"entanglement trajectory."To establish the trajectory's boundaries, we employ random matrix theory. Through the examination of examples such as quantum adiabatic computation, the Grover algorithm, and the Shor algorithm, we demonstrate that the entanglement trajectory remains within the established boundaries, exhibiting unique characteristics for each example. Moreover, we show that these boundaries and features can be extended to trajectories defined by alternative entropy measures. The entanglement trajectory serves as an invariant property of a quantum system, maintaining consistency across varying situations and definitions of entanglement. Numerical simulations accompanying this research are available via open access.


Introduction
An appropriate dose of entanglement can be seen as the resource of the quantum advantage for computation [1,2].It can not be too low since quantum states whose von Neumann entropy only scales logarithmically with the number of qubits can be described in terms of Matrix Product States [3].This means the amount of entropy in the system along the quantum algorithm should be large; otherwise, an efficient algorithm would exist to solve the same problem using Tensor Networks.It can not be too high either; otherwise, quantum resources can be replaced with a fair coin since maximum entanglement usually implies maximum randomness [4].
Therefore, all quantum states that present a potential for computational speedup lie in the separation between "not entangled enough" and "too entangled," which is very narrow since almost all quantum states are extremely entangled [4].This article introduces a new approach to studying this gap, referred to as the entanglement trajectory.Since genuine multipartite entanglement measures are not available for states with a large number of qubits [5], for simplicity, this article only focuses on bipartite pure quantum states.The result is already fruitful under these constraints.For such a quantum system, the value of the von Neumann entropy originates from the eigenvalues λ i of the reduced density matrix ρ A , where the subsystem A is partially traced over its complementary system.In our experience, these eigenvalues usually consist of an isolated dominant member λ 0 , while the rest are numerically approached to 0. This particular behaviour provides the initial motivation for this paper.We ask, what happens when we plot the entanglement trajectory for a quantum process performing a computation, with the dominant eigenvalue on the xaxis and the von Neumann entropy on the y-axis?How does the trajectory change when using other measures of entanglement on the y-axis?What additional information can be gained about the system?

Examples from quantum computation
We briefly summarize what we expect to see before studying different quantum systems.We use the von Neumann entropy defined on the natural logarithm E = − λ i ln λ i in this article, with λ i ≥ 0 and λ i = 1.When the size of the subsystem A is α and the size of its complementary system is β ≥ α, for the dominant eigenvalue, we should have 1 ≥ λ 0 ≥ 1 α and the maximum entropy is ln α.Also, when only two non-zero eigenvalues exist, the entropy reaches its lower bound This f 1 will be plotted in a gray dash-dotted line throughout this article.
When λ 0 = 1 w and w is an integer with α ≥ w ≥ 2, the entropy reach its lower bound f 2 when there is w non-zero eigenvalues and all of them are 1 w , with which will be plotted in a gray dashed line throughout this article.Notice that f 1 and f 2 do not depend on the size of the subsystem α and will be plotted in every figure in this article when the y-axis is the von Neumann entropy.Furthermore, the entropy will reach its upper bound f 3 when the rest of the eigenvalues are all equal and sum up to 1 − λ 0 , when α ∼ α − 1, we have which will be plotted in a gray solid line in every figure with one α and with the von Neumann entropy as the y-axis through this article.
The comparison between the boundary formed by f 1 , f 2 and f 3 and the actual numerical result is shown in FIG. 1.As we can see, the numerical boundary can be approximated with these three analytical lines while staying strictly within the area confined by them.At this point, the boundary is a direct result of the definition of von Neumann entropy.This boundary is tight and does not rely on the purity of quantum states.It is impossible to provide a counter-example.
In this section, we present three quantum algorithms and track their entanglement trajectories.They are all well-studied and with explicit numerical implementation on the open-access platform Qibo [6,7].Using them as examples is both convincing and convenient.The code which generates every figure in this article is on GitHub [8].We must highlight that in each figure, the lines that connected data points do not have any physical meaning and serve only the purpose of visualizing the sequential order between quantum states.The transitions between quantum states are not continuous and can never be presented as a path connecting two data points.Otherwise, it will be classical.For the result, we find out these trajectories present a unique feature for each quantum algorithm and can be used as "fingerprints."However, since the eigenvalues of these density matrices do not follow any analytical distribution, it will be difficult to guarantee a further mathematical explanation.

Adiabatic quantum computation
In Ref. [2], the authors have studied the entanglement properties of the adiabatic quantum algorithm for solving the NP-complete Exact Cover (EC) problem, which is a particular case of 3-SAT problem.An Exact Cover instance is characterized by a set of clauses containing 3 bits that are considered satisfied if one is in position 1 while the other remains at 0. The solution of this instance is one bit-string that fulfills all the clauses at the same time.To explain the adiabatic quantum algorithm in a nutshell [9], the quantum system is initialized with an "easy" Hamiltonian H 0 , which is a magnetic field in the x direction, where the ground state is an equal superposition of all possible computational states.The system is then adiabatically evolved to a "complicated" Hamiltonian H p that encodes an instance of EC problem, where the ground state is the solution.During the adiabatic evolution, the system can be described using a Hamiltonian H s , an interpolation between H 0 and H p , In the original paper, the authors tracked the von Neumann entropy of different equal bi-partition (using half of the qubits as the subsystem and α = 2 n/2 ) of the ground state of H s during the adiabatic evolution for solving EC problem.They found out that the entanglement always reaches its pic at s ∼ 0.7 while the pic scales O (n) for different sizes of systems.Their results imply that this pic corresponds to an undergo phase transition; therefore, a classical computer cannot efficiently simulate the adiabatic evolution.However, the difficulty of simulation does not signify a quantum advantage.Though it is almost impossible to prove, as pointed out by Ref. [10], the gap between the ground energy and the first excited energy of H s is believed to be exponentially small for an NP-complete problem [11].Evaluation requires an exponentially long time to reach the correct solution; even an optimal adiabatic scheduling function for s can only produce a quadratic speedup [9].The numerical simulation has been re-implemented in Qibo, and the code is on GitHub [12].In this section, we revisit their study and track the entanglement trajectory over random bi-partition of the ground state with the interval of s as 0.1.The results are shown in FIG. 2. As we can see, the data points collected are distinct for different instances and bi-partition, but they almost stay on the same path.

Grover algorithm
The quantum search algorithm for an unordered database was first introduced by Grover [13].The search space is initialized as the equal superposition of all possible solutions (usually the equal superposition of all computation basis).A Grover step consists of two parts: the oracle that encodes the description of the problem and inverts the sign of the solutions, followed by the diffusion operator that amplifies the probability of measuring the solutions.The repetition of Grover steps brings the probability of finding the right solutions close to 1. Then in Ref. [2], the authors have proved that the entropy in the register along the Grover algorithm is bounded by 1 after each Grover step.Consequently, the quantum advantage of the Grover algorithm should be expected from entanglement present in the middle of the oracle.In this section, we compare the entanglement trajectory of two search algorithms.The bi-partition splits between the search space (where the diffusion operator takes place) and its complementary space.
The first example is to use the Grover algorithm to solve the EC problem mentioned in the section 2.1.The oracle has the same number of ancillas as clauses to ensure that all clauses are satisfied.For each clause, a CNOT gate from each qubit to its ancilla, followed by a multi-Toffoli gate controlled by the three qubits targeting the ancilla, activates the ancilla only if the clause is satisfied.The circuit is shown in FIG. 3, and the implementation in Qibo is provided on GitHub [14].Due to its high entanglement, this quantum circuit can not be efficiently simulated on classical hardware.But it does not imply a quantum advantage.The quantum algorithm scales O 2 0.5n [15] while the best classical algorithm scales O 2 0.39n [16].The second example [17] is the application of the Grover algorithm to find the preimages of a toy Sponge Hash function [18].Their circuit consists of 18 qubits, where n = 8 are initialized as a uniform superposition of all states, referred to as the search space of Grover's algorithm.Like the previous example, this Grover circuit can not be efficiently simulated with a classical computer due to the high entanglement in the middle of the oracle, as indicated in their paper.However, the quantum algorithm does provide an advantage.For a given function, if there exists a classical algorithm that can output the pre-images of an arbitrary cipher-text more efficiently than blind search, it can not be qualified as a Hash function.And Grover algorithm performs blind search with a quadratic speed-up.Details are provided in their paper and their GitHub page [19].In this section, besides the entanglement, we have also tracked the entanglement trajectory.A simplified circuit is shown in FIG. 4.
The comparison of the entanglement trajectory of these two Grover algorithms can be found in FIG. 5.The sizes of subsystem α = 2 8 and the complementary system β = 2 10 are the same.As we can see, although the density of data points is distinct, two Grover algorithms that solve different problems share a similar path.In particular, the trajectory of the Grover algorithm that solves the Hash function is slightly longer and slightly above the one that solves the EC problem.

Semi-classical Shor algorithm
The original paper by Shor [20] stated that the factorization problem could be reduced to order finding.Then, Shor proposed a quantum order-finding algorithm that runs in polynomial time.To factorize a number N = p × q, a number 1 < a < N is randomly  chosen and the quantum algorithm finds a number r such that a r = 1 mod N .If r is neither odd nor verifying a r/2 = −1 mod N , p and q can be recovered from r.
Since the efficiency of the Shor algorithm is due to the efficiency of the quantum Fourier transform (QFT), previous research mainly focused on studying the entanglement entropy before the application of the QFT and in between the register of the QFT and the quantum modular exponentiation (QME) [2,21].It is concluded that with this partition, the reduced density matrix is diagonal with m non-zero elements, which all equals 1 m with positive integer m.This behaviour matches the scenario that f 2 is deduced; therefore, its entanglement trajectory (one data point for a given N and a) should locate on f 2 .And we verify this property with numerical simulation using Qibo, which is available on GitHub [22].They also showed that m is O (r), which indicates that the entropy scales O (ln r).
In this section, we investigate the entanglement trajectory of the Shor algorithm with the most up-to-date version for minimizing the number of qubits provided by Beauregard [23,24,25], which only needs 2n+3 qubits.The QFT is performed semi-classically with one controlling qubit, as shown in FIG. 6.The numerical simulation with Qibo is also available on GitHub [22].Since the space of QFT has been reduced to one qubit, the entanglement between this qubit and the rest of the system is upper bounded by 1, and it would be more convincing to study this circuit with another partition.We use the bi-partition between the register that stores the output of the QME (n qubits) and the rest of the qubits to track the entanglement trajectory and sample ρ A in the middle of each QME, as shown in FIG. 7. The size of the subsystem is α = 2 n , and the size of the complementary system is β = 2 n+3 .The result is shown in FIG. 8.The entanglement trajectories concentrate on a particular path connecting λ 0 = 1 2 with an interval of 1 2 ln 2 and despite the difference between system sizes and inputs.The maximum entanglement is located at λ 0 = 1 2 .We do not have an explicit interpretation of the 1 2 ln 2 interval nor f shor .As we can see from the numerical output, some data points are located quite far from this approximation.The interval is derived from the fact that the eigenvalues of ρ A after the controlled SWAP operator frequently consist of λ 0 = 1 2 and 2 m non-zero elements of 1  2 m , with positive integer m.And the formula f shor is derived from the fact that the eigenvalues of ρ A before the controlled SWAP operator frequently consist of λ 0 = m+1 2m and m − 1 non-zero elements of 1 2m , with positive integer m.

Mathematics Background
In this section, we briefly introduce two key ingredients for the asymptotic analysis of entanglement trajectory.As most quantum states are random, it is natural that the entanglement of a physical system should concentrate on the entanglement obtained with  random density matrices, which can be generated by partially tracing a larger random state [26,27].It is proven that the law of a random density matrix is the law of a Wishart matrix conditioned by trace 1 [28,29].Therefore, the first ingredient we present is one theorem in random matrix theory defined on the Wishart ensemble, the Marchenko Pastur distribution.
The second ingredient consists of analyzing the isolated dominant eigenvalue when the coefficients in the reduced density matrix do not centre at 0 using a basic theorem in linear algebra.This method mainly focuses on the expectation of the dominant eigenvalues rather than their distribution.It is sensitive to the condition that we generate random reduced density matrices (by setting the mean and the variance of the random variables).Also, it has a more straightforward analytical expression that helps us derive the relation between the entropy and the dominant eigenvalue.

Definition 1. Wishart matrix
Let X be an α×β complex matrix whose elements are complex normal N C (0, 1) random variables that are independently identically distributed (i.i.d).The α × α matrix XX † is referred to as a Wishart randomly generated matrix with parameters α and β.All Wishart matrices form the Wishart ensemble [36].
Then the ESD of Y β , µ α converges weakly, in probability, to the distribution with density function µ given by and with  Ref. [30] contains the proof for the MPD using the moment approach.In short, the MPD predicts the distribution of the eigenvalues of a Wishart matrix by giving the analytical outline of its normalized histogram.Let X be an α × β random matrix with variables N (0, 1), λ = α/β.Then the normalized histogram of the eigenvalues of as shown in FIG.
9. In this article, the λ + will be indicated with a gray dotted line (can be horizontal or vertical) in the figure when specific α and β are given.
We can rewrite X into and where H is an α × β with all elements γ, and X an α × β i.i.d random matrix with mean 0 and variance σ 2 .And we also have the relation and since the elements of X are centred at 0. We define the spectral norm [38,39] of a matrix P as which satisfies the inequality Therefore, we have which means In the same way, implies So we have when Then HH † is an α × α matrix with all elements equal to β|γ| 2 , which means all its eigenvalues are αβ|γ| 2 .
Finally we have as shown in FIG. 10, while the remaining eigenvalues still follow the MPD [40].We keep the notation γ 2 instead of |γ| 2 for simplicity.

Asymptotic Analysis
In this section, we first demonstrate how to use the section 3.1 to recover a famous theorem [49].After that, we introduce the element from the section 3.2 to give a new asymptotic result for the entanglement trajectory.Then we compare our analysis with previous examples in the section 2.

Average entropy of a subsystem
Here we consider random density matrices obtained by partially tracing larger random pure states.There is a strong connection between these random density matrices and the Wishart ensemble of random matrix theory, that can be used to determine the formula of entropy by Don Page in 1993 [28,49,5,50,51,52,53].

Theorem 2. Average entropy of a subsystem
If a quantum system of Hilbert space dimension αβ is in a random pure state, the average von Neumann entropy We can use the MPD to reproof it.
Proof.Consider the random reduced density matrix under the form ρ A = ZZ † / Tr ZZ † [28], where Z is an α × β random matrix for 1 ≤ α ≤ β.We have E = − α−1 i=0 λ i ln λ i where λ i are eigenvalues of ρ A .
is plotted in gray dotted line.The minor difference lies between λ 0 and λ + when γ ∼ 0 is due to the effect of the Tracy-Widom distribution [41,42,43,44] and can be generalized with a phase transition [45,46,47,48].Since this phenomenon is highly non-analytical also numerically negligible for our study, we will not elaborate in detail.
We can eliminate the trace and rewrite ρ A into a form that the MPD is applicable: , where X is an α × β random matrix with mean 0 and Therefore, X √ β can be considered a "normalized" matrix, with coefficients c i that correspond to the quantum state |ψ⟩ = αβ−1 i=0 c i |i⟩ with αβ−1 i=0 |c i | 2 = 1.We define λ = α/β.The von Neumann entropy of a subsystem of dimension α can be written as with and Then we have which can be simplified back into with the integral in Eq. ( 24)

Average entropy of a subsystem (general case)
In this section, we consider a slightly different variant when the mean of X is not 0 but γ.We still have λ = α/β and the eigenvalues of ρ A are still λ 0 , ..., λ α−1 with λ 0 ≥ ... ≥ λ α−1 .
The standard deviation σ γ should be evaluated with γ to maintain the sum of all eigenvalues as 1.We have The remaining eigenvalues λ 1 , ..., λ α−1 will follow the MPD, and their sum would be with the integral in Eq. ( 27) To have the sum of all eigenvalues of ρ A equal 1, we should have the relation This relation can also be obtained using the assumption that the quantum state is normalized, or ρ A has the trace 1.These three conditions are equivalent, and all originated from the fact that only pure quantum states are investigated in this article.
We can deduce that then use this expression to replace σ 0 in Eq. ( 21) and Eq. ( 22).We then have the entropy depending on λ 0 with After applying the integral in Eq. ( 24) and massive reduction, we finally have which is the average entropy for a given λ 0 .With this expression, we can recover Eq. ( 22) for λ 0 → 0. When β → ∞, we can recover the upper bound In FIG.11 we compare our analytical solution with the numerical results.The result given by asymptotic analysis reduces the tight boundary to a more flexible boundary, as shown in FIG.12.In contrast to "tight," which implies a barrier forbidden by mathematics, "flexible" indicates a barrier protected by randomness.The absolutely maximally entangled (AME) states [54,55] are rare and located at the extremity of the tight boundary if they exist [56].And almost all quantum states (average entropy of a subsystem) are located at the extremity of the flexible boundary.

Relation with previous examples
We re-plot FIG. 2, FIG. 5 and FIG. 8 with their flexible boundary, as shown in FIG. 13, FIG. 14 and FIG. 15.As we can see, the entanglement trajectories obtained from physical systems stay below the line deduced from MPD, which is a non-trivial result.Although it is difficult to theoretically or physically construct quantum states that are located above the line of average entropy over every possible bi-partition, such as AME states, preparing such a state for one bi-partition is easy.A possible explanation might be that reaching an entanglement higher than average will not help computation since the asymptotic analysis is derived from randomness.However, this argument is difficult to prove because counterexamples are not hard to formulate.In particular, we have In the case of the Shor algorithm using 2n + 3 qubits.If we postulate that when the entanglement trajectory reaches the maximum as λ 0 = 1 2 , and the entropy will be near and upper bounded by E λ 0 = 1 2 , then the maximum entropy that we expect from this particular partition of semi-classical Shor circuit would be n+2 2 ln 2 − 1 32 , by replacing α with 2 n and β with 2 n+3 .That means the semi-classical Shor algorithm can not be simulated efficiently on a traditional computer.This example can be seen as an application of the asymptotic analysis.

k-almost Prime states and their unions
In the previous section 4.3, we have demonstrated with examples that our asymptotic analysis defines a flexible upper bound for the entanglement trajectory.In this section, we provide another example to justify why this upper bound should not be lower.In this case, a series of quantum states are tracked to obtain the entanglement trajectory instead of those taken from different moments along a physical process.
In Ref. [57,58,59], the authors have studied the entanglement of the Prime state.The Prime state is defined as an equally weighted superposition of prime numbers in n qubits computational basis where each prime number p = p 0 2 0 + p 1 2 1 + ... + p n−1 2 n−1 is implemented as |p⟩ = |p n−1 , ..., p 0 ⟩, and π (2 n − 1) is the amount of prime numbers less than 2 n .For example, when n = 4, we have The natural bi-partition entanglement entropy of the Prime state for even n is investigated in Ref. [58] (the logarithm in base 2 was used).In this example, we explore its other variant, referred to as the k-almost Prime state.k-almost Prime is a number that has k prime factors, usually denoted as p k .Here we define k-almost Prime state |P k ⟩ as the superposition where π k (2 n − 1) is the amount of k-almost Prime numbers less than 2 n .And their union is defined as the equal superposition of k i=1 p k and are k-almost Prime state.Also are their union.
In particular, we notice the fact that of which the entanglement entropy is 0. On the other hand, |U n−1 ⟩ would be the equal superposition of all computational basis except |0⟩ and |1⟩, which also has 0 entanglement entropy if we take the natural bi-partition.As a consequence, before we draw the entanglement trajectories of |P k ⟩ and |U k ⟩ for k ∈ {1, ..., n − 1}, we can already anticipate that they will join on the point of the Prime state as k = 1, also on the point (1, 0) for k = n−1.
The actual numerical simulation is shown in FIG.16 and supports our prediction.The entanglement trajectory is extremely close to our asymptotic analysis while staying strictly below it.

General Entanglement Trajectory
This article mainly focuses on the entanglement trajectory defined with the von Neumann entropy since it is the most widely used entropy in other research.However, different variants of entropy are only different functions defined on the eigenvalues of ρ A .Their asymptotic analysis can be obtained with the same method in section 4.2, which allows us to define boundaries of different entanglement trajectories.Here we provide three examples.

Entanglement trajectory in the momentum space
If a quantum state |ψ⟩ is given in position state, its Fourier transform QF T |ψ⟩ gives the same state in the momentum space.Considering a random state |ψ⟩ = αβ−1 i=0 c i |i⟩.Assuming that the dominant eigenvalues of its reduced density matrix ρ A is still λ 0 , then the normalized state after performing the QFT on |ψ⟩ can be written as where |Q 0 | 2 ∼ λ 0 and Q i are random variables with mean 0 and variance σ 2 QF T = (1 − λ 0 ) /αβ, which according to the Central limit theorem [60], can be approximated by N C (0, σ QF T ).The reduced density matrix of QF T |ψ⟩ can be represented as q 0 q 2 ... q 2 q 2 q 2 q 1 ... q 3 q 3 . . . . . . . . . . . .q 2 q 3 ... q 1 q 3 q 2 q 3 ...
where 1) q 0 ∼ λ 0 , 2) The symbol O gives an order of approximation since the elements of QF T |ψ⟩ are correlated and can not be treated as i.d.d random variables.When α ≫ 0, for a given λ 0 , we have q 0 ≫ q 1 , q 2 , q 3 , which indicate that λ QF T 0 the dominant eigenvalues of ρ QF T is very likely to be near λ 0 .We can observe from FIG. 17 that QF T |ψ⟩ and |ψ⟩ share the exact dominant eigenvalue and entropy, which follows the analytical result in Eq. (32).Therefore, the boundary of entanglement entropy defined with von Neumann entropy remains valid in the momentum space.Also, as we can observe from FIG. 16, the entanglement trajectory after performing QFT almost overlaps with the original one.

Entanglement gap
In this section, we investigate the relation between λ 0 and the entanglement gap.The reduced density matrix ρ A can be written as  where H A is called the entanglement Hamiltonian [61,62,63,64,65,66].The entanglement spectrum levels ξ i = − ln λ i are the "energies" of H A .The entanglement gap is a natural quantity defined as With this definition and the approach that we deduced f 1 , f 2 and f 3 , we can have a tight boundary for entanglement gap confined by g 1 , g 2 and g 3 , such that is plotted with a gray dash-dotted line, and is plotted with a gray solid line.
To calculate the asymptotic value of δξ, we take λ + the right end of the MPD as λ 1 , which is , and δξ can be written as a function of α, β and λ 0 with Eq. ( 29) From FIG. 11 and FIG. 18, we can observe a monotone decreasing relation between λ 0 and E λ 0 , and a monotone increasing relation between λ 0 and δξ.These two relations imply that the von Neumann entropy decreases with the entanglement gap, as stated by Li and Haldane in 2008 with the ideal Moore-Read state [61].Furthermore, the trajectory boundaries defined with von Neumann entropy (FIG.12) and the entanglement gap share a similarity.Intuitively, we can obtain FIG.19 Also, Numerically, the maximum entanglement gap is becoming closer to the upper bound given by MPD as n increases.
Three non-trivial observations can be obtained from this comparison.Firstly, the trajectories still remain within the flexible boundary even though, mathematically speaking, they are allowed to step outside.Secondly, since there is a monotone decreasing relation with E λ 0 and δξ of the boundary, a smaller entanglement gap is expected with higher entropy.The simulation results do not fit this prediction.Instead, the figure of merit is the relative distance to the boundary.Data points that are closer to the boundary when the trajectory is defined with the von Neumann entropy remain closer to the corresponding boundary when the trajectory is defined with the entanglement gap.Furthermore, comparing FIG.14, FIG. 15 and (b), (c) of FIG. 20, trajectories closer to the boundary defined by MPD present more value for quantum computation.Thirdly, the trajectory preserves its particular "topology" in each example over different entanglement measures.Intuitively, we pull the trajectory down on the left side and push up on the right.These three conclusions are distinct, but they support each other.
Here we highlight that these three properties are observations which emerge naturally from simulations of quantum computing and are impossible to prove mathematically.A trajectory plotted with von Neumann entropy does not correspond to a unique path when plotted with the entanglement gap, and vice versa.Counterexamples of these properties can always be designed.In the abstract formulation, the eigenvalues (λ 0 , ..., λ α−1 ) locate on a compact α − 1 dimension surface in the space R α constrained by λ 0 ≥ ... ≥ λ α−1 ≥ 0 and α−1 i=0 λ i = 1.The von Neumann entropy and the entanglement gap are projections of this high dimensional surface on R + , and most information about ρ A is lost.The unquoted word topology should not be used without a proper definition of homomorphism or homotopy.
These observations grant the entanglement trajectory an intrinsic physical value.Instead of just a tool for visualization, it reveals invariant properties over several examples with a variant measurement of entanglement, which is not guaranteed by mathematics.

Different degrees of Rényi entropy
In this section, we study the entanglement trajectory defined with different degrees of Rényi entropy.To simplify the analysis, we keep α = β.
The Rényi entropy is defined as [67] where d ̸ = 1 is the degree, and λ i are the eigenvalues of ρ A , with λ 0 ≥ ... ≥ λ α−1 and We have and These three equations are verified with numerical results in FIG. 21.We can recover Eq. ( 56) from Eq. ( 32) by setting α = β.Rényi entropy with higher integer degrees can be obtained using integrals in Therefore, instead of a data point on a 2D figure, a quantum state that is partially traced over a subsystem can be represented by a continuous and differentiable path located on the plane sectioned by x = λ 0 .As shown in FIG.22, we can have a 3D version of the entanglement trajectory defined by Rényi entropy, which is a extension of the 2D version and preserve its properties.We remind our readers that in the 2D version, lines that connect data points do not have a physical meaning.For the same reason, the surface formed by paths in this 3D version does not process any physical meaning either.

Conclusion
We introduce a new way, the entanglement trajectory, to visualize the transition of entanglement along a quantum process or a series of quantum states.We set a tight boundary to this trajectory from the definition of entropy and reduce it to a flexible boundary with asymptotic analysis based on random matrix theory.This part of the research remains valid for fundamental quantum physics research outside the context of quantum computation.Then, for several well-known examples from quantum computation, we simulate their entanglement trajectory and compare them with analytical boundaries.We observe some invariant properties, such that each trajectory preserves its "topology" and stays within the flexible boundary over variant scenarios and measurements of entanglement.These properties cannot be derived directly from the mathematical definition of entropy.Hence instead of a simple geometric visualization, the entanglement trajectory processes an underlying physical meaning and can provide more information about the quantum system than only studying the entanglement.It has the potential to help us understand more about the advantages of quantum computation since interesting instances reach near the boundary.Furthermore, the study of entanglement trajectory is inevitable if the advantage is truly justifiable because it provides a more general description of entanglement, the most distinguishable feature between classic and quantum despite different computation models.

Figure 2 :
Figure 2: The entanglement trajectories for solving three instances of EC problem with n = 12, we took three random bi-partition (where the subsystem A is a random half of qubits) of each instance.Different instances are plotted with markers +, × and ⋆.

Figure 3 :
Figure 3: The quantum circuit for Grover algorithm to solve one instance of EC problem with n = 10 and 7 clauses.The number of bits and clauses are chosen such that two Grover examples share the same number of qubits and can be traced over subsystems of the same size.We sample ρ A in each iteration before and after the central multi-Toffoli gate and the diffusion operator as indicated with dashed lines.In this circuit, we have t = 25.

Figure 4 :Figure 5 :
Figure4: The quantum circuit of Grover algorithm to find the pre-image of a Hash function with a given cipher-text.The oracle consists of a unitary that encodes the Hash function, the inverse of this unitary, with a multi-Toffoli gate in the middle that encodes the cipher-text, which is 10100011 in this figure.We sample ρ A in each iteration before and after the multi-Toffoli gate and the diffusion operator as indicated with dashed lines.This cipher-text has 2 pre-images, and the number of Grover iterations needed is t = 8.

Figure 6 :
Figure6: The semi-classical circuit for the Shor algorithm using 2n + 3 qubits.Details are provided in Ref.[25].|c⟩ is the controlling qubit that performs the QFT.Register |x⟩ consists of n qubits and encodes the output of the QME, and is initialized as |1⟩ (it means the least significant qubit is |1⟩ while others in |0⟩).Register |0⟩ consists of n + 2 qubits and serves as ancillas.

Figure 7 :
Figure 7: The decomposition of the controlled QME.CM is the circuit for the controlled multiplication.Details are provided in Ref. [25].The reduced density matrices ρ A are sampled before and after each time the controlled SWAP operator in the middle is applied, and the subsystem A is |x⟩, as shown with dashed lines.

Figure 8 : 1 2
Figure8: The entanglement trajectories of different semi-classical Shor circuits.For each n, three different N and three different a were chosen (nine combinations).Data that were taken after the controlled SWAP gate cluster on λ 0 = 1 2 (plotted as a black vertical dotted line) with an interval of 1 2 ln 2 (plotted with gray horizontal lines).And data that were taken before the controlled SWAP gate accumulate on the curve f shor (plotted as a black dashed line).Data points are plotted with high transparency to show the concentration of the trajectory.

Figure 9 :
Figure 9: The MPD (in violet) versus the numerical results (in blue) for different λ with α = 10000.The gray dotted line indicates λ + .

Figure 11 :
Figure 11: The numerical results (in data points) versus the analytical solution (in line).The von Neumann entropy of three different β ≥ α (30 data points for each β) for the same α = 128.Each of the 90 data points in the plot is the entropy of a randomly generated normalized reduced density matrix under the form ρ A = ZZ † / Tr ZZ † .

Figure 12 :Figure 13 :Figure 14 :Figure 15 :
Figure 12: The comparison between the tight and flexible boundaries for the entanglement trajectory defined with the von Neumann entropy.The vertical gray dotted line is λ + in Eq. (21) to indicate the lower bound of λ 0 for the flexible boundary, when the variance reaches the theoretical maximum, representing the case of average entropy.

Figure 16 :
Figure 16: The entanglement trajectories for k-almost Prime states and their unions when n = 26.The Prime state |P⟩ is plotted with black + and located on the intersection between |P k ⟩ and |U k ⟩.QF T |P k ⟩ and QF T |U k ⟩ are plotted with slightly bigger gray markers.They are very close to the data points of |P k ⟩ and |U k ⟩, therefore almost invisible.

Figure 17 :
Figure 17: The von Neumann entropy before and after applying the QFT of different β ≥ α for the same α = 128, 10 random states for each β.Notice that data points after applying the QFT are plotted with black markers and overlap completely with the original state.

)Figure 19 :
Figure 19: The comparison between the tight and flexible boundaries entanglement trajectory defined with the entanglement gap.
by pulling down the left side of FIG. 12 to 0 and pushing up the right side to ∞.The trajectories of previous examples are re-plotted with the entanglement gap in FIG.20.In the example of semi-classical Shor algorithm, f shor is reformulated into

α− 1 i=0 λ i = 1 .
When d = 0, E d = ln α and called Hartley entropy or max-entropy.Notice E d behaves very differently from d → 0 to d = 0.When d → 1, E d is the von Neumann entropy − α−1 i=0 λ i ln λ i [68], when d = 2, E d is the Rényi entropy or the Collision entropy.Furthermore, when d → ∞, E d is the min-entropy − ln λ 0 , which is another reason why λ 0 should be the x-axis, that it naturally gives the asymptotic analysis of Rényi entropy before we define entanglement trajectory.Since the series of Rényi entropy are continuous extensions of the von Neumann entropy of different d, we will no elaborate on the corresponding f 1 , f 2 and f 3 in detail and focus on the asymptotic value, which for d → 1 and d = 2 can be calculated analytically by integrating over the MPD using Eq.(55) b 0 ln x x (b − x)dx = π 16 2b 2 ln b − b 2 (4 ln 2 − 1) .
Corresponding to FIG. 13, adiabatic quantum computation to solve EC problem.Corresponding to FIG. 14, Grover algorithm to solve EC problem and the Hash function.Corresponding to FIG. 15, Shor algorithm.Gray horizontal lines are plotted with interval ln 2. Corresponding to FIG. 16, k-almost Prime states and their unions.