Concurrence percolation threshold of large-scale quantum networks

Quantum networks describe communication networks that are based on quantum entanglement. A concurrence percolation theory has been recently developed to determine the required entanglement to enable communication between two distant stations in an arbitrary quantum network. Unfortunately, concurrence percolation has been calculated only for very small networks or large networks without loops. Here, we develop a set of mathematical tools for approximating the concurrence percolation threshold for unprecedented large-scale quantum networks by estimating the path-length distribution, under the assumption that all paths between a given pair of nodes have no overlap. We show that our approximate method agrees closely with analytical results from concurrence percolation theory. The numerical results we present include 2D square lattices of 2002 nodes and complex networks of up to 104 nodes. The entanglement percolation threshold of a quantum network is a crucial parameter for constructing a real-world communication network based on entanglement, and our method offers a significant speed-up for the intensive computations involved. The move to the quantum internet demands developments in communication networks that are based on quantum entanglement. The authors discuss the phenomenon of entanglement percolation in a quantum network presenting solutions to significantly accelerate the intensive computation effort involved in the process.

T he application of network science to problems in quantum physics is a relatively new and rapidly developing field 1 .Quantum networks, where the links between nodes represent entangled qubits [2][3][4][5][6] , are expected to form the basis of the quantum internet.Recent advances in quantum repeater technology have made long-distance, noise-resilient quantum communication possible [7][8][9][10][11] .These networks have quantum correlations that can be exploited by performing specific local measurements on any node.
Protocols for quantum communication rely on the conservation of correlations in entangled states, and the generation and distribution of entanglement are necessary for quantum networks 12 .For a given network topology, we want to determine the minimum amount of entanglement necessary between qubits to maintain a giant component in the network, which is a problem analogous to percolation on classical networks [13][14][15][16] .However, there are crucial differences between classical and quantum networks, limiting the extent to which we can map a classical percolation theory to quantum networks.For example, in a classical random network with N nodes, if an edge between nodes exists with probability p, a subgraph with n nodes and l edges exists above a critical threshold of p given by p c ∝ N −n/l17 .For a quantum network, however, p c ~N−2 for all subgraphs for large N 18 .We can also use measurement strategies to alter the topology of a quantum network, meaning that the optimal entanglement percolation threshold needs to be minimized over all possible measurement strategies 3 .
Quantum networks 6 are used to model a scalable communication network based on quantum teleportation.It consists of nodes that denote a local set of qubits and edges, which represent a bipartite and entangled state of qubits shared between the two connected nodes (Fig. 1).The simplest practical quantum network can be built from quantum repeaters 19,20 which share only one entangled pair between nodes 5 .Quantum communication networks are expected to have several advantages over classical communication networks, including the ability to use quantum cryptography and send "quantum software" 2,6 .
The bipartite state of any entangled qubits in the network can be defined as up to a unitary transformation, where θ can change from 0 to π/4.Each entangled pair can be converted to a maximally entangled pair with a certain probability p, known as the singlet-conversion probability (SCP), given by p ¼ 2sin 2 θ 21 .This represents the probability of establishing a perfect communications channel 18 .
Converting every entangled pair in the network to a singlet is equivalent to a bond-percolation process and this measurement strategy is called classical entanglement percolation (CEP) 3,13 .
Having established the mapping to percolation, it is natural to ask what is the minimum level of entanglement necessary for the formation of a perfect communication channel between any two stations.Using classical-percolation arguments, we can establish the minimum level of entanglement necessary for establishing an infinite cluster under CEP for some simple cases.For example, it is θ th = π/4 for 1D chains and θ th = π/6 for 2D square lattices 3 .
Unfortunately, CEP does not give us the lowest possible percolation threshold value for a quantum network because it is possible to lower the entanglement threshold necessary for creating an infinite cluster by changing the network's topology.This is done through a process known as entanglement swapping, shown in Fig. 1, where two previously unentangled qubits are entangled using local operations and classical communication (LOCC) 22 .
The network topology may be altered to lower the percolation threshold before converting every link to a singlet by performing a series of entanglement swapping operations.This strategy is known as quantum entanglement percolation (QEP) 3,23 .The limitation of QEP is that it is not generally adaptable to arbitrary network topology as CEP is.For most network topologies, QEP cannot improve the percolation threshold in general.Note also that neither CEP nor QEP are optimal, meaning that it is impossible to determine if any given measurement strategy results in the lowest potential value of the percolation threshold 1 .A new local statistical theory, concurrence percolation theory (ConPT), has recently been proposed to explain the observed quantum advantage in quantum networks over the prediction of classical percolation theory 24 .This theory is analogous to classical percolation theory but fundamentally different, as it is not built on probability measure p but a new variable c for each link, which stands for concurrence-a key measure of bipartite entanglement 25 .However, ConPT can be computationally very expensive 24 .
Here, we present a fast and tangible solution for calculating the ConPT threshold.Our method relies on two approximations.The first is the parallel approximation, which treats all paths in the network as non-overlapping.The second is what we call the S m approximation, where we calculate the total concurrence between nodes using a subset of paths consisting of the m-th shortest paths on the network, with m = 1 referring to the shortest paths.We find that our approximate method agrees closely with the analytical results provided by Meng et al. 24 .Depending on the choice of S m , the computation based on this method can be several orders of magnitude faster than the analytical approach.By combining our method with combinatorial expressions for shortest and second-shortest paths we can calculate, for the first time, an approximation for the concurrence for much larger networks than would be analytically possible.Here, we calculate the sponge-crossing concurrence for 2D lattices with up to 200 2 nodes.We also extend the notion of concurrence to networks without boundaries and present results for Erdős-Rényi and Barabási-Albert networks of up to 10 4 nodes.Our results are summarized in Table 1.Each link in the network represents a pair of entangled qubits shared between the two connected nodes.The original network is represented by the black solid and dashed lines.Node D performs a local measurement (represented by the red rectangle) on the qubits it shares with nodes B and C, causing the qubits in B and C to become entangled, represented by the dotted red line.The previous entanglement, represented by the dashed black line, is lost.Entanglement swapping can modify the topology of a quantum network, changing its percolation threshold.In this example, a fully connected network is split into two disjoint components, consisting of nodes A-C (in turqoise) and nodes D-F (in gray).

Results and discussion
Concurrence percolation theory.For pure bipartite states, the concurrence is defined as 26 where ρ r is the reduced density matrix of one party (subsystem) of the bipartite state.For qubits, Eq. ( 1), the concurrence is simply c ¼ 2 cos θ sin θ: Meng et al. then use this quantity in place of probability to construct a concurrence percolation theory (ConPT) on arbitrary network topologies 24 .
To be more specific, recall that for classical bond percolation on a lattice, any link in the lattice is active with probability p-a variable that should be considered as the link weight.We may then define a "sponge-crossing" quantity, P SC , as the probability that at least one path connecting the two distant boundaries is fully active.As a function of p, P SC can be calculated by summing up all paths that connect the two boundaries of the lattice, following basic addition and multiplication rules of probability measures 13 .Essentially, we treat P SC as a "weighted sum of all paths".Now, given an n-node quantum network, G θ ðnÞ, where all the link weights are θ, by the CEP/QEP schemes we have the mapping p 2sin 2 θ (i.e., the SCP).From classical percolation theory, it is known that a minimum value of p exists, below which the sponge-crossing probability, P SC , becomes zero in the thermodynamic limit n → ∞: This minimum value, p th , is known as the percolation threshold.The ConPT is constructed differently, using the mapping c sin 2θ instead 24 .Still, an analogous quantity, C SC , referred to as the sponge-crossing concurrence can be defined as the weighted sum of all paths in terms of this new weight c 24 .It is believed that a nontrivial threshold on c also exists: such that c th is the minimum value of the concurrence c per link, below which C SC becomes zero when n → ∞.
It remains to show how the "weighted sum of paths" is calculated for ConPT.As a problem of path connectivity, the calculation turns out to closely resemble the analysis of an electrical resistor network, where a set of series and parallel rules are needed as the basic connectivity rules 24,27 .Fundamental quantum communication theorems demand that, for k links connected in series, the total concurrence must be given by The rule for k links connected parallel to each other is more involved, given by ; which yields . A caveat lies in the fact that if the network topology is not series-parallel 28 but has nontrivial loops (e.g., a bridge-circuit topology), then C SC cannot be calculated using only series and parallel rules.Higher-order connectivity rules are needed, of which general forms are unknown.There is, however, a heuristic way to approximate these higher-order rules: by employing the so-called star-mesh transform, all possible higher-order rules can be approximated using only the series and parallel rules 24 .
Equations ( 4), ( 5), together with the star-mesh transform, allow us to calculate the "weighted sum of paths" between arbitrary two nodes in a quantum network of arbitrary topology.Formally, we denote the two nodes as the source node (s) and the target node (t), respectively, and we define the final concurrence between them as the s-t concurrence, C st .Note that although s and t are named differently, they are symmetric and exchangeable.Hence, between any two nodes, a C st can be calculated by the connectivity rules mentioned above (see Fig. 2 for example).
On regular lattices, the sponge-crossing concurrence C SC can be calculated by contracting two separate boundaries into two "mega" nodes 24 and calculate the s-t concurrence between them.As we increase the network size n, a threshold c th will emerge, accompanied with a sudden jump of C SC as soon as the concurrence c per link becomes larger than c th .This observation supports the existence of ConPT.Also, the observed c th is significantly smaller than all previously known classicalpercolation-theory-based schemes 24 , exhibiting in large-scale quantum networks a quantum advantage that is purely structural.
Despite the fresh insights the ConPT has offered, it has two main limitations: 1.The heuristic approximation (star-mesh transform) used for higher-order connectivity rules is a double-recursive process that is computationally intensive, thus only feasible for networks of very small size.2. Although an s-t concurrence can be calculated between any two nodes in any network topology, the sponge-crossing concurrence C SC is only defined for regular lattices that have apparent boundaries, and thus so is c th .It is unknown how to define c th on complex network topology where we cannot define a boundary, and, provided a proper definition, how (non)trivial the numerical result of c th would be.Bethe Lattice (Cayley Tree) (L = 100, k = 3) 0.5 0.5 Bethe Lattice (Cayley Tree) The concurrence percolation threshold θ th determines the fundamental long-distance entanglement transmission capability in large-scale quantum networks, where each link is a bipartite state of qubits, Ψ j i ¼ cos θ 00 j iþ sin θ 11 j i, θ ∈ [0, π/4].We compare our results to those provided in Meng et al. 24 for the Bethe lattice and 2D square lattice.We also report results on Erdős-Rényi (ER) and Barabási-Albert (BA) networks.
Unlike cluster-based percolation theories, ConPT is based on path connectivity, which is arguably more general 24 .This is why c th simply cannot be defined by clusters like in classical percolation theory for complex network topology.A proper definition and feasible calculation of c th will be of great interest for the theory itself as well as for its applications.Below, we will show our suggested solution that can satisfactorily handle these two limitations of ConPT.
A fast and tangible solution for concurrence percolation.We start by generalizing C SC from being defined between two apparent boundaries to two arbitrary sets of nodes, denoted S and T. It is reasonable that Eq. (3) will yield a nontrivial ConPT threshold c th for this generalized C SC , as long as the lengths of all paths connecting S and T increase with the network size n.We contract the two sets S and T into two "mega" nodes, which amounts to erasing the internal network topologies of S and T, and then calculate the s-t concurrence between them.This provides us a definitive way of calculating C SC for arbitrary network topology and inferring c th from Eq. (3).
Our numerical computation of c th ("Fast ConPT computation") on large-scale quantum networks is further made possible by introducing two key simplifying approximations: the parallel approximation and the S m approximation Parallel approximation.In this Section, we introduce the parallel approximation, where we treat all paths connecting nodes of interest to be parallel, i.e., we treat them as if they have no shared edges.For an arbitrary network with n nodes and uniform edge-weights c, the parallel approximation C 0 SC of the true sponge-crossing concurrence between two sets of nodes, S and T, is given by where N l is the total number of self-avoiding paths of length l that connect s and t for all s ∈ S and t ∈ T, respectively.Equation ( 6) is the mathematical statement of the parallel approximation, indicating that we are taking each of the N l paths to be parallel (Fig. 3).We illustrate the approximation with a simple example, and then show that on series-parallel networks 28 the concurrence calculated under the parallel approximation forms an upper bound to the true concurrence.First, we consider the case where our network is essentially parallel, i.e., it can be expressed as the parallel combination of k subnetworks each with concurrence c i .In this case, the parallel approximation is exact: The more interesting case is that of an essentially series network, i.e., a network that can be decomposed as a combination of subnetworks in series.We consider an exemplary network that splits into k branches, each with concurrence c p i (Fig. 3).The concurrence of the segment before branching is c s .Following the series and parallel rules (Eqs.4, 5), the sponge-crossing concurrence from the left of this network segment to the right is Under the parallel approximation, the network is transformed so that the concurrence of the segment is given by There are three cases: We would like to show that the above expression is no less than unity by considering some limiting cases.When c p i ¼ 0 for all i it results in C SC ¼ C 0 SC ¼ 0. If we increase the concurrence of a single branch i to be greater than zero while holding the other branches to be zero, then the expression becomes or gðc s c p i Þð1 À gðc s c p i ÞÞ ¼ c 2 s gðc p i Þð1 À gðc p i ÞÞ.This is simply to state that, given only one branch, the parallel approximation is exact (since there is one and only one path).We now add a second branch, j, while leaving the remaining k − 2 branches with zero concurrence.Now we have For c s → 0, the expression above reduces to lim and for c s = 1, the expression is equal to unity.Hence, since the derivative of the expression w.r.t.c s , is nonpositive between 0 ≤ c s ≤ 1, we derive that C 0 SC is no less than C SC .This tells us that the addition of another branch results in the approximated concurrence being an upper bound for the true concurrence calculated through series-parallel rules.As more branches are added to the expression, the result can be easily generalized, so that As before, we consider For c s → 0, the expression reduces to lim In particular, the inequality holds because f ðc p 1 ; c p k Þ is a concave function defined on the simplex where it satisfies f ðc p 1 ; c p k Þ > 1=2 except at the vertices where all but one c 2 p i are equal to zero.Thus, due to the monotonicity of f, one requires 1/2 and we will revert to case 1.Therefore, there is an upper bound on c s .As c s approaches its upper bound, f ðc s c p 1 c s c p k Þ !1=2, and thus lim is nonpositive between c s → 0 and the presumed upper bound of c s , again we derive that Finally, since every series-parallel network can be decomposed into essentially series or parallel configurations, taken together we have shown that C 0 SC is an upper bound for C SC on series-parallel networks.Interestingly, as we will see, this upper bound seemingly becomes tighter as the network becomes larger.We hence expect that a new concurrence threshold on C 0 SC can emerge, which should numerically approach the true c th from below and match c th in the thermodynamic limit n → ∞. S m approximation.For most regular lattices and complex networks, however, the distribution of N l (Eq.6) is not trivial.When we look at arbitrary networks, the calculation for the spongecrossing concurrence is essentially a path-counting problem which may require approximation as well.
Although the literature of path counting on graphs is rich and well studied, we are not aware of any closed-form solutions for enumeration of self-avoiding walks of arbitrary length for even the simplest network (like 2D lattices) 29 .While approximate path enumerations exist for both 2D lattices 30 and random networks 31 , we find them impractical, since the concurrence calculation is very sensitive to N l for small l.Indeed, observation of Eq. ( 6) implies that a single path's contribution to the total concurrence decreases with increasing l and increases with increasing N l .Even though longer paths (l ≃ n) will outnumber shorter paths by several orders of magnitude, shorter paths will still contribute significantly more to the concurrence.
Based on this, if we define S m as the set which contains up to the m-th shortest paths (i.e., the shortest paths, the 2nd shortest paths, and so on up to the m-th shortest paths) between s and t for all s ∈ S and t ∈ T, then it is possible to approximate the sponge-crossing concurrence between S and T using only these paths.When m ¼ m max , S m becomes the set of all sponge-crossing paths.
In this Section, using our Fast ConPT computation, we present numerical results for different networks of large size n.We numerically estimate the finite-size ConPT threshold in terms of θ th 1 2 sin À1 c th , determining its position on the curve by matching the corresponding sponge-crossing concurrence at the half point, C SC = 1/2.Fig. 3 A schematic representation of the parallel approximation.The parallel approximation preserves the number of paths and their lengths, but ignores overlap between the paths.In this example, the source node (in green) is connected to target nodes (in blue) by k branches of the network (represented by the orange boxes), each with concurrence c p i .Each branch joins to a segment with concurrence c s (represented by the blue box).Under the parallel approximation the network is transformed so that the k branches have no overlap with each other.
Bethe Lattice (Cayley Tree).Given a finite Bethe lattice (i.e., a Cayley tree) with L layers and coordination number k 32,33 , all paths from the root node to any one of the boundary nodes have the same length, L. Since only one path exists from the root node to any node on the boundary, the number of paths of length L is There is no need to employ the S m approximation since all paths are exactly known.Only the parallel approximation C 0 SC of the sponge-crossing concurrence C SC is to be calculated, which is given by (following Eq. 6) To solve for c th , near given an arbitrarily small positive ϵ.This gives rise to and thus in the limit of large L.This is identical to the exact ConPT threshold calculated by Meng et al. 24 using a recursive renormalization trick on the series and parallel rules (Eqs.4, 5).Interestingly, it is known that a saturation point c sat < 1 also exists in ConPT 24 , namely, before c reaches unity, C SC will already reach unity at c = c sat .This is because of the maximum function appearing in the parallel rule (Eq.5).It is also obvious that c sat ≥ c th , given the monotonicity of the series and parallel rules.
To see if we can solve for c sat using the parallel approximation too, let We see that the saturation point calculated using the parallel approximation is equal to c th , which is underestimated, since the true saturation point is given by c sat ¼ ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ð1=2Þ 1=k À ð1=4Þ 1=k q = ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi ð1=2Þ ðkÀ1Þ=k À ð1=4Þ ðkÀ1Þ=k q 24 , which can be calculated similarly using a recursive renormalization trick on the series and parallel rules (Eqs.4, 5).
For validation purposes, numerical results of the spongecrossing concurrence on the Bethe lattice using the parallel approximation versus the true ConPT results are shown in Fig. 4. We see that as L increases, both c th and c sat (where c ¼ 2 cos θ sin θ) approach 1= ffiffiffiffiffiffiffiffiffiffi ffi k À 1 p from below and above, respectively, consistent with our theoretical result.Hence, it is highly suggested that the parallel approximation can correctly estimate the true ConPT threshold c th in the thermodynamic limit.
2D square lattices.In a 2D square lattice of n nodes ( ffiffiffi n p 2 Z), the length of the mth shortest self-avoiding path, between source and target nodes of coordinates s = (x s , y s ) and t = (x t , y t ) ), is simply Now, let S and T denote the left (x s = 1) and right (x t ¼ ffiffiffi n p ) boundaries.Let s = (1, y s ) ∈ S and t ¼ ð ffiffiffi n p ; y t Þ 2 T.Under the S m approximation, the total number of self-avoiding paths of length l between S and T is given by where δ ij is the Kronecker delta.This approximation of N l is then substituted into the parallel approximation (Eq.6) to calculate C SC between S and T.
For m ≤ 2, it is possible to directly enumerate the 1st and 2nd shortest self-avoiding paths between every pair of s and t.The general expressions are given by (w.l.o.g., x s ≤ x t and y s ≤ y t ) and where the boundary effect is also taken into account.In particular, Eq. ( 17) is obtained due to the fact that every 2nd shortest self-avoiding path in the square lattice, having length l 2 = |x s − x t | + |y s − y t | + 2, must contain one and only one of the configurations as shown in Fig. 5.The first and second terms in Eq. ( 17) account for the two "Z"-shape configurations (Fig. 5a, b), respectively; the third term for the two "L"-shape configurations (Fig. 5c, d); and the last term for the other two "L"-shape configurations (Fig. 5e, f).
Piecewise path enumeration algorithm.For m > 2, it becomes difficult to write down a closed-form combinatorial expression like Eqs. (16, 17) for N l m ðs !tÞ.A path enumeration algorithm is thus needed.We treat paths of length l m with m > 2 as deviations from the 1st and 2nd shortest paths.For a given m, these deviations can only take a finite number of shapes.Once we have identified these primitive deviations, we must next identify positions in the lattice where these deviations can be placed.Finally, we count the total number of paths by counting the number of shortest paths between deviations using Eqs.(16, 17).For example, given source and target nodes s and t, all 3rdshortest paths (m = 3) have either two single-step deviations or one double-step deviation from the 1st shortest path.For the case where we have two single-step deviations, we first identify two sets of points, D 1 and D 2 , where the first and second deviations can happen respectively.Then we calculate N s;D 1 (the number of shortest paths from s to every point in D 1 ), N D 1 ;D 2 (the number of shortest paths from every point in D 1 to every point in D 2 ), and N D 2 ;t (the number of shortest paths from every point in D 2 to t).The total number of 3rd-shortest paths is then given by N l 2 ðs !tÞ ¼ N s;D 1 N D 1 ;D 2 N D 2 ;t .This algorithm, while significantly faster than a brute-force path enumeration, is still too involved for large m.We use this algorithm to calculate S 3 exclusively.
Numerical calculations.The final numerical results of C SC , calculated using the exact combinatorial expressions (S 1 , S 2 ) and/or the piecewise path enumeration algorithm (S 3 ), are shown in Fig. 6.We can see from Fig. 6a that the transition in the value of the sponge-crossing concurrence becomes sharper as the network size increases.From Fig. 6b-d we see that for large enough m or n, the numerical threshold θ th levels out at constant values that are very close to those calculated using the star-mesh transform.For example, for n 2 = 8 the Fast ConPT method yields θ th = 0.4, compared to the value of θ th = 0.416 calculated using the starmesh transform 24 .This suggests that our Fast ConPT calculation can yield a good approximation of the ConPT threshold.We can also see from Fig. 7 that the Fast ConPT computation is over 100 times faster than the star-mesh transform method.
Complex network topologies.Unlike 2D square lattices, we cannot write down any analytical expressions for the path-length distribution of complex networks.While techniques to enumerate paths give a good estimate of the total number of paths 31 , they approximate the path-length distribution poorly.This means that we must enumerate paths through brute-force methods and this restricts our analysis to sparse graphs.
For complex networks, we simply define the sponge-crossing concurrence as the C st between two nodes s and t which means that S = {s} and T = {t}.We pick s and t such that the shortest path between them is equal to the diameter of the network.In general there might be multiple choices for s and t that meet this criterion, and we randomly choose one of these pairs.We randomly generate 100 network realizations of a given size and degree distribution and average the concurrence percolation threshold of all networks.These results are reported in Table 1 along with the standard error, σ= ffiffiffiffi N p , where σ is the standard deviation and N = 100 is the number of realizations of each random graph.
Erdős-Rényi Network.Results for Erdős-Rényi (ER) networks 34 are shown in Figs. 8, 9a, 10a, b.The concurrence is calculated under the S 6 approximation for different settings of network sizes and average degrees.The results are averaged over 100 network realizations for each setting.For small values of m, the behavior of the concurrence for ER networks can be approximated with a power-law fit, as shown in Fig. 8. Figure 10a shows that the value of θ th converges with increasing network size for smaller values of the average degree, e.g., k h i ¼ 4. For larger values of the average a b Fig. 5 Characterization of second-shortest paths using primitive steps.Every 2nd shortest self-avoiding path must contain one and only one of the configurations (solid line): either "Z"-shape (a, b) or "L"-shape (c-f), then the rest connected by shortest paths (dashed lines).
degree, such as k h i ¼ 8, we do not see the value of the threshold converging for the same network size.
The calculation of the number of paths becomes increasingly computationally intensive for larger values of the average degree and we must restrict our analysis to small values of m. Figure 10b shows how the concurrence percolation threshold changes as a function of average degree under the S 1 approximation.
Barabási-Albert network.Many real-world networks show power-law degree distribution, such as the Internet, WWW, scientific collaboration networks, protein-protein interaction networks, and actor networks 14,35 .Barabási-Albert (BA) model 36 is the first model to describe the structure property of such networks, using preferential attachment.In this model, every new node in the graph is assigned z edges, where z is known as the coordination number, and nodes with higher degrees are more likely to be selected.The classical bond-percolation threshold for a BA network with z > 1 and n → ∞ is p c = 0 37,38 .
Results for BA networks are shown in Figs.9b, 10c, d.For z = 1 there are no loops in the network and the relatively small number of paths connecting any two nodes allows us to calculate the concurrence for up to 10 4 nodes, shown in Fig. 10c.We also look at smaller networks with higher coordination numbers, up to z = 25, shown in Fig. 10d.Unlike ER networks, the value of θ th decreases with the increasing network size.
Comparison with Classical Entanglement Percolation (CEP).As a baseline comparison, we numerically calculate θ CEP th , the percolation threshold associated with classical entanglement percolation 3,13 (CEP) for ER networks.As before we define the percolation threshold on random networks as the minimum entanglement necessary for the existence of a path between two nodes s and t, where s and t are a randomly selected pair of nodes with the property that their distance is the diameter of the network.We generate 100 random networks and eliminate edges with probability 1 − p, where p ¼ 2sin 2 θ is the singlet-conversion probability.For each network we perform 1000 simulations and calculate θ CEP th as the average minimum value of θ such that the probability of the existence of a path Fig. 7 The speed-up obtained by the Fast ConPT algorithm over starmesh reduction.The figure shows the computing time (in seconds) to calculate C st , the s-t concurrence between two nodes s and t, on 2D square lattices with n nodes using the Fast ConPT method with the S 1 and S 2 approximation, as well as using the star-mesh reduction method of Meng et al. 24 .We can see that the Fast ConPT method speeds up the calculation over the star-mesh transform method by two orders of magnitude.
Fig. 8 The concurrence percolation threshold as a function of m for Erdős-Rényi (ER) networks.The relationship can be approximated with a power law, θ th (m) ∝ m ϕ , which is represented by the dashed black lines.The networks have size n = 100.between s and t is 0.5.We also calculate θ th (m) for these networks.The results are shown in Fig. 11.We can see that for all our samples θ th ðm ¼ 1Þ < θ CEP th .Since CEP represents the naive, baseline measurement strategy it is encouraging that our approximate ConPT threshold always lies below it even when we restrict ourselves to m = 1.Therefore θ th heuristically approaches a lower bound on the true concurrence and even for low values of m it predicts a lower concurrence percolation threshold than that predicted by CEP.

Conclusion
Table 1 summarizes the numerical results of our fast concurrence percolation theory (Fast ConPT) computation compared with previously known results 24 .The algorithm we have presented in this report utilizes two approximations to allow for numerical calculations of ConPT.Where available, our results are in good agreement with the analytical values of the concurrence percolation.We have also extended the analysis of the ConPT threshold to complex networks and demonstrated that our method could be applied to square lattices of 200 2 nodes and complex networks of 300 nodes.Combining our method with more efficient path-counting algorithms would allow us to probe a more significant fraction of the total paths of a network for the ConPT calculation and provide a more robust estimate for the ConPT threshold.We believe that this work is an important step towards understanding the structural and communication properties of large-scale quantum networks.
We believe that ConPT is a promising tool for practically designing and analysing quantum networks.It offers crucial insights into how entanglement strength, viewed as a costly resource, should be distributed throughout a network to ensure resilient communication.In full-optical quantum communication a b Fig. 9 The s-t concurrence as a function of entanglement strength for complex networks.The solid lines show the curves averaged over 100 network realizations and the bands represent one standard deviation.We calculate C st , the s-t concurrence, between two maximally separated nodes as a function of link weight, θ, for (a) Erdős--Rényi (ER) with k h i ¼ networks, for example, entanglement strength is usually expressed as a function of the number of entangled photons shared between nodes 39 .Current methods used in simulations to determine the entanglement strength necessary for the emergence of a giant component in quantum networks or the connectivity of two random nodes implicitly assume the CEP measurement strategy, i.e., they assume that the topology of the quantum network is immutable like that of a classical network 39,40 .The numerical methods we have presented in this paper allow concurrence percolation theory to be practically useful in the analysis of large complex networks, providing a lower bound on the entanglement strength necessary for communication between distant nodes, therefore allowing the cost associated with establishing communication channels of a certain strength in these networks to be lowered.As we showed in our comparison with CEP, even for m = 1 the Fast ConPT method already predicts a lower entanglement percolation than CEP, demonstrating its effectiveness for determining how close any given measurement strategy is to being optimal.
Still, the critical behaviors of ConPT near θ th remain an interesting and open question.Previous studies have indicated that some critical phenomena, such as the emergence of subgraphs, are drastically different in quantum random networks than in classical networks 6 .While Meng et al. has provided a finite-size analysis of the critical behaviors of Bethe lattices, their analysis of 2D lattices is limited by the size of the lattices they could investigate 24

Fig. 1
Fig. 1 Schematic representation of a quantum network.The larger circles represent nodes.The smaller, enclosed circles represent entangled qubits.Each link in the network represents a pair of entangled qubits shared between the two connected nodes.The original network is represented by the black solid and dashed lines.Node D performs a local measurement (represented by the red rectangle) on the qubits it shares with nodes B and C, causing the qubits in B and C to become entangled, represented by the dotted red line.The previous entanglement, represented by the dashed black line, is lost.Entanglement swapping can modify the topology of a quantum network, changing its percolation threshold.In this example, a fully connected network is split into two disjoint components, consisting of nodes A-C (in turqoise) and nodes D-F (in gray).

Fig. 2
Fig.2Demonstration of calculating the s-t concurrence C st .The s-t concurrence between nodes 1 and 6 on a 2D rectangular lattice can be calculated using ConPT, in the order from (a) to (e).The edge weight θ characterizes the entanglement strength between a pair of entangled particles shared between nodes (Eq. 1) and c sin 2θ refers to the corresponding concurrence.a Original lattice.b, c Series rules.d Star-mesh transform on the star-graph (edges 4 ↔ 1, 4 ↔ 3, 4 ↔ 6), then parallel rule for edges 1 ↔ 3 and 3 ↔ 6. e Series rule for edges 1 ↔ 3 and 3 ↔ 6, then parallel rule for edge 1 ↔ 6.

Fig. 4
Fig. 4 The sponge-crossing concurrence for the Bethe lattice under the parallel approximation.Results are shown for coordination numbers (a) k = 3 and (b) k = 4.As the number of layers, L, in the network become larger the numerical values of θ th approaches the analytical value.The solid black lines represent the true concurrence values for the Bethe lattice 24 .

Fig. 6 Fast
Fig.6Fast ConPT calculation on 2D square lattices.a Sponge-crossing concurrence C SC as a function of link weight θ, calculated under the S 1 -S 3 approximations.Only the result of S 3 is plotted.The results of S 1 and S 2 are nearly identical to S 3 .b Numerical ConPT threshold θ th under the S m approximation.As m increases, θ th approaches a constant value.c θ th for different size n.d Same as (c) but for larger n. S 3 becomes too computationally intensive to calculate for n > 202 .As n increases, θ th also approaches a constant value.

Fig. 10 The
Fig.9The s-t concurrence as a function of entanglement strength for complex networks.The solid lines show the curves averaged over 100 network realizations and the bands represent one standard deviation.We calculate C st , the s-t concurrence, between two maximally separated nodes as a function of link weight, θ, for (a) Erdős--Rényi (ER) with k h i ¼ 4 and (b) Barabási--Albert (BA) networks with coordination number z = 2, calculated under the S 1 approximation.

Fig. 11 A
Fig. 11 A comparison of the Fast ConPT algorithm and classical entanglement percolation (CEP).We plot θ CEP th , the entanglement percolation threshold calculated under CEP, against θ th (m), the approximate concurrence percolation threshold under the S m approximation for Erdős--Rényi (ER) networks with k h i ¼ 3.Each dot represents a network realization, with a total of 100 network realizations.The black line represents x = y and is added for comparison.a, b For small network sizes (n = 100) we can choose high enough value of m such that we expect θ th (m) to be a good approximation to θ th (∞), the threshold calculated with all paths.a Shows results for m = 1 while (b) shows results for m = 10.c, d For larger networks (n = 1000) we are restricted to smaller values of m. c Shows results for m = 1 while (d) shows results for m = 6.

Table 1
Concurrence percolation thresholds of different network topologies.