Heuristic for estimation of multiqubit genuine multipartite entanglement

# Heuristic for estimation of multiqubit genuine multipartite entanglement

Paulo E. M. F. Mendonça111Permanent address: Academia da Força Aérea, C.P. 970, 13.643-970 Pirassununga, SP, Brazil Marcelo A. Marchiolli Gerard J. Milburn ARC Centre for Engineered Quantum Systems, School of Mathematics and Physics, The University of Queensland, St. Lucia, Queensland 4072, Australia Avenida General Osório 414, centro, 14.870-100 Jaboticabal, SP, Brazil
July 12, 2019
###### Abstract

For every -qubit density matrix written in the computational basis, an associated “X-density matrix” can be obtained by vanishing all entries out of the main- and anti-diagonals. It is very simple to compute the genuine multipartite (GM) concurrence of this associated -qubit X-state, which, moreover, lower bounds the GM-concurrence of the original (non-X) state. In this paper, we rely on these facts to introduce and benchmark a heuristic for estimating the GM-concurrence of an arbitrary multiqubit mixed state. By explicitly considering two classes of mixed states, we illustrate that our estimates are usually very close to the standard lower bound on the GM-concurrence, being significantly easier to compute. In addition, while evaluating the performance of our proposed heuristic, we provide the first characterization of GM-entanglement in the steady states of the driven Dicke model at zero temperature.

###### keywords:
Multipartite Entanglement, Genuine Multipartite Concurrence, Multiqubit X-states
###### Pacs:
03.65.Ud, 03.67.Mn, 03.65.Aa
journal: Annals of Physics\biboptions

sort&compress

## 1 Introduction

During the last one and a half decades, increasing interest has been manifested in the topic of multipartite entanglement. Just as we learnt, in the early days of quantum information, that bipartite entanglement is an invaluable resource for quantum cryptography 91Ekert661 (), communication 92Bennett2881 () and speed-up of classical algorithms 92Deutsch553 (); 94Simon116 (); 96Grover212 (); 97Shor1484 () (to name but a few), we have now a growing awareness of the role that multipartite entanglement plays in quantum computing 03Jozsa2011 (); 11Bruss052313 (); 01Raussendorf5188 (); 09Briegel19 (), high-precision metrology 04Giovannetti1330 (); 12Hyllus022321 (); 13Marchiolli1330001 (); 14Lucke155304 (); 14Vitagliano032307 (), quantum phase transitions 02Osborne032110 (); 05Bruss014301 (); 06Oliveira010305R (); 06Oliveira039902E1 (); 07Oliveira039901E2 (); 10Montakhab062313 (); 13Giampaolo052305 (); 14Hofmann134101 (); 14Giampaolo93033 (); 14Stasinska032330 () and even biology 10Sarovar462 (); 10Caruso062346 ().

Despite significant interest and progresses, many limitations from the theory of bipartite entanglement are naturally inherited in the multipartite setting. For example, since the entanglement detection problem was shown to be NP-hard already in the bipartite case 03Gurvits10 (); 07Ioannou335 (), there is very little hope that computable measures of (multipartite) entanglement can be found for an arbitrary quantum state. As a result, much of the effort in this field is centered around devising sufficient criteria for (multipartite) entanglement detection 96Peres1413 (); 96Horodecki1 (); 00Terhal319 (); 10Huber210501 (), which are then turned into computable (multipartite) entanglement estimators 02Vidal032314 (); 05Brandao022310 (); 11Jungnitsch190502 (); 11Ma062325 ().

Amongst all types of multipartite entanglement, special interest is devoted to the detection and quantification of genuine multipartite entanglement, a type of multipartite entanglement established collectively between all parties of an -partite system. More precisely, a GM-entangled state is any state that cannot be written as a mixture of biseparable states, which, in turn, are those states that are separable with respect to some bipartition of the relevant Hilbert space 09Guhne1 (); 09Horodecki865 (). Besides, it is GM-entanglement that occupies the highest position in the hierarchy of multipartite entanglement 13Levi150402 ().

In this paper, we consider a particular set of sufficient conditions for detecting GM-entanglement introduced by Huber et al. 10Huber210501 () — and turned into an entanglement estimator by Ma et al. 11Ma062325 () — to obtain our main result: a heuristic approach for approximating these estimates with a considerably lower computational cost than that required by the optimization problem appearing in Refs. 10Huber210501 (); 11Ma062325 (). In practice, this result opens the way for systematical and analytical estimation of GM-entanglement in symmetric multiqubit mixed states and, more generally, enables a significant reduction in the computational time required by standard numerical solvers, usually accompanied by some small impact on the accuracy of the estimate.

Our heuristic relies upon three basic facts: (i) the estimates proposed by Ma et al. 11Ma062325 () are actually lower bounds on the GM-concurrence — a GM-entanglement monotone introduced by the authors to generalize the Wootters concurrence 98Wootters2245 () to the multipartite setting; (ii) this lower bound is saturated for the family of -qubit X-states, in which case it is also remarkably easy to compute 12Rafsanjani062303 (); and (iii) in general, it is computationally easier to determine a local-unitary (LU) transformation that approximates a generic density matrix to the X-form, than it is to solve the optimization problem in Refs. 10Huber210501 (); 11Ma062325 (). Accordingly, our so-called X-heuristic, estimates the GM-concurrence of an arbitrary -qubit density matrix by applying a LU-transformation that minimizes the entries out of its main- and anti-diagonals and, subsequently, outputting the GM-concurrence of the associated X-density matrix (obtained by neglecting any entry that may have remained out of the main- and anti-diagonals of the LU-transformed density matrix).

In order to illustrate the strengths and drawbacks of the X-heuristic, we explicitly apply it to the family of diagonal symmetric states 54Dicke99 (); 14Quesada052319 (); 14Wolfe140402 () and steady states of the driven Dicke model 79Puri200 (); 80Drummond1179 (); 81Lawande4171 (); 02Schneider042107 (). In both cases, we discuss the trade-off between computational efficiency and accuracy vis-à-vis the standard scheme for GM-concurrence estimation 10Huber210501 (); 11Ma062325 (). Our comparative analysis reveals the X-heuristic as a useful computational tool for GM-concurrence estimation in multiqubit systems.

As side results, while assessing the X-heuristic with the steady states of the driven Dicke model at zero temperature, we provide the first characterization of multipartite entanglement in this model. In particular, we note that our estimates agree on that GM-entanglement is maximal for the parameter values corresponding roughly to a bifurcation of a fixed point in the corresponding semiclassical dynamics, reinforcing the long-standing (and widely-accepted) conjecture that multipartite entanglement must be maximal at quantum phase transitions 02Osborne032110 (). While the driven Dicke model is a zero-dimensional many-body problem, the bifurcation is controlled by a parameter in the Hamiltonian and the correlations in the quantum steady state change at a particular value of that parameter in analogy to the change in the character of a ground state in a spatially extended many-body quantum phase transition. Moreover, we show that the two-qubit steady states of this model are examples of non-X states whose concurrence can be exactly computed with the simple concurrence formula for X-states, settling an open problem posed in Ref 12Rafsanjani12043912 ().

Our paper is structured as follows. In Sec. 2 we review some key concepts relating GM-entanglement and its quantification. In Sec. 3 our heuristic for GM-concurrence estimation is presented and, in Sec. 4, its accuracy is benchmarked against the standard scheme for the two aforemetioned families of -qubit mixed states. In Sec. 5 we focus on the computational advantage of our heuristic over the standard scheme. Finally, in Sec. 6, we summarize our results and discuss some possible avenues for future work.

## 2 Preliminaries

In this section we briefly review the concept of GM-entanglement and some key results regarding its detection and quantification. In particular, we outline the scheme introduced in Refs. 10Huber210501 (); 11Ma062325 () (here referred to as the “optimal-” method), a current standard for GM-concurrence estimation in multipartite mixed states. Although we are only interested in multiqubit states, throughout most of this section we leave the dimensionality of each party arbitrary.

### 2.1 Genuine Multipartite Entanglement

An -partite pure state (with ), is said to be GM-entangled if it is not biseparable; is said to be biseparable, in turn, if there is a bipartition of and a pair of states and such that

 |ψ⟩=|ψA⟩⊗|ψB⟩. (1)

In general, there are possible bipartitions of and, if a decomposition of the form (1) is found with respect to any one of them, then can be immediately declared biseparable. On the other hand, can only be declared GM-entangled after ruling out the existence of such a decomposition with respect to every bipartition of .

The notions of GM-entanglement and biseparability are extended to mixed states via a convex roof construction: An -partite mixed state acting on is GM-entangled if it is not biseparable; is biseparable if it can be expressed as a convex combination of biseparable pure states, that is,

 ρ=∑ipi|ψi⟩⟨ψi| (2)

where every is biseparable. Remarkably, although the biseparability of requires every in Eq. (2) to be biseparable, each one of them can be so with respect to a different bipartition of , meaning that a biseparable mixed state need not to be separable with respect to any particular bipartition of .

So far we have only established the notion of GM-entanglement. In what follows we review some recent results on (i) how to tell if a given density matrix is biseparable (detection) and (ii) if not, how to estimate the amount of GM-entanglement that it contains (quantification). Needless to say, these are hard problems even in the case  03Gurvits10 (); 07Ioannou335 (); 09Horodecki865 (), let alone .

### 2.2 GM-entanglement detection

Detection of GM-entanglement has been intensely studied (see, e.g., Ref. 09Guhne1 () for a review). To date, one of the most effective detection schemes was proposed by Huber et al. 10Huber210501 (), where it was shown that every biseparable state satisfies every inequality of the -parametrized family

 I|Φ⟩(ρ)≤0, (3)

where

 I|Φ⟩(ρ):=√⟨Φ|ρ⊗2Π|Φ⟩−n∑i=1√⟨Φ|(ΠAi⊗1Bi)†ρ⊗2(ΠAi⊗1Bi)|Φ⟩. (4)

In the above, can be chosen as any product state of , i.e.,

 |Φ⟩=|μ1μ2μ3⋯μN⟩⊗|ν1ν2ν3⋯νN⟩. (5)

The symbol denotes the global permutation operator that peforms simultaneous permutations of all subsystems across the two copies of ,

 Π(|μ1μ2μ3⋯μN⟩⊗|ν1ν2ν3⋯νN⟩)=|ν1ν2ν3⋯νN⟩⊗|μ1μ2μ3⋯μN⟩, (6)

whereas only permutes those subsystems whose labels are in (of a given Hilbert space bipartition ). For example, consider the following bipartitions of the -partite Hilbert space :

 {A1|B1}={1|2,3,…,N},{A2|B2}={2|1,3,…,N}and{AN+1|BN+1}={1,2|3,…,N}. (7)

Then,

 (ΠA1⊗1B1)|μ1μ2μ3⋯μN⟩⊗|ν1ν2ν3⋯νN⟩ =|ν1μ2μ3⋯μN⟩⊗|μ1ν2ν3⋯νN⟩, (8) (ΠA2⊗1B2)|μ1μ2μ3⋯μN⟩⊗|ν1ν2ν3⋯νN⟩ =|μ1ν2μ3⋯μN⟩⊗|ν1μ2ν3⋯νN⟩, (9) (ΠAN+1⊗1BN+1)|μ1μ2μ3⋯μN⟩⊗|ν1ν2ν3⋯νN⟩ =|ν1ν2μ3⋯μN⟩⊗|μ1μ2ν3⋯νN⟩. (10)

Thus, if one can find a state such that , then can be promptly declared to be GM-entangled. Unfortunately, such a criterion does not detect every GM-entanglement [there are GM-entangled density matrices that satisfy the entire family of inequalities (3)], but it is stronger than many commonly used criteria (see, e.g., 08Seevinck032101 () and references therein).

### 2.3 GM-concurrence estimation: The optimal-|Φ⟩ scheme

The relevance of transcends its application in the detection problem, manifesting itself also in the context of GM-entaglement quantification. Indeed, it was shown by Ma et al. 11Ma062325 () that, for every choice of , the following inequality holds

 CGME(ρ)≥max[0,2I|Φ⟩(ρ)], (11)

where denotes the concurrence of GM-entanglement of ; a well-defined (but difficult to compute) entanglement monotone.

Unlike , the lower bounds are easy to measure and to compute, emerging thus as natural candidates for GM-concurrence estimation. However, the accuracy of such estimates are strongly dependent on the choice of for a given , highlighting the importance of devising efficient ways of making such a choice.

Formally, the strongest inequality of the family (11) is given by

 CGME(ρ)≥max[0,2max|Φ⟩I|Φ⟩(ρ)]≡C|Φ⟩(ρ), (12)

where the optimization runs over all product states of . Since any product state can be constructed from LU transformations on a reference product state , we can set

 |Φ⟩=(N⨂n=1Vn⊗N⨂n=1Wn)|Φ0⟩, (13)

where and are elements. In this case, the optimization is performed over independent (special) unitary transformations ( and ), in terms of which the objective function takes the form

 I|Φ⟩(ρ)=|⟨0|¯¯¯¯U†0ρU0|0⟩|−n∑i=1√⟨0|U†iρUi|0⟩⟨0|¯¯¯¯U†iρ¯¯¯¯Ui|0⟩, (14)

where we have chosen and defined the LU transformations

 U0:=N⨂n=1Vn,¯U0:=N⨂n=1Wn,Ui:=N⨂n=1Si,n,¯¯¯¯Ui:=N⨂n=1¯¯¯¯Si,n. (15)

In the above, the single particle unitary matrices and have the index referring to a particular bipartition of and, in terms of this, are given by:

 Si,n:={Wnif% n∈AiVnotherwiseand¯¯¯¯Si,n:={Vnif n∈AiWnotherwise. (16)

Thus, we end up with the optimization problem

 maxVn,Wn∈SU(dn)|⟨0|¯¯¯¯U†0ρU0|0⟩|−n∑i=1√⟨0|U†iρUi|0⟩⟨0|¯¯¯¯U†iρ¯¯¯¯Ui|0⟩. (17)

Let us now restrict the analysis of this optimization problem to the case of -qubits (i.e., for every ). Since single-qubit unitary transformation require only real parameters and our objective function involves unitaries, the number of real variables to be optimized is . Although the linear scaling with is certainly an appealing feature of this approach, the resulting objective function is not everywhere differentiable. This fact compromises the performance of numerical algorithms that rely on differentiating the objective function to approach the optima.

## 3 X-Heuristic

We now introduce a heuristic approach for estimating the GM-concurrence of a density matrix of -qubits. To motivate our heuristic, we start by reviewing the concept of X-density matrices and how to quantify their GM-concurrence. An -qubit -density matrix is any density matrix that, in the computational basis, takes the form

 (18)

with and, for every , , , (non-negativity of ) and (normalization of ). As shown by Hashemi Rafsanjani et al. 12Rafsanjani062303 (), such states can have their GM-concurrence easily — and exactly — quantified by the formula:

 CX(ρX)=max[0,cX(ρX)] (19)

where

 cX(ρX):=2maxk∈Sn⎡⎣rk−∑j∈Sn∖{k}√ajbj⎤⎦%andSn:={1,…,n}. (20)

The proposed heuristic takes advantage of the simplicity of Eq. (19) to produce computable estimates of GM-concurrence for non-X density matrices. In fact, as noted in Ref. 12Rafsanjani062303 () (and explictly demonstrated in the A), the evaluation of for an arbitrary -qubit density matrix (written in the computational basis), leads to a lower bound on the GM-concurrence of . In principle, depending on the particular form of , this lower bound can be very loose and, consequently, a very poor estimate of [see, e.g., Eqs. (44) and (45)]. However, since this lower bound is saturated for X-density matrices, one should expect that the “closest” is to a X-density matrix, the better that estimate will turn out to be.

The X-heuristic materializes this reasoning by “approximating” to a X-matrix with a LU transformation and, subsequently, estimating its GM-concurrence as , where

 ~ρ:=UlρU†l. (21)

Of course, the constraint of LU transformation is sufficient to imply that and have the same amount of GM-concurrence. However, whereas is invariant under LU transformations of , is not.

Once again, we are left with the problem of choosing a LU transformation, but now aiming to implement the “approximation to the X-form”. In order to formalize the meaning of this approximation, we introduce the penalty function

 f(~ρ)=∑j>ij+i≠2N+1|~ρi,j|2, (22)

which (i) vanishes if and only if is a X-density matrix and (ii) smoothly increases as the magnitudes of the off-X entries of increase. Under this “-metric”222Here and throughout, the term “metric” is used in a loose way to suggest a notion of distance between a given density matrix and the set of X-density matrices., the optimal choice of is given by the solution of the optimization problem

 minUlf(UlρU†l), (23)

which is a non-linear least square problem on real variables with residuals given by . We can think of it as finding the LU related density matrix that most accurately “fits” the X-density matrix model.

Besides halving the number of variables in the optimization problem (17), the objective function in (23) is smooth, which implies that numerical algorithms can take advantage of well-defined derivatives to more efficiently approach a minimum. On the analytical side, for density matrices with enough symmetry, the objective function will be sufficiently simple to enable educated guesses on the location of the minima, which can be certified a posteriori by evaluating the gradient and the spectrum of the Hessian matrix at the candidate point. In any case, even if the gradient vanishes and the Hessian turns out to be positive definite, we cannot be sure that we have attained a global minimum. Nevertheless, at this point we may have already obtained a non-trivial lower bound on the amount of GM-concurrence of at a considerably low computational cost.

## 4 Applications

In this section we apply the X-heuristic to estimate the amount of GM-concurrence in two families of -qubit density matrices: the diagonal symmetric states and the steady states of the driven Dicke Model at zero temperature. For each family, we benchmark the quality of the resulting estimates against those produced by the optimal- scheme.

### 4.1 Diagonal symmetric states

Diagonal symmetric states are a natural extension of W-states to the domain of multiqubit mixed states and have recently attracted much attention in the fields of GM-entanglement detection and quantification 09Toth170503 (); 12Tura060302R (); 13Novo012305 (); 14Quesada052319 (); 14Wolfe140402 (). Formally, they can be defined as

 ρds,N:=N∑k=0pk|DNk⟩⟨DNk|, (24)

where and . Here, the states are the (totally symmetric) -qubit Dicke states of excitations 54Dicke99 (); 03Stockton022112 (), defined in the computational basis as

 |DNk⟩:=1√CNk∑σ|σ(1,\lx@stackrelk⋯,1,0,\lx@stackrelN−k⋯,0)⟩ (25)

with the summation running over every distinct permutation of the sequence of ones and zeros. Of course, there are such permutations, justifying the normalization constant upfront. Moreover, note that are the -qubit W-states 00Dur062314 (), which sets the context in which can be considered a generalization thereof. In the B, we give explicit matrix forms of in the computational basis for .

For , the resulting density matrix is already an X-state [cf. Eq. (46)]. As a result, , and . For , the corresponding density matrices are no longer of the X-form [see, e.g., Eqs. (47) and (48)] and thus, in order to optimally estimate with , we must first determine a LU transformation

 Ul=N⨂j=1[cosϑjsinϑjeiφj−sinϑje−iφjcosϑj] (26)

that minimizes . Thanks to the symmetry of , some simple analysis of this function reveals a minimum at333A complementary numerical analysis strongly suggests that Eq. (27) gives, actually, a global minimum.

 ϑj=π4andφj=0for % everyj∈SN, (27)

for which the corresponding values of can be promptly computed with the aid of Eqs. (21) and (22). In Fig. 1, these are plotted for a family of populations specified by the single parameter as follows:

 (28)

Intuitively, the plots of Fig. 1 can be interpreted as some error measure associated with the process of estimating as . Indeed, if (e.g., case ), then is exactly equal to . As increases, deviates from the X-form and, hence, the estimate is expected to be poorer. Nevertheless, recall that is always a lower bound on , so the “-error bars” are not centered in , but lay strictly above it.

In Fig. 2 we present the resulting values of (solid line) along with the points corresponding to for one hundred values of uniformly distributed in the range . While the line was analytically constructed, each point was obtained by numerically solving the optimization problem (17). Overall, we notice that both estimates follow very similar trends and coincide in a significant portion of the domain. Next, their similarities and differences are discussed in greater detail.

For , owing to the X-form of , is precisely equal to the Wootters concurrence and, therefore, exactly matches for every — see Plot 2a. For , this exact matching no longer occurs for () and (). Nevertheless, as the insets in Plots 2b,c show, the difference between and is orders of magnitude smaller than the actual value of the estimates, being thus negligible for most practical applications. For odd , a more significant mismatch between and occurs at the extremes of the domain (cf. Plots 2d,e,f), but it rapidly decreases and vanishes in a large central region of the domain, where neither estimates detect GM-entanglement. An important difference between and becomes apparent in the case (Plot 2f), where we note that each estimate arrives to (and departs from) zero at slightly different values of , in such a way that there are a few states whose GM-entanglement can be detected by the optimal- scheme, but not by the X-heuristic.

### 4.2 Steady states of the driven Dicke model

From the theoretical point of view, the Dicke model 54Dicke99 () represents a special quantum mechanical model whose unique mathematical and physical properties allow us to describe important cooperative phenomena (such as, for example, resonance fluorescence and cooperative emission) far from the thermodynamic limit. Their virtues are also reflected in works involving the exact steady state solutions of its master equation that include (or not) the effects of detuning between a collective driving field and the atomic resonant frequency 79Puri200 (); 80Drummond1179 (); 81Lawande4171 (); 02Schneider042107 (). In this section, we focus on the steady state solution of the driven Dicke model at zero temperature

 ρs,N=1DNN∑m,n=0(J−g∗)m(J+g)n, (29)

where depicts the number of two-level atoms (ions), is a normalization constant, denotes the ratio between the Einstein -coefficient of each atom (ion) and the Rabi frequency , and . Note that for correspond to the collective raising and lowering operators expressed in terms of the Pauli matrices which satisfy, together with the collective inversion operator , the usual angular momentum commutation relations.

In the B, explicit matrix forms of in the computational basis are given for . Noticeably, the resulting matrices are not of the X-form and, as shown in Fig. 3, they cannot be exactly brought to the X-form via LU transformations (i.e., they do not assume the X-form in any orthonormal basis of product states). These observations qualify as a good testbed for the X-heuristic and so, in Fig. 4, we present the corresponding estimates obtained from numerical implementations of the X-heuristic and optimal- scheme.

At a glance of Fig. 4, both and are seen to follow very similar trends: a zero plateau that extends from up until , followed by a sudden growth and, finally, an asymptotic decay. In what follows we offer a more in-depth analysis of these results by considering the cases and separately.

In the case , on top of and , we have also plotted the Wootters concurrence of which, in this case, is the “genuine multipartite” concurrence . Despite not being of the X-form, we note that the three measures coincide up to the numerical precision of , with the zero plateu ranging up until and the maximal concurrence () occurring at . Interestingly, the observed coincidence between the three measures signals that the X-concurrence formula may also be an exact concurrence formula for some non-X states. In A, this point is further explored with a demonstration that matches exactly the Wootters concurrence () of . This observation provides a (constructive) negative answer to an open question in Ref. 12Rafsanjani12043912 () on whether the saturation of requires to be of the X-form in some orthonormal product basis.

For , and no longer coincide for every value of . Nevertheless, the difference between them asymptoticaly decreases and nearly vanishes at relatively small values of . For example, the ratio is smaller than for (), (), (), () and (). Of course, this behavior was to be expected since for such values of we have already seen, in Fig. 3, that is nearly X-formed. In fact, thanks to that, we can also infer that must well aproximate the (unknown) value of in such values of .

As it ocurred with the diagonal symmetric states, the main drawback in replacing with is that there are certain states that have its GM-entangled detected by the former but not by the latter. In Fig. 4, this is expressed by the fact that the zero plateau generated by the X-heuristic is longer than that obtained with the optimal- scheme. However, as the plots show, the difference is usually very small, meaning that failure to detect GM-entanglement with the X-heuristic will only occur for very few and specific values of . Other than that, it is also worth noticing that each method foresees a maximal value of at slightly different values of , in such a way that . However, the difference is again very small and, ultimately, it is not clear which one more accurately describes the actual location of the maximum of .

## 5 Computational Efficiency

The GM-concurrence estimates produced by the X-heuristic cannot improve on the corresponding estimates produced by the optimal- scheme; a fact that follows directly from the derivation of the X-heuristic (cf. Sec. 3). There is, however, an important trade-off between accuracy and efficiency that must be considered before disregarding the X-heuristic in favor of the optimal- scheme. On the accuracy side, we have already seen that the X-estimates are usually pretty good approximations of the optimal- estimates. On the efficiency side, we now explicitly demonstrate that the X-estimates are significantly easier to compute both for the diagonal symmetric states and for the steady states of the driven Dicke model.

As mentioned before, the computational advantage of the X-heuristic over the optimal- scheme relies upon two aspects of the objective function associated with the former: (i) it is smooth and (ii) it depends on a number of variables that halves the number of variables in the (non-smooth) objective function of the optimal- scheme. While the relevance of (ii) is obvious, the importance of (i) cannot be overestimated. Thanks to (i), the X-heuristic optimization problem can be regarded as ideal for the application of a quasi-Newton method 06Nocedal (), a family of nonlinear optimization algorithms that take advantage of well-defined derivatives of smooth objective functions to more efficiently converge to a minimum.

In this work, the numerical computations of the X-heuristic and optimal- estimates were performed with the MATLAB function fminunc, which implements a popular quasi-Newton algorithm known as BFGS (after its discoverers Broyden, Fletcher, Goldfarb and Shanno). As with any quasi-Newton algorithm, the BFGS performs successive evaluations of the gradient to build a quadratic model of the objective function that is sufficiently good to attain a superlinear rate of convergence. It contrasts with Newton’s method 06Nocedal () in the sense that it does not require (nor attempts to compute) the Hessian matrix, typically a time-consuming and error-prone task. Instead, at each step, it gains information about the second derivative along the search direction by considering changes in the gradient. Although the BFGS algorithm requires repeated computations of the gradient, it has been noted to perform well also in nonsmooth optimization problems, as long as it does not run into a nonsmooth point (see, e.g., 13Lewis135 () and references therein). Consistently with this observation, we have observed a better performance of the optimal- scheme with the BFGS algorithm than with algorithms that do not rely on evaluations of the gradient (e.g., Nelder Mead algorithm 65Nelder308 (), implemented by the MATLAB function fminsearch), especially for larger values of . For this reason, aiming to build the most efficient implementation of the two schemes, we employed the BFGS algorithm for both.

In order to compare the performance of the two schemes in producing reasonable estimates of GM-concurrence, we set a threshold equal to the X-estimates presented in Figs. 2 and 4 and timed how long it took for each scheme to reach that threshold. As it is usually the case with search algorithms, a starting value had to be provided, and so we resorted to random initial guesses (uniformly sampled) to avoid biasing the optimizer toward (or against) a satisfactory minimum. By doing so, we could gauge how hard it is for each scheme to reach the threshold from a zero-knowledge initial condition.

Each optimization was repeated times444Except for the time measurements of the optimal- scheme for states with , which were repeated only times due to the long time-frame necessary to converge to a satisfactory value., and the five-number summary of the resulting “wall-clock time” distributions is shown in Fig. 5 as a box-and-whiskers plot. Plot 5a refers to the time measurements obtained for the diagonal symmetric state with for every , whereas Plot 5b refers to the time measurements for the steady states of the driven Dicke model that, for each , have the value of corresponding to the peak of the X-estimates in Fig. 4. We note that to obtain each box-and-whisker appearing in Fig. 5, several initial guesses were typically made because either (i) the optimizer converged to a unsatisfactory minimum or (ii) the optimizer failed to converge after iterations. Of course, the resulting box-and-whiskers take into account every (if any) unsuccessful attempt as well as the successful one.

Remarkably, Fig. 5 reveals a significant advantage of the X-heuristic over the optimal- scheme, especially for the larger values of . For , for example, the X-heuristic estimates were produced in roughly  s, whereas the optimal- scheme was unable to reach the threshold once, even after a week (approximately  s) of computation. Most importantly, the plots suggest that the computational time for each scheme scales differently with , with the X-heuristic resembling an exponential growth and the optimal- some worse scaling. Naturally, this observation cannot be explained by the difference in the number of variables associated to each scheme (since they both scale linearly with ), and we assign it to the smoothness (or lack thereof) of the objective functions. Because of this scaling, the X-heuristic may be the only viable option for GM-concurrence estimation in generic -qubit states with intermediate values of .

The results of this section indicate that, despite the theoretical superiority of the optimal- scheme over the X-heuristic, the latter may actually overperform the former in practical situations where a time horizon has to be considered. For example, if we allow one second to obtain a GM-concurrence estimate for a state with , then Fig. 5 implies that we are much more likely to obtain a better estimate with the X-heuristic than with the optimal- scheme, which may have converged to a lesser estimate (or even have not converged at all). Naturally, time constraints like this become virtually more important as increases.

## 6 Concluding Remarks

In this paper we introduced a heuristic method for estimating the GM-concurrence of an -qubit density matrix. It consists of evaluating the -qubit X-state GM-concurrence formula for density matrices that are not of the X-form, but were previously brought “as close as possible” to it with a LU transformation. We have shown that the estimates thus produced are lower bounds on the GM-concurrence by demonstrating that they actually lower bound a standard lower bound on the GM-concurrence (attainable via a numerical procedure which we refer to as the “optimal-” method). Most importantly, by examining two prominent families of mixed -qubit states, we found that our estimates are usually very close to those produced by the optimal- scheme and are significantly easier to be numerically computed and, in certain cases of high symmetry, can be analytically obtained.

As for future directions, let us mention two possible extensions of the present work. First, it would be interesting to evaluate the performance of the X-heuristic in other well-established families of -qubit mixed states such as reduced ground states of spin chain models. As a matter of fact, it is fair to say that the X-heuristic has already been successfully applied in this context, in Ref. 14Giampaolo93033 (), where the genuine tripartite concurrence of the reduced ground state of three spins symmetrically distributed in the cluster-Ising model was exactly calculated thanks to the fact that such states can be turned into X-states via LU transformations. Nevertheless, the application of our heuristic in alike models for which the X-form is not exactly attainable is yet to be explored and is likely to provide an efficient tool for GM-entanglement estimation in condensed matter systems.

Second, the requirement of approximating the X-form via LU transformations is unnecessarily strong for the goal of preserving GM-concurrence. In fact, it was recently conjectured that a generic two-qubit density matrix can always be turned into a X-density matrix of same entanglement via a unitary transformation that is not necessarily local 13Hedemann (). In Ref. 14Mendonca79 (), this conjecture was rigourously proved for three entanglement measures and, moreover, a semi-analitic prescription for constructing the corresponding “entanglement-preserving-unitaries” was delineated. Along the same lines, it would be interesting to characterize the most general family of “GM-concurrence-preserving-unitaries” and attempt to approximate the X-form via conjugation with elements of this less constrained family. In that case, better X-form approximations should be expected, presumably improving the quality of the resulting GM-concurrence estimates. Of course, this new optimization problem may be less suitable for numerical implementation, in which case the trade-off between efficiency and accuracy should be carefully reconsidered.

## Acknowledgements

P.E.M.F.M. acknowledges the financial support of the Brazilian Air Force and the program “Ciência sem Fronteiras”, Project No. 200024/2014-0.

## Appendix A X-concurrence is a lower bound for GM-concurrence

In this appendix we demonstrate that the X-concurrence formula (19) applied to a generic -qubit density matrix is equal to for a certain family of states (recall that ), being thus a lower bound on the GM-concurrence.

Consider the following (decimal) representation for the -qubit computational basis

 |0⟩≡|00⋯00⟩,|1⟩≡|00⋯01⟩,…,|N⟩≡|10⋯00⟩, |N+1⟩≡|10⋯01⟩,…,|2N−1⟩≡|11⋯11⟩, (30)

and, in terms of this, define

 |Φμ⟩:=|μ⟩⊗|2N−1−μ⟩for everyμ∈{0,…,n−1}. (31)

Substitution of this into Eq. (4) gives, after some straightforward algebra,

 I|Φμ⟩(ρ)=|⟨μ|ρ|2N−1−μ⟩|−∑ν∈{0,…,n−1}∖{μ}√⟨ν|ρ|ν⟩⟨2N−1−ν|ρ|2N−1−ν⟩. (32)

Now, let the density matrix of be decomposed (in the computational basis) as

 (33)

where is some Hermitian matrix with zeros along the main- and anti-diagonals. In terms of this, Eq. (32) clearly takes the form

 I|Φμ⟩(ρ)=rμ+1−∑ν∈{0,…,n−1}∖{μ}√aν+1bν+1. (34)

Making and (in such a way that ), it is clear that

 I|Φk−1⟩(ρ)=rk−∑j∈Sn∖{k}√ajbj (35)

and, thus [cf. Eq. (20)]

 cX(ρ)=2maxk∈SnI|Φk−1⟩(ρ). (36)

This fact, combined with inequality (12) and identity (19), leads to the desired result:

 CGME(ρ)≥max[0,2max|Φ⟩I|Φ⟩(ρ)]≥max[0,2maxk∈SnI|Φk−1⟩(ρ)]=max[0,cX(ρ)]=CX(ρ). (37)

As suggested in Ref. 12Rafsanjani12043912 () (within the framework of ), it is interesting to look for density matrices that saturate the above inequalities. Of course, saturation occurs if is an X-state, other than that (and to the best of our knowledge), no examples of density matrices saturating inequality (37) have been identified so far. In what follows, we introduce a family of two-qubit density matrices which are not of the X-form in any product state basis and, nonetheless, satisfy , where denotes the Wootters concurrence (i.e., the “genuine multipartite” concurrence of two-qubit states).

Consider the two-qubit steady states of the driven Dicke model at zero temperature, whose density matrix in the computational basis is presented in Eq. (53). Noticeably, it is not an X-state for and, as Plot 3a shows, there is no LU transformation capable of simulatenously vanishing all of its non-X entries (although this can be remarkably well aproximated for sufficiently large values of ). Nevertheless, as the following computation shows, the Wootters concurrence exactly matches the X-concurrence in this case.

We first apply the X-concurrence formula to Eq. (53), which yields

 CX(ρs,2)=2max[0,γ2−1D2,γ2−√1+2γ2+4γ4D2]=max[0,2(γ2−1)D2]. (38)

Then, in order to compute the Wootters concurrence of , we start by evaluating

 (39)

whose eigenvalues are

 {1+2γ4+2γ2√1+γ4D22,1D22,1D22,1+2γ4−2γ2√1+γ4D22}. (40)

Since the first element in the above set is clearly the largest eigenvalue, one has

 CW(ρs,2) =1D2max[0,√1+2γ4+2γ2√1+γ4−2−√1+2γ4−2γ2√1+γ4] (41) =max[0,2(γ2−1)D2] (42)

where, in order to obtain the second row, we used the denesting identities

 √1+2γ4±2γ2√1+γ4=√1+γ4±γ2. (43)

Comparing Eqs. (38) and (42), it is clear that .

Unfortunately, saturation of (37) does not occur for the zero-temperature steady states of the driven Dicke model with . For , for example, density matrices in the computational basis are explicitly presented in Eqs. (54) and (55), respectively. Although Plots 4b,c reveal that these states are GM-entangled for , a straightforward computation of theirs X-concurrence yields

 CX(ρs,3)=2D3max[0,6γ3−3√Γ1(Γ2+4γ4),2γ3−2√Γ1(Γ2+4γ4)−√Γ3(1+12γ4)]=0 (44)

and

 (45)

which shows that is only a trivial lower bound for in these cases.

## Appendix B Density Matrices in the Computational Basis

In this appendix we present the density matrices, written in the computational basis, of the -qubit diagonal symmetric states and steady states of the driven Dicke model, for . For ease of visualuzation, we restrict to show the lower triangular entries, replace zeros with dots and boldify the anti-diagonal entries.

### b.1 Diagonal Symmetric States

 ρds,2 =12⎡⎢ ⎢ ⎢ ⎢⎣2p0⋅p1⋅p1p1⋅⋅⋅2p2⎤⎥ ⎥ ⎥ ⎥⎦ (46)
 ρds,3 =13⎡⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢⎣3p0⋅p1⋅p1p1⋅⋅⋅p2⋅p1p1⋅p1⋅⋅⋅p2⋅p2⋅⋅⋅p2⋅p2p2⋅⋅⋅⋅⋅⋅⋅3p3⎤⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥⎦ (47)
 ρds,4 =14⎡⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢⎣4p0⋅p1⋅p1p1⋅⋅⋅p2⋅p1p1⋅p1⋅⋅⋅p2⋅p2⋅⋅⋅p2⋅p2p2⋅⋅⋅⋅⋅⋅⋅p3⋅p1p1⋅p1⋅⋅⋅p1⋅⋅⋅p2⋅p2p2⋅⋅p2⋅⋅⋅p2⋅p2p2⋅⋅p2p2⋅⋅⋅⋅⋅⋅⋅p3⋅⋅⋅p3⋅⋅⋅p2⋅p2p2⋅⋅p2p2⋅p2⋅⋅⋅⋅⋅⋅⋅p3⋅⋅⋅p3⋅p3⋅⋅⋅⋅⋅⋅⋅p3⋅⋅⋅p3⋅p3p3⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅⋅4p4⎤⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥⎦ (48)

### b.2 Steady states of the driven Dicke model

For brevity, we define the following family (parametrized by ) of quadratic functions of :

 Γk:=1+kγ2. (49)

The normalization appearing as a multiplicative factor in front of the matrices forms are given below:

 D2 :=4(1+γ2+γ4), (50) D3 :=4(2+3γ2+6γ4+9γ6), (51) D4 :=16(1+2γ2+6γ4+18γ6+36γ8). (52)
 ρs,2 =1D2⎡⎢ ⎢ ⎢ ⎢⎣1iγ1+γ2iγγ21+γ2−2γ2iγ(1+2γ2)iγ(1+2γ2)1+2γ2+4γ4⎤⎥ ⎥ ⎥ ⎥⎦ (53)