text:
# Reduction of Spin Glasses applied to the Migdal-Kadanoff Hierarchical Lattice

## I Introduction

We propose a set of reduction rules applicable to spin glasses at $`T=0`$ on any sparse graph and arbitrary bond distribution. These reductions strip graphs of all variables that are connected to at most two neighbors while accounting exactly for the ground state energy, entropy, and approximately for the overlap distribution of the system. In this paper, we introduce these reduction rules, and test them on the hierarchical lattice (see Fig. 1) obtained from the Migdal-Kadanoff bond-moving scheme MK . The recursive nature of the hierarchical lattice permits us to quickly reduce the equivalent of $`10^9`$ graphs with lengths corresponding to $`L=2^{100}`$ in $`d=3`$ and $`L=2^{35}`$ in $`d=4`$ dimensions, limited only by accumulating rounding errors. Throughout this paper we focus exclusively on $`d=3`$ and 4, and use only a discrete bond distribution,

$`P(J)=p\delta (J^21)+(1p)\delta (J).`$ (1)
But our procedure is equally applicable to any continuous bond distribution, such as a Gaussian.

The reduction produces high-accuracy results for the scaling of the defect energy width $`\sigma (\mathrm{\Delta }E)`$ with $`L`$. We also study the diluted lattice, where $`p`$ shall refer to the fraction of occupied bonds, and show that this scaling emerges for all $`p>p^{}`$. The critical point $`p^{}=0.31032`$ BF  is particular to the bond distribution in Eq. (1) and located just above the percolation point $`p_c=0.281837`$ of the lattice. Since $`p_c`$ is a purely topological property of the lattice itself, having $`p>p_c`$ is merely a necessary condition for long-range correlated behavior in a spin glass. Cooperative effects from the bond disorder suppress correlations even for $`p_c<p<p^{}`$, as was already discussed in Ref. BF . Correspondingly, the moment $`|q|`$ of the overlap distribution MPV  becomes non-zero only at $`p^{}`$, while the ground state energy and entropy remain smooth for all $`0p1`$, even at the transition. Both, entropy and overlap, exhibit an extremum at some $`p>p^{}`$ which can be explained in terms of the peculiar lattice hierarchy.

This study suggests that the reduction method provides a powerful means to explore the (potential) onset of replica symmetry breaking (RSB) MPV  on finite-dimensional spin glasses for arbitrary bond distributions at $`T=0`$ but variable $`p`$. The onset of glassy behavior (albeit not RSB) in the Migdal-Kadanoff model at $`p^{}`$, close to and intrinsically linked with $`p_c`$, may suggest such a close link between percolation and the possible onset of RSB in realistic lattice models, which we are currently exploring BoPe4 . In fact, we already found a transition for the stiffness exponent in real three- and four-dimensional Edwards-Anderson models, resembling the one observed here, at a $`p^{}>p_c`$ eo\_y .

Considered in terms of the average connectivity $`\alpha =2dpO(1)`$ of spins, dilute lattice glasses as well as the hierarchical lattice at $`T=0`$ behave strikingly similar to many mean-field models, which are currently studied in the context of combinatorial optimization Monasson ; Leone . Those models also enter into a glassy state (exhibiting RSB) at some finite connectivity $`\alpha ^{}`$, above a finite-connectivity percolation transition $`\alpha _c<\alpha ^{}`$. But in mean-field models the onset of RSB is often studied as a precursor to yet another critical point, the satisfiability transition AI , whereas geometrically defined models are unsatisfiable at any finite connectivity. These features should be independent of the bond distribution, assuming zero mean and unit variance.

While the hierarchical lattice is completely reducible, diluted realistic lattices will become irreducible soon above percolation, in which case our algorithm eo\_y  is used to reduce any graph as much as possible first, followed by a complete exploration of the ground states of the far more compact remainder graph using the extremal optimization heuristic eo\_prl . With this hybrid approach, we are typically able to compute ground state energies, entropies and overlaps for graphs 10-100 times larger than previously recorded, extending well beyond the percolation point.

## II Migdal-Kadanoff Hierarchical Lattice

To illustrate the reduction method, and to test our algorithm, we consider here the hierarchical lattice (see Fig. 1), obtained in the Migdal-Kadanoff real-space renormalization scheme MK  for low-dimensional spin glasses. These lattices have a simple recursive, yet geometric, structure and are well-studied Kirkpatrick ; SY ; BF ; BM ; DM ; BKM . Most importantly, they are completely reducible, and we can discuss the method independently of any subsequent optimization that may be required for more complicated models. The most interesting property of these lattices is the curious fact that the scaling of its defect energy distribution behaves very similar to that measured for actual two- and three-dimensional lattices BM ; F+H ; eo\_y .

As described in Fig. 1, starting from generation $`I=0`$ with a single link, at each subsequent generation $`I+1`$, all links from $`I`$ are replace with a new subgraph. The structure of the subgraph arises from the bond-moving scheme in $`d`$ dimensions, and has $`2^d=8`$ links for $`d=3`$ here. Thus, a hierarchical lattice of generation $`I`$ has $`l_I=(2^d)^I`$ links, thus corresponding to a $`d`$-dimensional lattice of “length” $`L=2^I`$ but $`n_I=2+2^{d1}(L^d1)/(2^d1)=O(L^d)`$ vertices. While the average connectivity is $`2l_I/n_I42^{2d}`$, the two root-vertices from generation $`I=0`$ themselves obtain in generation $`I1`$ a connectivity of $`2^{(d1)I}`$, and $`2^{dI1}`$ vertices are only two-connected, i. e. 7 in 8 for $`d=3`$.

The diluted hierarchical lattice percolates when there is a path between the two root-vertices. This notion leads to a simple recursion relation for the percolation threshold by counting the weights of all diluted subgraphs from Fig. 1 that percolate:

$`p_{I+1}=4p_I^26p_I^4+4p_I^6p_I^8,`$ (2)
which has a non-trivial stationary point at $`p_c=0.2818376366`$. It has been pointed out by Ref. BF  that a spin glass with the discrete $`\pm J`$-bond distribution in Eq. (1) exhibits instead a critical transition between a paramagnetic and a spin glass phase at $`p^{}=0.31032`$, which is closely related to the percolation transition. While below $`p_c`$ disconnected bonds ($`J=0`$) clearly dominate and prevent long-range correlations, even for $`p_c<p<p^{}`$ such correlations remain suppressed due to the cooperative behavior of parallel bond structure pervasive in the lattice that lead to many cancellations (see “double rule” in Sec. III) and additionally disconnects subgraphs at some higher level of the hierarchy.

## III Reduction Rules for Spin Glasses

In the Ising spin glass problem we assign to each vertex $`i`$ of a graph or lattice, $`1in`$, a spin variable $`x_i\{1,+1\}`$, while each link between two connected vertices $`i`$ and $`j`$ obtains a bond variable $`J_{i,j}`$ that is drawn at random from some distribution $`P(J)`$, which may be discrete or continuous. The problem consists of finding spin configurations $`(x_1,\mathrm{},x_n)`$ which minimize the Hamiltonian, or total energy,

$`H(x_1\mathrm{},x_n)={\displaystyle \underset{<i,j>}{}}J_{i,j}x_ix_j`$ (3)
for a fixed (quenched) set of bonds $`J_{i,j}`$. Of course, for $`J_{i,j}>0`$ (only ferromagnetic bonds), the solution is simply to make all connected spins aligned with each other. But if some or all bonds are negative, and spins are sufficiently connected, the problem can be frustrated Toulouse  and solutions become nontrivial.

In the following, we explain the rules by which to eliminate recursively all two-connected spins while accounting for the energy, entropy, and approximate overlap of ground-state configurations. To this end we have to generalize bonds by adding internal degrees of freedom, which evolve during the reduction process. Each bond between two spins $`i`$ and $`j`$ now consists of a tuple

$`𝒥_{i,j}=(J_{i,j},m_{i,j}^+,m_{i,j}^{},s_{i,j}^+,s_{i,j}^{}),`$ (4)
where $`J_{i,j}`$ is the weight of the bond, $`m_{i,j}^\pm `$ is the entropy accumulated by the bond, and $`s_{i,j}^\pm `$ counts the number of previously reduced spins whose state is completely determined (“entrained”) by the bond for $`x_ix_j=\pm 1`$, respectively. In addition, there is an energy offset $`E_{\mathrm{offset}}0`$ accounting for the energy difference between the original and the reduced graph. In general, determining entropy and overlap from the $`m`$\- and $`s`$-values can be quite complicated for the interesting case of an irreducible remainder graph BoPe4 . For the completely reducible hierarchical lattice, it is a lot simpler and we will describe it below. Initially, for the unreduced graph, all $`m`$\- and $`s`$-values and $`E_{\mathrm{offset}}`$ are zero, and the bonds are drawn from one of the usual distributions, for example, Eq. (1).

These rules are elaborate but easily enumerated, independent of the previous or future structure of the graph or the spin configuration. A more complete set of rules including zero-, one-, and three-connected vertices required for arbitrary graphs will be presented elsewhere eo\_y ; BoPe4 . There, each spin $`i`$ obtains one additional number $`s_i`$ (aside from its state $`x_i`$) that counts the previously eliminated spins whose state is entrained to $`i`$.

To demonstrate the use of these extra degrees of freedom, consider, for instance, the reduction of a one-connected spin $`j`$: It and its bond $`𝒥_{i,j}`$ to its sole neighbor $`i`$ disappears from the graph, and the connectivity of $`i`$ is decremented. The new graph will have itself a ground state energy offset by $`|J_{i,j}|`$ relative to the unreduced graph, i. e. $`E_{\mathrm{offset}}E_{\mathrm{offset}}|J_{i,j}|`$. For $`J_{i,j}=0`$ (a possibility we eventually have to allow for, even if our unreduced bond distribution did not), all information contained in $`𝒥_{i,j}`$ is discarded and $`s_i`$ doesn’t change, since $`x_j`$ is free to take on any value $`\pm 1`$ in the ground state; only the entropy of the old graph is offset by $`\mathrm{ln}2`$ with respect to the new graph. In turn, for $`J_{i,j}>0`$ ($`J_{i,j}<0`$), $`j`$ is entrained to the state of $`i`$ in any ground state configuration, and $`s_i`$ is incremented. If the bond $`𝒥_{i,j}`$ had previously inherited nonzero entropies $`m_{i,j}^\pm `$ and entrainments $`s_{i,j}^\pm `$, there would be an offset in entropy by $`m_{i,j}^+`$ ($`m_{i,j}^{}`$) between the reduced and the original graph, and $`s_i`$ would be further increased by $`s_{i,j}^+`$ ($`s_{i,j}^{}`$); $`m_{i,j}^{}`$ ($`m_{i,j}^+`$) and $`s_{i,j}^{}`$ ($`s_{i,j}^+`$) are discarded because the ground state must have $`x_ix_j=+1`$ ($`x_ix_j=1`$).

With the two-point rule we use exclusively in this paper, a spin variable $`x_i`$ and its two bonds $`𝒥_1`$ and $`𝒥_2`$ are replaced by a new bond $`𝒥`$. It is easy to show that

$`J={\displaystyle \frac{1}{2}}\left(\left|J_1+J_2\right|\left|J_1J_2\right|\right),`$ (5)
which assigns to $`J_{j,k}`$ the bond with the smaller modulus and the sign of $`J_1J_2`$. The energy offset $`E_{\mathrm{offset}}`$, accounting for the energy of the original graph, will have to be lowered by

$`{\displaystyle \frac{1}{2}}\left(\left|J_1+J_2\right|+\left|J_1J_2\right|\right).`$ (6)
At this level, our reduction rules correspond to the familiar traces (at $`T=0`$) used many times before to study the hierarchical lattice (see, e. g., Refs. DM ; F+H ).

Clearly, for $`|J_1||J_2|`$, the reduced spin $`i`$ is entrained to the neighbor with the bond of larger modulus, while the weaker bond may be violated. The entropy and entrainment variables that $`𝒥`$ inherits from the two reduced bonds are controlled by the sign of the stronger one. We can summarize these rules by defining $`𝒥_>`$ and $`𝒥_<`$ as the bond with the stronger and weaker (in modulus) weight, respectively. Then, we obtain for the new bond:

$`\underset{¯}{|J_>|>|J_<|}`$ $`:`$
$`m^\pm `$ $`=`$ $`m_>^{sign(J_>)}+m_<^{\pm sign(J_>)},`$
$`s^\pm `$ $`=`$ $`s_>^{sign(J_>)}+s_<^{\pm sign(J_>)}.`$ (7)
The case $`|J_1|=|J_2|`$ leads to further distinctions:

$`m^\pm `$ $`=`$ $`\{\begin{array}{cc}m_1^{sign(J_1)}+m_2^{sign(J_2)}\hfill & ,J_1=\pm J_20\hfill \\ \mathrm{ln}\left(e^{m_1^++m_2^\pm }+e^{m_1^{}+m_2^{}}\right)\hfill & ,J_1=J_2\hfill \end{array}`$
$`s^\pm `$ $`=`$ $`\{\begin{array}{cc}s_1^{sign(J_1)}+s_2^{sign(J_2)}\hfill & ,J_1=\pm J_20\hfill \\ 0\hfill & ,J_1=J_2\hfill \end{array}`$ (8)
Finally, we need a rule to combine double bonds, which will be created eventually. When two bonds, $`𝒥_1`$ and $`𝒥_2`$, connect the same two spins, they are simply replaced with a new bond

$`𝒥=𝒥_1+𝒥_2.`$ (9)
It is important to note that the two-point and double rules only effect the newly created bond and the energy offset, but not the two neighboring spins (or the entropy offset). Unlike for the one-point rule, we do not need to consider an internal entrainment number for spins here.

## IV Algorithm and Implementation

With the reduction rules from Sec. III, it is now straightforward to evaluate hierarchical lattices of any size. In effect, these rules allow us to reverse the generating mechanism depicted in Fig. 1. Recursively, from one generation to the next, we reduce all two-connected spins in the middle of the subgraph first, then use the double-bond rule sequentially, until a single, reduced bond is left. Once an entire graph has been reduced into such a single bond, the weight $`J`$ of that bond is in fact half the defect energy, $`\mathrm{\Delta }E=2J`$, of that graph, i. e. the energy difference between having both root spins aligned and anti-aligned. The entropy of that graph is simply given by $`m^{sign(J)}`$, or by $`\mathrm{ln}[\mathrm{exp}(m^+)+\mathrm{exp}(m^{})]`$ for $`J=0`$. The energy of the graph can be obtained from the running energy-offset $`E_{\mathrm{offset}}`$, accounting for the reduction of bonds along the way, and $`J`$. For our purposes, it is more useful to extract the (non-negative) “cost” of a graph given by the absolute weight of all violated bonds for any bond distribution. In case of the discrete $`\pm J`$ distribution in Eq. (1), this sum reduces to a count of all violate bonds (all of absolute weight 1) $`(l+E_{\mathrm{offset}}|J|)/2`$, where $`l=pL^d`$ is the average number of non-vanishing bonds in a graph. As an approximate measure of the overlap, we merely store the largest entrainment number $`s_{\mathrm{max}}`$ observed at any point during the reduction process, which is a measure of the largest correlated cluster of spins within the lattice. The size of this cluster will either dominate the overlap (with smaller clusters adding some negligible fluctuations) in the spin glass state or be itself a vanishing fluctuation in the paramagnetic state.

Although these measurements are not difficult, the reduction of ever larger graphs soon becomes untenable, since the cost of reducing true graphs would grow by a factor of $`2^d`$ for each generation. Instead, we exploit the fact that each subgraph within the hierarchy is independent of any other parallel subgraph. As above, we reduce all graphs down to a single bond. To assemble a graph of generation $`I+1`$ we only need a large enough pool of independent graphs of generation $`I`$. In fact, many nearly-independent graphs of generation $`I+1`$ can be assembled by repeatedly drawing at random on $`2^d`$ graphs of generation $`I`$ from that pool.

We proceed as follows: Assume that we already possess a pool of $`A_I`$ reduced graphs at every generation $`I`$ (up to a maximal generation $`I_{\mathrm{max}}`$), each graph represented by its bond $`𝒥`$, $`E_{\mathrm{offset}}`$, and $`s_{\mathrm{max}}`$. Each cycle of the algorithm generates an elementary subgraph (generation $`I=1`$) with bonds $`J`$ drawn from the bond distribution in Eq. (1) and empty entropy, entrainment, $`E_{\mathrm{offset}}`$, and $`s_{\mathrm{max}}`$ variables. The elementary subgraph is reduced to obtain a bond $`𝒥`$ of the first generation, replacing the oldest entry in pool $`A_1`$. After every $`k`$th new reduced bond entering into pool $`A_{I1}`$, we randomly assemble $`2^d`$ such bonds into a new subgraph to create a new reduced bond to be entered into pool $`A_I`$, as long as $`II_{\mathrm{max}}`$.

If we choose $`A_I=k=2^d`$ and assemble bonds sequentially instead of at random, this algorithm produces exact hierarchical lattices, where each subgraph is certain to be independent. But to obtain just one graph of, say, generation $`I=10`$, we have to assemble $`2^{30}`$ elementary subgraphs for $`d=3`$. Instead, for all $`p`$ not too close to $`p^{}`$, we simply choose large pools, $`A_I=2048`$, and as small as $`k=1`$. That is, each time we assemble an elementary graphs at $`I=0`$, we enter its value to pool $`A_1`$, assemble a subgraph out of $`A_1`$, add the reduced bond to pool $`A_2`$, etc, all the way up to $`I_{\mathrm{max}}`$ in each cycle. Hence, our implementation becomes independent of $`L`$, and for instance at $`p=1`$ we have assembled graphs up to generation $`I_{\mathrm{max}}=100`$ in $`d=3`$ and $`I_{\mathrm{max}}=35`$ in $`d=4`$, corresponding to $`L=2^{100}`$ or $`2^35`$, respectively. $`I_{\mathrm{max}}`$ is only limited by apparent limits on floating point accuracy due to repeated cancellation between almost equal bonds. Note that each pool gets completely refreshed after $`A_I`$ creations of new graphs from the pool just below it, while there are about $`O(2^{dA_I})`$ different subgraphs that could be created out of each pool at any one time.

Every time a graph gets added to a pool, its properties are evaluated from its $`𝒥`$, $`E_{\mathrm{offset}}`$, and $`s_{\mathrm{max}}`$. While ground state energy, entropy, and (approximate) overlap converge quickly with $`I`$, and throughout are quite independent of the choice of $`k`$ and $`A_I`$, the measurement of the defect energy becomes quite sensitive to the algorithmic parameters near $`p^{}`$. In those cases, we increase $`k`$ (and lower $`I_{\mathrm{max}}`$ correspondingly), until stable results were obtained. Only at $`pp^{}`$, the exact algorithm was necessary to observe the true defect energy distribution.

## V Numerical Results

Using the algorithm from Sec. IV, we have studied various properties of a Ising spin glass with bond distribution as given in Eq. (1) on the diluted and undiluted hierarchical lattice. In fact, using the easiest implementation with $`k=1`$, pool sizes of $`A_I=2048`$, and $`I_{\mathrm{max}}=20`$, we first scanned lattices for $`0<p1`$ in steps of $`\mathrm{\Delta }p=0.1`$. At each $`p`$ and $`I`$, we averaged over $`10^7`$ values. We found that for each $`p`$, the ground state energies per spin, entropies per spin, and the moment of the overlap, $`|q|`$ were well-converged at $`I_{\mathrm{max}}=20`$, and insensitive to the choice of $`k`$ using trials with $`k=2`$.

We also found the width of the defect energy distribution, $`\sigma (\mathrm{\Delta }E)=\sqrt{\mathrm{\Delta }E^2\mathrm{\Delta }E^2}`$, to be stable except for $`p=0.3`$. For large $`I=\mathrm{log}_2(L)`$ the distribution of defect energies either collapses into a $`\delta `$-function at zero for small $`p`$, or eventually converges to the distribution of the undiluted lattice ($`p=1`$). Near $`pp^{}`$, the evolution of that distribution becomes highly sensitive to small fluctuations, and the minute statistical dependencies inherent in our algorithm for small $`k`$ escalate. Thus, we ran another scan for $`0.28p0.34`$ using $`k=2`$ and $`10^4`$ graphs at $`I_{\mathrm{max}}=20`$ (and hence gaining a factor $`k`$ more graphs at each lower generation $`I`$). Now, only for $`p=0.31p^{}`$ did $`\sigma (\mathrm{\Delta }E)`$ become unstable (by comparison with trial runs at $`k=4`$). Finally, for $`p=0.31`$, we used an exact algorithm and indeed found $`\sigma (\mathrm{\Delta }E)`$ to be horizontal for many generations until it eventually decays, indicating that $`p^{}0.31`$, consistent with the expected value of $`0.31032`$ from Ref. BF . The results for $`\sigma (\mathrm{\Delta }E)`$ as a function of $`I=\mathrm{log}_2(L)`$ for all $`p`$ are shown in Fig. 2.

In Figs. 3-5 we plot the number of violated bonds (or cost) per spin, entropy per spin, and a characteristic moment of the (approximate) overlap distribution for $`L\mathrm{}`$ and $`d=3`$. The cost is already non-zero for small $`p`$, rising consistent with $`p^4`$, proportional to the probability of creating a small, frustrated loop that requires four bonds. Just beyond $`p^{}`$, the cost per spin rises more slowly, but without drastic changes for any $`p`$. Similarly, the entropy density in Fig. 4 remains unaffected by the transition at $`p^{}`$, as is expected Zecchina . Instead, it reaches a minimum at about $`p=0.7`$ that is an artifact of the hierarchical lattice. The large majority of spins in the graph (a fraction $`7/8`$ for $`d=3`$) are at most two-connected spins whose contribution dominates the entropy density. In a fully bonded graph ($`p=1`$), about half of these spins are free, the other half is constraint, with an energy density of $`\mathrm{ln}(2)/2`$. For small dilutions $`p1`$, about every removed bond turns a two-connected spin into a one-connected spin, which is always constrained, and the entropy declines. Eventually, at large dilution (small $`p`$), an increasing number of one-connected spins become totally disconnected and free again, raising to entropy density to its ultimate value, $`\mathrm{ln}(2)`$, for the unconnected ($`p=0`$) lattice. By this argument the entropy density should be given simply by

$`{\displaystyle \frac{S}{n}}\left[(1p)^2+{\displaystyle \frac{1}{2}}p^2\right]\mathrm{ln}(2),`$ (10)
which fits the data exceedingly well.

In contrast with energy and entropy, the overlap distribution is clearly sensitive to the transition from paramagnetic to glassy behavior at $`p^{}`$, as shown in Fig. 5. There, we plot only its moment, since the distribution itself quickly converges to a structureless $`\delta `$-peak as expected for a model that is replica symmetric (RS) throughout. In fact, we find that $`|q|/\sqrt{q^2}1`$ to within 1% for all $`p`$. (For one, a broad overlap distribution indicative of widely separated ground state configurations and RSB MPV  would be antithetical to the reducibility of the hierarchical lattice.) About at $`p^{}`$, the location of the peak becomes non-zero and splits into two symmetric peaks, leading to a non-zero moment $`|q|`$.

The moment experiences a surprising maximum very close to $`p=1`$, varying sharply for small dilutions. Again, the maximum may be explained in terms of the abundance of two-connected spins: Whatever the size of strongly correlated spin clusters may be in the fully connected lattice ($`p=1`$), the overwhelming conversion of two-connected to one-connected spins at small dilution can only increase the size of those clusters at its edge, since one-connected spins are always entrained. But more importantly, the elimination of two-connected spins, say, within a subgraph as shown in Fig. 1 will enhance the probability that the two root-spins of that subgraph become correlated, which would merge the clusters that those spins are part of. While the latter effect could suddenly increase the size of typical clusters by a factor, it is not clear whether the steep slope observed at $`p=1`$ is not an artifact of our approximation to the overlap (see Sec. IV). Further numerical studies suggest that

$`|q|_p|q|_{p=1}(1p)^{\frac{1}{5}}(p1),`$ (11)
as seen in the inset of Fig. 5. We have not found an explanation for the scaling in Eq. (11) in terms of the lattice yet, since it is also a property of the entrainment in the spin glass state (and hopefully, the overlap).

As shown in Fig. 2, for all $`p>p^{}`$ we obtain the expected scaling of the width of the defect energy distribution with system size $`L`$,

$`\sigma (\mathrm{\Delta }E)aL^y\left[1+ϵ(L)\right],`$ (12)
for large enough $`L`$, including a next-to-leading correction $`ϵ(L)1`$. Most conveniently, it appears that our $`k=1`$ algorithm (for large enough pools $`A_I`$) described in Sec. IV is least likely to create spurious dependencies for $`p=1`$, where we would expect the quickest and most stable convergence toward the renormalized bond distribution. Thus, we disabled the sampling of energy-offsets, entropies, and entrainments, which already are well-converged at $`I_{\mathrm{max}}=20`$, and only focused on the reduction of bonds. Setting $`I_{\mathrm{max}}=100`$ ($`d=3`$) or $`I_{\mathrm{max}}=40`$ ($`d=4`$) , $`k=1`$, and $`A_I=1024`$ for each $`I`$, we measured more than $`10^9`$ reduced bonds at each $`I`$ to obtain the distributions of defect energies $`\mathrm{\Delta }E=2J`$ in $`d=3`$ and 4. (In fact, we repeated the same calculation with a different random number generator, without markable difference in results. A previous run with $`I_{\mathrm{max}}=150`$ in $`d=3`$ showed divergence due to limited floating point precision at $`I120`$; the same problem emerged already at $`I_{\mathrm{max}}=40`$ in $`d=4`$.) At that level, we obtained a ratio between its absolute mean and its width $`\sigma (\mathrm{\Delta }E)`$ of $`<10^4`$ in both dimensions, even for the largest $`I`$, which is a measure of the statistical error (the distribution is supposed to be symmetric).

Surprisingly, although we extended the scaling regime in $`d=3`$ by a power of 10 compared to previous calculations ($`L=2^{100}`$ compared to $`L=2^{10}`$ DM ), the extrapolation for the scaling exponent $`y`$ in Eq. (12) only yields 5 digits accuracy, see top of Fig. 6. First, we approximated $`y0.2555`$ from the intercept of $`\mathrm{log}_2(\sigma )/\mathrm{log}_2(L)`$ with respect to $`1/\mathrm{log}_2(L)`$ at $`L\mathrm{}`$, which displays large corrections to scaling. We are then able to refine our estimate of $`y`$ by converting Eq. (12) to $`\mathrm{ln}(\sigma )y\mathrm{ln}(L)\mathrm{ln}(a)+ϵ(L)`$, using $`\mathrm{ln}(1+ϵ)ϵ`$. Plotting that expression again vs. $`1/\mathrm{log}_2(L)`$ (note that $`\mathrm{log}_2(L)=I`$) suggest a straight-line asymptotic behavior that would imply $`ϵb/\mathrm{log}_2(L)`$ for $`L\mathrm{}`$. The top of Fig. 6 indicates that a value of $`y=0.25546(3)`$ in $`d=3`$ leads to an apparently converging result with $`\mathrm{ln}(a)3/4`$ at the intercept for $`L\mathrm{}`$. We have used the same procedure to find $`y=0.76382(4)`$ in $`d=4`$, see bottom of Fig. 6. Although the floating point accuracy here deteriorated already at $`L=2^{35}`$, the asymptotic scaling regime is reached for much small $`L`$ than in $`d=3`$, leading to comparable accuracy in the fit.

Our value for the stiffness exponent $`y`$ in $`d=3`$ is 6% below $`y=[d1+\mathrm{log}_2(12/\pi )]/2=0.2698`$ obtained analytically in Ref. SY  with a Gaussian approximation, and consistent with $`y0.25`$ obtained numerically in Ref. BM  by evolving a “pool” of bonds through several generations, very similar to our algorithm in Sec. IV. Since both of these references use an initial Gaussian bond distribution, in contrast to our discrete distribution in Eq. (1), the results suggest that the stiffness exponent is independent of the bond distribution. In Ref. BKM , an expansion in the limit of large dimension $`d\mathrm{}`$ is used to demonstrate the independence of $`y`$ of the initial bond distribution. A recursion in the generation $`I`$ is constructed which for $`I\mathrm{}`$ leads invariably to the Gaussian derived in Ref. SY . With finite-dimensional corrections, this expansion improves on the Gaussian result above to yield $`y_30.255`$, in excellent agreement with our result here. In contrast, we expect that energy, entropy, overlap, and even $`p^{}`$ would depend on the distribution at least quantitatively.

## VI Conclusions

We have presented a new method to reduce spin glasses on arbitrary graphical structures. The reductions are most appropriate for sparse graphs near their percolation threshold, which may be either completely (as in this paper) or substantially reduced, facilitating further optimization using standard methods eo\_y ; BoPe4 . In either case, the reductions allow to extract exact ground state energies and entropies, and even provide a good approximation to collective properties such as the overlap distribution.

With the reduction method we have reproduced a number of ground-state properties of the Migdal-Kadanoff hierarchical lattice at any dilution, and discovered various new features, which can be mostly explained in terms of the lattice structure. The simplicity of this model further permitted us to determine the exponent $`y`$ in $`d=3`$ and 4 for the growth of the defect energy with system size to unprecedented accuracy. We have recently applied the reduction method to diluted Edwards-Anderson models for $`d=3`$ and 4, obtaining $`y_3=0.240(5)`$ and $`y_4=0.60(1)`$ for the respective stiffness exponents eo\_y . These results demonstrate that the Migdal-Kadanoff approximation in $`d=3`$ is even closer to the Edwards-Anderson value than had been previously expected BM ; F+H , although they are definitely distinct. The differing values in $`d=4`$ reflect the fact that Migdal-Kadanoff is a low-dimensional approximation.

The results obtained here and in Ref. eo\_y  suggest that the reduction method may be a useful tool in the exploration of Ising spin glass systems on any graph or lattice near percolation that may reveal the onset of non-trivial configurations of ground states.

## Acknowledgments

I would like to thank A. Percus and M. Paczuski for helpful discussions.


format: text

text:
# A Modified Wigner’s Inequality for Secure Quantum Key Distribution

Permanent address: \]Istituto Elettrotecnico Nazionale G. Ferraris, Strada delle Cacce 91, 10135 Torino (Italy)

## Abstract

In this report we discuss the insecurity with present implementations of the Ekert protocol for quantum-key distribution based on the Wigner Inequality. We propose a modified version of this inequality which guarantees safe quantum-key distribution.

preprint: wigner PRA 04 2.tex

Following the first proposal by Bennett and Brassard bennet&brassard  and the later introduction of the Ekert protocol invoking entangled states ekert , various systems of quantum key distribution (QKD) have been implemented and tested by groups around the world, turning QKD into the most advanced application in quantum information science.

QKD offers the possibility that two remote parties, conventionally called Alice and Bob, exchange a secret random key to implement a secure encryption-decryption algorithm, without meeting bennet&brassard ; ekert ; bb92 . QKD provides a significant advantage over the public-key cryptography because the security of the distributed key relies on the laws of quantum physics bennet&brassard ; ekert ; bb92 , i.e. the wave packet collapse prohibits gaining information from a quantum channel without disturbing it. Indeed, any attempt by a third party (Eve) to obtain information about the key is detected.

Two main goals underly the implementation of QKD schemes. One is to create and preserve authentic quantum channels against decoherence effects induced by any interaction with the environment gisinrevmod , eventually reducing or even destroying invulnerability of quantum channels against Eve’s attack. The other is to provide a true guarantee of absolute security against any possible eavesdropping attack i.e. the security is not simply based on technological feasibility.

Scientists are currently using one of two photon sources for practical QKD, faint lasers faint1 ; faint2 ; faint3 ; faint4 ; faint5 ; faint6  or spontaneous parametric down-conversion to generate entangled photon states (SPDC) ekertrarity ; sasha ; qk2 ; qk3 ; qk1 .

Both schemes have disadvantages, although Brassard et al. brassard  proved theoretically that QKD schemes based on entangled photons offer enhanced performance and security relative to schemes based on weak coherent pulses. In fact, entangled states allow a further test of security based on the completeness of quantum mechanics; in other words, the Ekert protocol exploits either the Clauser-Horne-Shimony-Holt (CHSH) ekert  or the Wigner inequality qk2  as an ultimate test of eavesdropping attack.

Wigner’s inequality was originally intended qk2  to provide an easier but as equally reliable eavesdropping check as the CHSH test in implementing the Ekert protocol. In this report we disclose the weakness of the Wigner inequality as a security test when Eve has the total control of quantum channels.

In Fig. 1 we depict a typical scheme for QKD by means of entangled photons. QKD performed by pure entangled states relies on the realization of two quantum-correlated optical channels yielding quasi single-photon polarization states. Alice’s measurement of the polarization of a photon of the pair, automatically projects Bob’s photon into a defined polarization state, and reversibility. Alice’s and Bob’s detection apparatuses are built using polarization analyzer systems (PAS), polarizing beam splitters (PBS), detectors ($`+_A,_A,+_B,_B`$), data storage systems (computers), and synchronization systems.

Consider, as an example, type II parametric down-conversion entangled states kwiat2 , where the output two-photon states are a quantum superposition of orthogonally polarized photons, i.e. in the singlet state

$$|\psi ^{}=\frac{1}{\sqrt{2}}\left(|H_A|V_B|V_A|H_B\right).$$

Without any loss of generality, a particular QKD procedure can be described by the transformations, $`\widehat{T}_z`$,with $`z=A,B`$, on the polarization state of a photon, e.g. a rotation of angles $`\alpha _A`$ and $`\alpha _B`$ respectively according to

$`\widehat{T}_z|H_z`$ $`=`$ $`\mathrm{cos}\alpha _z|H_z+\mathrm{sin}\alpha _z|V_z,`$
$`\widehat{T}_z|V_z`$ $`=`$ $`\mathrm{sin}\alpha _z|H_z\mathrm{cos}\alpha _z|V_z.`$
We assume ideal PBS in channel $`z`$ transmitting photons in the state $`|H_z`$ towards detectors $`+_z`$ and reflecting photons in the state $`|V_z`$ towards detector $`_z`$.

The joint probabilities of each possible pair of detectors firing due to the arrival of a photon pair with PAS angles of $`\alpha _A`$ and $`\alpha _B`$ are given by

$`p_{\alpha _A,\alpha _B}(+_A,+_B)`$ $`=`$ $`\left|H_B|H_A\left|\widehat{T}_A\widehat{T}_B\right|\psi ^{}\right|^2`$
$`=`$ $`{\displaystyle \frac{1}{2}}\mathrm{sin}^2(\alpha _A\alpha _B)`$
$`p_{\alpha _A,\alpha _B}(+_A,_B)`$ $`=`$ $`\left|V_B|H_A\left|\widehat{T}_A\widehat{T}_B\right|\psi ^{}\right|^2`$
$`=`$ $`{\displaystyle \frac{1}{2}}\mathrm{cos}^2(\alpha _A\alpha _B)`$
$`p_{\alpha _A,\alpha _B}(_A,+_B)`$ $`=`$ $`\left|H_B|V_A\left|\widehat{T}_A\widehat{T}_B\right|\psi ^{}\right|^2`$
$`=`$ $`{\displaystyle \frac{1}{2}}\mathrm{cos}^2(\alpha _A\alpha _B)`$
$`p_{\alpha _A,\alpha _B}(_A,_B)`$ $`=`$ $`\left|V_B|V_A\left|\widehat{T}_A\widehat{T}_B\right|\psi ^{}\right|^2`$
$`=`$ $`{\displaystyle \frac{1}{2}}\mathrm{sin}^2(\alpha _A\alpha _B).`$
Here, we consider the variant of the Ekert protocol based on the Wigner inequality originally proposed in ref. qk2 . For the analyzer settings the possible combined choices by Alice and Bob can be split into two groups: the first for key distribution and the second testing security. We assume Alice and Bob measure randomly among four combined analyzer settings. Alice’s possible choices are $`\alpha _A=(\alpha _{A,1}=\pi /6,\alpha _{A,2}=0)`$ and Bob’s possible choices are $`\alpha _B=(\alpha _{B,2}=0,\alpha _{B,3}=\pi /6).`$ The key is obtained from the subset of measurements corresponding to parallel PAS settings (i.e. $`\alpha _{A,i}=\alpha _{B,i}=0`$) and the firing detectors $`+_A_B`$ and $`_A+_B`$.

The Wigner inequality is based on the Wigner parameter

$`W=p_{\alpha _{A,1},\alpha _{B,2}}(+_A,+_B)+p_{\alpha _{A,2},\alpha _{B,3}}(+_A,+_B)`$
$`p_{\alpha _{A,1},\alpha _{B,3}}(+_A,+_B).`$ (1)
Note that for the maximally entangled states $`W=1/8`$ while for any local realistic theory $`W0`$.

To obtain the Wigner inequality $`W0`$ it is necessary to summarize the Wigner argument wigner . Two assumptions are stipulated in the proofs of the Wigner inequality: locality and realism. Locality means that Alice’s measurement results do not influence Bob’s results, and vice versa. Realism means that, given any physical property its value exists independently of its observation or measurement.

In the present case this is translated in terms of a classical probability distribution, $`𝒫(x_1,x_2;y_2,y_3)`$, where $`x_1`$ and $`x_2`$ are the hidden variables associated with the physical property inducing Alice’s measurement outcome in the presence of the PAS rotation angle $`\alpha _{A,i}`$. Similarly $`y_2`$ and $`y_3`$ correspond to the physical property inducing Bob’s outcomes when PAS angle is set as $`\alpha _{B,2}`$ and $`\alpha _{B,3}`$. Thus, we can identify the possible values of $`x_{1,2}`$ and $`y_{2,3}`$ with Alice and Bob’s measurement outcomes, in other words $`x_{1,2}=+_A,_A`$ and $`y_{2,3}=+_B,_B`$. Following the Wigner approach we write,

$`p_{\alpha _{A,1},\alpha _{B,3}}(+_A,+_B)={\displaystyle \underset{x_2,y_2}{}}𝒫(+_A,x_2;y_2,+_B)`$
$`=𝒫(+_A,+_A;+_B,+_B)+𝒫(+_A,_A;_B,+_B)`$
$`+𝒫(+_A,+_A;_B,+_B)+𝒫(+_A,_A;+_B,+_B).`$ (2)
In Wigner’s original paper wigner , $`\alpha _{A,2}=\alpha _{B,2}`$ and $`𝒫(x_1,+_A;+_B,y_3)=𝒫(x_1,_A;_B,y_3)=0`$. The assumption of perfect anticorrelation is obviously reasonable in the test of realism and locality of a physical theory, because it reflects the classical counterpart of a quantum system prepared in the singlet state, i.e. $`p_{\alpha _{A,2},\alpha _{B,2}}(+_A,+_B)=p_{\alpha _{A,2},\alpha _{B,2}}(_A,_B)=0`$. Thus, the inequality $`W0`$ is obtained from Eq. (2) by simply observing that $`𝒫(+_A,_A;+_B,+_B)_{x_2,y_3}𝒫(+_A,x_2;+_B,y_3)=p_{\alpha _{A,1},\alpha _{B,2}}(+_A,+_B)`$ and $`𝒫(+_A,+_A;_B,+_B)_{x_1,y_2}𝒫(x_1,+_A;y_2,+_B)=p_{\alpha _{A,2},\alpha _{B,3}}(+_A,+_B)`$.

Ref.s qk2 ; xuo  suggest to be cautious in applying Wigner’s inequality to test the security of cryptography schemes, since authors are aware that the Wigner’s inequality is derived assuming perfect anticorrelations, which are only approximately realized in practical situations. Nevertheless a punctual quantification of the insecurity induced by this assumption in the Ekert protocol based on Wigner’s inequality has never been pointed out.

In particular, after stating that the violation of the Wigner’s inequality ascertains the security of the quantum channels, Ref. qk2  suggests to replace the Wigner’s inequality by the generalization of the Bell’s inequality presented in Ref. ryff  for a substantial deviation from perfect anticorrelations. This last point has been overlooked by most readers probably because it is not clear what a ”substantial deviation from the perfect anticorrelations” is.

When the eavesdropper, Eve, measures photons on either one or both of Alice’s and Bob’s channels, her presence might be expected to be revealed by a higher value of $`W`$ than the local realistic theory limit, as it happens for the CHSH inequality ekert . Unfortunately this is not the case. Even if it can be proved that, for Eve adopting intercept-resend strategy and detecting only one photon of the pair, the limit becomes $`W_{\text{eve}}1/16`$, this is not for eavesdropping on both channels, because in this case there is no bound.

A proof of this last assertion, according to ref. ekert , can be obtained by considering Eve preparing each particle of the pairs separately, so that each individual particle has a well defined polarization direction, possibly varying from pair to pair. We write the probability $`P(\mathrm{\Phi }_A,\mathrm{\Phi }_B)`$ of having Alice’s particle in state $`|\mathrm{\Phi }_A`$, and Bob’s particle in state $`|\mathrm{\Phi }_B`$, where $`\mathrm{\Phi }_A`$ and $`\mathrm{\Phi }_B`$ are two angles of Eve’s polarization preparation i.e. $`|\mathrm{\Phi }_z=\mathrm{cos}\mathrm{\Phi }_z|H_z+\mathrm{sin}\mathrm{\Phi }_z|V_z`$. If the density operator associated with Eve’s pairs of photons is $`\widehat{\rho }_{\text{eve}}=P(\mathrm{\Phi }_A,\mathrm{\Phi }_B)|\mathrm{\Phi }_A|\mathrm{\Phi }_B\mathrm{\Phi }_B|\mathrm{\Phi }_A|`$d$`\mathrm{\Phi }_A`$d$`\mathrm{\Phi }_B,`$ then

$$W_{eve}=P(\mathrm{\Phi }_A,\mathrm{\Phi }_B)\left[\begin{array}{c}\mathrm{cos}^2(\mathrm{\Phi }_A+\frac{\pi }{6})\mathrm{cos}^2\mathrm{\Phi }_B+\mathrm{cos}^2\mathrm{\Phi }_A\mathrm{cos}^2(\mathrm{\Phi }_B\frac{\pi }{6})\\ \mathrm{cos}^2(\mathrm{\Phi }_A+\frac{\pi }{6})\mathrm{cos}^2(\mathrm{\Phi }_B\frac{\pi }{6})\end{array}\right]\text{d}\mathrm{\Phi }_A\text{d}\mathrm{\Phi }_B.$$
(3)
In the case of single-channel eavesdropping with intercept-resend strategy, the anticorrelation between the two photons in the Eve measurement base is preserved because $`\mathrm{\Phi }_B=\mathrm{\Phi }_A`$ $`\pi /2`$, and thus $`W_{eve}1/16.`$ In the most general case, when Eve has the total control over the state of individual particles, there is no physical bound, as $`W_{eve}`$ results in some cases below the limit of local realistic theory, i.e. $`W_{eve}<0`$ and surprisingly also below the quantum limit $`W_{eve}<1/8`$. For instance, if Eve gains total control of the source of photons, she can send photons in proper polarization states to avoid disclosure by a security test by Alice and Bob, and for $`P(\mathrm{\Phi }_A,\mathrm{\Phi }_B)=\delta (\mathrm{\Phi }_A0.6\pi )\delta (\mathrm{\Phi }_B0.4\pi )`$ the value $`W_{eve}=0.1995`$ is below the quantum limit.

In other words, the guarantee of Wigner’s security test against eavesdropping strategies is limited only to the detection of one photon of the pairs while the CHSH security is independent from the adopted strategy ekert .

Therefore, we end up with a modified Wigner’s parameter $`\stackrel{~}{W}`$ suitable for QKD security test by removing from Wigner’s original argument the anticorrelation assumption when $`\alpha _{A,2}=\alpha _{B,2}`$. Starting again from Eq. (2) but rejecting the assumption $`p_{\alpha _{A,2},\alpha _{B,2}}(+_A,+_B)=p_{\alpha _{A,2},\alpha _{B,2}}(_A,_B)=0`$ and observing that $`𝒫(+_A,+_A;+_B,+_B)+𝒫(+_A,_A;+_B,+_B)p_{\alpha _{A,1},\alpha _{B,2}}(+_A,+_B)`$, $`𝒫(+_A,+_A;_B,+_B)p_{\alpha _{A,2},\alpha _{B,3}}(+_A,+_B)`$ and $`𝒫(+_A,_A;_B,+_B)p_{\alpha _{A,2},\alpha _{B,2}}(_A,_B),`$ we obtain

$`\stackrel{~}{W}`$ $`=`$ $`p_{\alpha _{A,1},\alpha _{B,2}}(+_A,+_B)+p_{\alpha _{A,2},\alpha _{B,3}}(+_A,+_B)`$
$`+p_{\alpha _{A,2},\alpha _{B,2}}(_A,_B)p_{\alpha _{A,1},\alpha _{B,3}}(+_A,+_B)0.`$
for any local realistic theory.

We finally observe that for the singlet state we still obtain $`\stackrel{~}{W}=1/8`$ also in this case. This result happens because the modified parameter $`\stackrel{~}{W}`$ equals the original $`W`$ except for the additional term $`p_{\alpha _{A,2},\alpha _{B,2}}(_A,_B)`$, which is zero in the case of singlet state.

To demonstrate the robustness of this security test against an eavesdropping attack we perform a calculation analogous to Eq. (3)

$$\stackrel{~}{W}_{eve}=P(\mathrm{\Phi }_A,\mathrm{\Phi }_B)\left[\begin{array}{c}\mathrm{cos}^2(\mathrm{\Phi }_A+\frac{\pi }{6})\mathrm{cos}^2\mathrm{\Phi }_B+\mathrm{cos}^2\mathrm{\Phi }_A\mathrm{cos}^2(\mathrm{\Phi }_B\frac{\pi }{6})\\ +\mathrm{sin}^2(\mathrm{\Phi }_A)\mathrm{sin}^2(\mathrm{\Phi }_B)\mathrm{cos}^2(\mathrm{\Phi }_A+\frac{\pi }{6})\mathrm{cos}^2(\mathrm{\Phi }_B\frac{\pi }{6})\end{array}\right]\text{d}\mathrm{\Phi }_A\text{d}\mathrm{\Phi }_B.$$

Even in the general case when Eve has total control of the polarization states of photons in the two channels, we obtain that the minimum of $`\stackrel{~}{W}_{eve}=0.04428`$, well above the limit for local realistic theory. This result is completely equivalent to the one obtained with the security test based on the CHSH inequality ekert .

Some further analysis of $`\stackrel{~}{W}`$ must be considered for the practical implementation of the Ekert protocol based on Wigner’s inequality. According to qk2 , we highlight that the Ekert’s protocol based on modified Wigner’s inequality still guarantees a simplification with respect to the one based on the CHSH inequality, because Alice and Bob randomly choose between two bases rather than three. Though the necessity of an experimental evaluation of the term $`p_{\alpha _{A,2},\alpha _{B,2}}(_A,_B)`$ forces Alice and Bob to sacrifice part of the key for the sake of security. We note that in any practical implementation of QKD protocols, Alice and Bob distill from the noisy sifted key a nearly noise-free corrected key by means of error correction procedures subjected to the constraint of knowing the quantum bit error rate (QBER). Also, the QBER is estimated at the cost of losing part of the key. Thus, we suggest using the same sacrificed part of the key to estimate both $`\stackrel{~}{W}`$ and QBER.

To compare the performances of the Ekert protocol based on Wigner’s inequality versus the one CHSH-based ekert  we consider a protocol (Table I) where Alice and Bob measure randomly using three analyzer settings (as in the case of CHSH), where $`\alpha _{A,1}=\alpha _{B,1}=\pi /6`$, $`\alpha _{A,2}=\alpha _{B,2}=0`$ and $`\alpha _{A,3}=\alpha _{B,3}=\pi /6`$. $`\stackrel{~}{W}`$ and $`\stackrel{~}{W}^{}`$ correspond to two distinct test parameters and $`K`$ to the key distribution. This protocol is more efficient than the protocol based on CHSH. In terms of key generation we observe that in the case of CHSH only 2/9 of the qubits exchanged are devoted to the key generation ekert  while in the case of our protocol this quantity is a number between 2/9 and 1/3, depending on the security needs.

Even if, for some strict security request, all the qubits exchanged by the two parties with analyzer settings $`\alpha _{A,2}`$, $`\alpha _{B,2}`$ are devoted to the evaluation of $`\stackrel{~}{W}`$ ($`\stackrel{~}{W}^{}`$), and also for this protocol only 2/9 of the data still contribute to the key, in this protocol none of the qubits exchanged are discarded (Table I) while in the case of CHSH 1/3 of the qubits are discarded.

The results so far obtained allow us to quantify the so called ”substantial deviation from prefect anticorrelations” as that value of the QBER such that the quantum key distribution fails, because the key is insecure even if the Wigner’s inequality is violated.

For the original protocol introduced in qk2  the QBER of the key is the probability of correlated results when analyzer settings are $`\alpha _{A,2},\alpha _{B,2}`$, i.e.

$$\mathrm{QBER}=p_{\alpha _{A,2},\alpha _{B,2}}(_A,_B)+p_{\alpha _{A,2},\alpha _{B,2}}(+_A,+_B).$$

Thus, this leads to

$$0\stackrel{~}{W}W+p_{\alpha _{A,2},\alpha _{B,2}}(_A,_B)+p_{\alpha _{A,2},\alpha _{B,2}}(_A,_B),$$

which can be easily rewritten in terms of QBER as

$$\mathrm{QBER}+W0.$$
(5)
Thanks to the modified Wigner’s inequality a limit is found for both $`W`$ and the anticorrelation check. Thus, for the Ekert’s protocol based on the original Wigner’s inequality a secure key distribution is provided when, instead of the violation of the Wigner’s inequality ($`W0`$), the Wigner’s parameter $`W`$ satisfies $`W`$QBER. Eq. (5) and the related considerations obviously still hold for the enlarged protocol presented in Table I. Here one has to consider only the QBER of that part of the key obtained with analyzers settings $`\alpha _{A,2},\alpha _{B,2}`$, and not the QBER of the whole key.

As a final remark we highlight that our approach not only solves the problem of estimating the ”substantial deviation from the perfect anticorrelations” but completely overtakes it. Following our approach the experimentalist has simply to check the modified Wigner’s inequality from his experimental data to guarantee non-locality, i.e. security of the quantum channels as in the original Ekert’s idea ekert , without any further anticorrelation check.

Furthermore, the experimentalist in a noisy environment can discard a secure key when Eq. (5) is used instead of the modified Wigner’s inequality, since $`\stackrel{~}{W}W+`$QBER.

In conclusion, this paper discusses the security of Ekert’s protocol based on the Wigner inequality. We showed that the QKD Ekert protocol based on Wigner’s inequality presents a serious lack of security against some eavesdropping strategies other than intercept-resend. We emphasized the motivation beneath the missing security and propose a modified test ultimately guaranteeing secure QKD.

We would like to thank M. Rasetti and I. Ruo Berchera, G. Di Giuseppe, A. M Colla, F. Bovino, P. Varisco for useful discussion and helpful suggestions. This work was developed in collaboration with Elsag S.p.A., Genova (Italy), within a project entitled ”Quantum Cryptographic Key Distribution” co-funded by the Italian Ministry of Education, University and Research (MIUR) - grant n. 67679/ L. 488. In addition S. C. acknowledges the partial support of the DARPA QuIST program and M. L. R. acknowledges the partial support by INFM.


text:
# Slow dynamics of a confined supercooled binary mixture: direct space analysis

## I  Introduction

Liquids under confinement represent a field of growing interest in science due to the connection with relevant technological and biophysical problems grenoble . In this field modifications of both the phase diagram and dynamics of the confined liquid with respect to the bulk represent a key point, in particular as far as the possibility of supercooling is concerned pablo . Experiments exhibit a very diversified phenomenology for supercooled confined liquids. In particular the calorimetric glass transition temperature, $`T_g`$, can increase melni , decrease pissis ; kremer2  or be unaffected richert  with respect to the the bulk value depending on the liquid, the confining geometry and the nature of the substrate.

Phenomenological arguments predict the existence of domains of cooperative dynamics in the supercooled liquid adam-gibbs . The size of these domains is expected to grow upon supercooling. The existence of an upper bound for the domain size in confined media implies a decrease of $`T_g`$ upon confinement. In fact if this size remains small on supercooling then the probability of a cooperative motion, and therefore of a structural relaxation, is large, and dynamics is faster.

The mode coupling theory of the glass transition (MCT) goetze  is able to describe the dynamics of bulk liquids in the supercooled region on approaching a crossover temperature $`T_C`$. At $`T_C`$ the system passes from a regime where ergodicity is attained through structural relaxations to a regime where this mechanism is frozen and only activated processes permit the exploration of the configurational space. $`T_C`$ can be estimated through experiments and computer simulations or predicted by the theory.

Above $`T_C`$ the relaxation mechanism of the supercooled liquid can be described as mastered by the cage effect. Nearest neighbors surround and trap the tagged particle forming a cage around it. When the cage relaxes, due to cooperative motions, the particle moves. MCT translates this phenomenological description into a precise mathematical framework where, starting from the classical equation of motion for the density correlator $`\varphi _q(t)`$, a retarded memory function is introduced. In the idealized version of MCT hopping processes are neglected and the retarded non linear set of integro-differential equations can be solved analytically to the leading order in $`ϵ=(TT_C)/T_C`$, deriving universal results for the density correlator. With these approximations $`T_C`$ is the temperature of structural arrest of the ideal system. The success of this theory is due to the fact that on approaching $`T_C`$ from above hopping processes can be neglected in many liquids and the predictions of the idealized version of MCT hold.

Molecular dynamics studies of model liquids in restricted geometries intended to assess the applicability of MCT can give an important contribution to the characterization of the glass transition scenario in confinement. In this field recent progresses have been done for water confined in a silica nano pore gallo1 ; gallo2  and for confined thin polymer films baschnagel .

We performed a Molecular Dynamic (MD) study of the single particle dynamical properties of a Lennard-Jones binary mixture (LJBM) embedded in an off lattice matrix of soft spheres. The model considered represents a situation of strong confinement analogous to real silica xerogels rosinberg ; monson ; page . The bulk phase of the chosen mixture is known to test most of the MCT predicted features kob . A brief report on some of our findings has been recently published euro-noi .

We present in this paper direct space quantities evaluated from MD trajectories in order to define the range of validity of MCT with respect to the bulk phase. In the next section we describe the model adopted and the details of the MD simulations. In the third section we illustrate the radial pair distribution functions and the coordination numbers of the system. The fourth section is devoted to the study of the dynamic quantities in the direct space: the mean square displacement (MSD), the van Hove self correlation function and the non-Gaussian parameter. Last section contains the concluding remarks.

## II Model and Molecular dynamics

The confining medium consists in a rigid disordered array of soft spheres. The simulation box contains 16 soft spheres labeled in the following with the letter M. The interspaces of this structure host a liquid LJBM of 1000 particles, composed of 800 particles of type A and 200 particles of type B. The interaction pair potential can be written in the following general form:

$$V_{\mu \nu }(r)=4ϵ_{\mu \nu }\left[\left(\frac{\sigma _{\mu \nu }}{r}\right)^{12}\eta _{\mu \nu }\left(\frac{\sigma _{\mu \nu }}{r}\right)^6\right]$$
(1)
where indexes $`\mu ,\nu `$ run on the particle types A,B,M. The parameters for the three components are listed in Table 1. The parameters of the LJBM have been chosen as in ref. kob . The shifted potential technique has been used with a cutoff of $`2.5`$ $`\sigma _{\mu \nu }`$. Periodic boundary conditions are applied. In the following Lennard-Jones units will be used, therefore energy will be given in units of $`ϵ_{AA}`$, temperature in units of $`ϵ_{AA}/K_B`$, length in units $`\sigma _{AA}`$ and time in units of $`(m\sigma _{AA}^2/(48ϵ_{AA}))^{1/2}`$. The box length is $`L=12.6`$. The simulation box of the system is represented in Fig. 1. As it can be seen from the picture the LJBM is in a strongly confined environment as only few layers of particles can accommodate among the soft spheres. Fig. 2 shows a cell made of four simulation boxes. Only equipotential surfaces of the soft spheres are displayed in order to offer a clearer image of the disordered porous structure hosting the liquid and to highlight the analogies with silica xerogels.

We have conducted MD simulations of this system in the NVE ensemble. The equations of motion were solved by the velocity Verlet algorithm. The system was equilibrated at different reduced temperatures via a velocity rescaling procedure. The temperatures investigated are $`T=5.0,2.0,0.8,0.58,0.538,0.48,0.465,0.43,0.41,0.39`$. The timestep used for $`T1`$ was $`0.01`$ and for $`T<1`$ was $`0.02`$. For the lowest temperature investigated, $`T=0.37`$, a production run of $`14\times 10^6`$ timesteps was performed. We verified that the results obtained do not depend on the specific choice of the disordered matrix by running MD simulations for two other different configurations of the disordered matrix of soft spheres. Both thermodynamics and dynamics appeared the same for the three systems, in particular deviations among the diffusion coefficients extracted at given temperatures from the MSD of different configurations are within $`2`$%.

Values of thermodynamics quantities of the equilibrated systems are shown in Fig 3. Averaged pressure, potential energy and total energy per particle of the mixture behave similar to the bulk. A stability of the system throughout all the isochore explored is evident from the smoothness of the curves.

Due to the soft spheres potential the free volume accessible for the mixture strongly depends on temperature. In Fig. 4 we show an estimate of the free volume for A and B particles in our confined system as a function of temperature. The estimate has been carried out through the Voronoi tessellation. This procedure has been applied only to non-interfacial particles. The free volume has been calculated by averaging over $`160`$ configurations of the system for each temperature.

## III Radial pair correlation functions

The radial pair correlation functions $`g_{\mu \nu }(r)`$ are shown in Fig.5. For all the temperatures investigated the system is in the liquid phase as no abrupt sharpening of the peaks is evident. No indications of phase separation appear upon supercooling either. The latter phenomenon would in fact result in a substantial decrease of the area of the first peak, $`ab_1`$, of the $`g_{AB}(r)`$, which we do not observe in our system. On lowering temperature the peaks sharpen and enhance more than in the bulk due to the strong packing induced by the soft spheres. A clear element of distinction with respect to the bulk is the enhancement of intensity of the first $`g_{BB}(r)`$ peak, $`bb_1`$, as supercooling progresses. The same peak appeared instead reduced to a small shoulder for supercooled temperatures in the bulk.

The functions $`g_{MA}(r)`$ and $`g_{MB}(r)`$ are characteristic of this confined system and are representative of the interaction between the confining soft spheres and the particles of the mixture. We note the gradual disappearance of the first peak, $`mb_1`$, of the $`g_{MB}(r)`$ that begins at $`T=2.0`$ while the second peak, $`mb_2`$, enlarges its area. Correspondingly the first peak of the $`g_{MA}(r)`$ shifts to larger $`r`$. For lower temperatures, namely $`T<0.8`$, the structure of the mixture stabilizes so that MA and MB peaks alternate in space.

In order to quantify the above considerations for the radial distribution functions we calculated the coordination number defined as:

$$N_{\mu \nu }(r_1,r_2)=4\pi \rho x_\nu _{r_1}^{r_2}r^2g_{\mu \nu }(r)𝑑r$$
(2)
where $`\rho `$ is the total number density $`N/V`$, where $`N=1016`$, $`V=L^3`$ and $`x_\nu `$ is the relative fraction of the $`\nu `$ species. The quantity of Eq. 2, evaluated for each peak of the $`g_{\mu \nu }(r)`$, is shown in Fig. 6. The number of B particles around the soft sphere in the first shell lowers upon decreasing temperature in favor of an enhancement in the second shell. Correspondingly the number of A particles in the first shell enhances.

The above considerations on the radial distribution functions and the coordination numbers lead to the conclusion that B particles avoid the interface with the soft spheres for low temperatures. In Fig. 7 a snapshot of the system for $`T=0.37`$ evidences this behavior. The fluctuations of the error bars in Fig.4 also confirm this phenomenon. In fact since interfacial particles have been excluded from the calculation of the free volume a larger fluctuation corresponds to a frequent exchange of the particles with the interface.

## IV Dynamic quantities

In the following we shall describe the temperature variations of the time dependent direct space quantities analyzed for our confined system upon supercooling, namely the MSD, the van Hove self correlation function and the non-Gaussian parameter.

The asymptotic behaviour of the MSD can be used as a test of an MCT prediction. The MSD is in fact expected, as supercooling progresses, to develop a plateau at intermediate time scale, corresponding to the rattling of the particle in the cage. After the initial short time ballistic motion the particle enters the plateau whose extension in time enhances upon lowering temperature. After leaving the plateau region the Brownian behaviour predicted for a simple liquid in a stable state is restored. From the slope of the late MSD the diffusion constant $`D`$ can be extracted through Einstein relation $`<r^2(t)>=6Dt`$. The temperature behaviour of diffusion coefficients can give an estimate of the mode coupling crossover temperature through the MCT power law behaviour

$$D(TT_C)^\gamma $$
(3)
In Fig.8 we show the MSD calculated for both A and B particles as a function of time for all the temperatures investigated. It can be seen that also in this confined system the binary mixture displays the typical features of a glass former. From the height of the plateau the cage radius can be extracted and the value obtained for both A and B particles is $`r_c^2=0.04`$, analogous to the bulk result kob .

The lowest temperature reached in our study is $`T=0.370`$ against $`T=0.466`$ for the bulk. The elongation of the plateau of the MSD for our lowest temperature is comparable to that of the lowest temperature of the bulk. The diffusion coefficients extracted from the slopes are reported in Fig. 9. At higher temperatures $`T>1.0`$ the larger A particles diffuse slower than B particles while they move faster at lower temperatures. As shown in the previous section, in the confined case for lower temperatures B particles tend to avoid the interface of the soft spheres. As a consequence they remain confined in the inner part of the channels, surrounded by A particles, and this diminishes the accessible free volume lowering mobility below the one of A particles, see also Fig.4.

In Fig. 9 the fit to the MCT power law of Eq. 3 is also reported. The estimate of the crossover temperature given by the fit is $`T_C=0.356`$ for both A and B particles. The exponents are $`\gamma =1.86`$ and $`1.89`$ respectively for A and B particles. We do observe deviation from a power law both at high and low temperatures. The temperatures used for the fit are in the range $`0.410T0.58`$ corresponding to $`0.153<ϵ<0.631`$. The range of validity of MCT reported for the bulk is much wider, $`0.07<ϵ<1.30`$. The upper limit is generally due in bulk to the non-validity of the asymptotic expansion for too large values of $`ϵ`$. Here we infer that a further decrease of this limit is due to the progressive B particles avoidance of the soft sphere interfaces upon supercooling. For the temperatures below the upper limit the peaks of g(r) for MA and MB have in fact already stabilized, see Fig.5. The deviation observed at the two lowest temperatures is to be connected to the presence of hopping mechanisms, analogously to bulk liquids, as we will see in more detail when we illustrate the behaviour of the VHSCF.

The Van Hove self correlation function, VHSCF, is defined, for a system of N particles, as

$$G_S(r,t)=\frac{1}{N}\underset{i=1}{\overset{N}{}}\delta [𝐫+𝐫_i(0)𝐫_i(t)]$$
(4)
$`4\pi r^2G_S(r,t)dr`$ is the probability of finding a particle at distance $`r`$ after a time $`t`$ if the same particle was at the origin $`r=0`$ at the initial time $`t=0`$.

In Fig.10 we show $`4\pi r^2G_S(r,t)dr`$ at fixed times for three different temperatures. On this scale the correlators of A and B particles do not differ substantially. At high temperature the VHSCF decays regularly as time increases as expected for a simple liquid. As supercooling progresses the system begins to display a new feature that is best evident for the lowest temperature. This feature consists in the clustering of the curves at intermediate times which is the signature of the caging of the particles. The zone where the peak value clusters is $`r=0.150.25`$ and corresponds to the zone in which the MSD flattens. For the lowest temperature investigated the maximum and minimum times for the clustering are $`2.6t82.0`$. These values define the so called $`\beta `$ relaxation temporal region. In the comparison between the VHSCF of A and B particles it can be observed the crossing between the velocity of the two species. At high temperatures, in fact, the VHSCF of B particles extends further than the one of A particles for long times, while for low temperatures behaviors are swapped.

A behavior not present in the bulk can be observed for the lower temperatures by looking at the peak of the VHSCF for times longer than the plateau region. In this temporal range the peak of the VHSCF is expected, in the framework of MCT, to shift to larger $`r`$ for large times. The absence of this shift in our confined mixture, as we shall see in the next picture, must be connected with an important hopping phenomenon able to restore diffusive motions when the cage is not jet dissolved.

In Fig.11 we show a blow up of the tails of the curves of Fig.10. Aside the previously cited cage peak, we note the presence of a hopping peak around $`r=1.0`$. This peak enhances for late time and it is more pronounced for A particles than for B. In the bulk on approaching $`T_C`$ from above a hopping peak appears only for B particles. From the MSD it is evident that our hopping phenomena can not be neglected at the two lowest temperatures, namely at $`T=0.37,0.39`$. As a consequence these temperatures can not be accounted for a test of the MCT in the region around the cage relaxation time, $`\alpha `$ relaxation region. Nevertheless these temperatures can be used to test the theory in the $`\beta `$-time region. In fact, see also ref.kob , hopping peaks do not affect curves in the time range of the clustering region.

For the A particle a second broader peak is visible in Fig.11. This peak is also connected to the presence of hopping. The two hopping peaks correspond to the first and second peak of the $`g_{AA}(r)`$, see Fig.6, meaning that there are some privileged position for the atoms to be reached and that those atoms that reached these positions have traveled more than the average. Similar, the hopping peak of the B VHSCF corresponds to the first peak of the $`g_{BB}(r)`$. Hopping phenomena are connected with the existence of dynamical heterogeneities in supercooled liquids schober . In confinement they appear to have above $`T_C`$ a more important role than in the bulk.

Dynamical heterogeneities are also known to cause a non gaussianity of atomic displacements. We evaluated the first non-Gaussian parameter (NGP), that works as a quantifier of the degree of non-Gaussianity. This parameter is defined yip  in three dimension as

$$\alpha _2(t)=\frac{3<r^4(t)>}{5<r^2(t)>^2}1$$
(5)
This quantity calculated for our system is reported in fig.12. Also for this parameter we have marked differences with respect to the bulk that reflect the different weight of the hopping processes in the two cases. In the bulk the NGP of A and B particles were different and never exceeded the values of 2 for A particles (no hopping above $`T_C`$) and 3.2 for B particles. Here instead we have similar NGP for A and B particles and the value grows enormously as we supercool. For the last temperature investigated the peak is around 5. We found no relation between the peak positions and MCT power law of Eq.3 as it was found for example in bulk water fs  and the curves of the confined NGP do not appear to follow a master curve before the maximum as found in bulk liquids kob ; fs . The trend of the NGP of bulk and confined mixture depicts a scenario in which MCT behaviour causes a mild deviation from gaussianity. The peak is around 2 when close to $`T_C`$ kob ; fs , its position shifts substantially upon supercooling and the NGP follows a master curve in the temporal region before the peaks. Hopping appears to cause a marked deviation from gaussianity instead, leading to a disappearance both of the MCT master curve and of the shift of the peaks upon supercooling.

## V Summary and Conclusions

We presented an MD study of a LJBM confined in a disordered array of soft spheres. The direct space quantities have been analyzed upon supercooling for a test of the MCT behaviour.

Two main differences with respect to the bulk mixture appear due to confinement. The $`g(r)`$ have evidenced that smaller B particles tend to avoid interfaces at low T and correspondingly their diffusion coefficient becomes lower than that of the A particles. Hopping processes are present also for A particle above $`T_C`$ as shown in the VHSCF.

Not obviously, these differences do not prevent the MCT asymptotic predictions from holding for this mixture also in this kind of confinement. The range of validity of MCT in the region of the $`\alpha `$ relaxation is however much more limited. The range of the bulk kob  is $`0.07<ϵ<1.30`$ against the $`0.153<ϵ<0.631`$ found in our confined model. Disturbances in the dynamical behaviour of the bulk introduced by a different type of confinement are also reported to lead to a complete disappearance of MCT kobconf .

Power law fit to the diffusion coefficient extracted from the slope of the MSD give an estimate of $`T_C=0.356`$ for both species and $`\gamma =1.86`$ for A particles and $`\gamma =1.89`$ for B particles in agreement with MCT that states that these quantities should be the same for both A and B particles. The values of the bulk extracted from D in literature are kob  $`T_C=0.435`$ for both species and $`\gamma =2.0`$ for A particles and $`\gamma =1.7`$ for B particles. The exponents $`\gamma `$ appear similar while we have a reduction of $`T_C`$ upon confinement. A complete test of the correlators in the Q space will be presented in a subsequent paper noi .

Acceleration of dynamics in presence of repulsive confinement is expected baschnagel  based on the fact that in confinement the cooperative rearranging region cannot grow beyond the size of the system. In the case of an attractive interaction this effect might superpose to the often severe slowing down caused by the substrate to the closest layers of the liquidgallo1 ; gallo2 .

In the present work we have considered a much more complex system than the bulk. Therefore the fact the MCT could be used in this context as an unifying theoretical approach is highly relevant as a guideline for the systematic study of the important phenomenology of confined and interfacial liquids.

## VI Acknowledgments

P.G. wishes to thank J. Baschnagel,W. Götze and F. Varnik for very useful discussions.