IsSeparable

From QETLAB
Jump to: navigation, search
IsSeparable
Determines whether or not a bipartite operator is separable

Other toolboxes required cvx
Related functions InSeparableBall
IsBlockPositive
IsPPT
IsPSD
OperatorSchmidtDecomposition
OperatorSchmidtRank
SymmetricExtension
SymmetricInnerExtension
Function category Entanglement and separability

IsSeparable is a function that determines whether a bipartite operator is separable or entangled. A value of 0 indicates that the operator is entangled, a value of 1 indicates that the operator is separable, and a value of -1 indicates that the script was unable to determine whether the operator is separable or entangled.

Syntax

  • SEP = IsSeparable(X)
  • SEP = IsSeparable(X,DIM)
  • SEP = IsSeparable(X,DIM,STR)
  • SEP = IsSeparable(X,DIM,STR,VERBOSE)
  • SEP = IsSeparable(X,DIM,STR,VERBOSE,TOL)

Argument descriptions

  • X: A bipartite positive semidefinite operator.
  • DIM (optional, by default has both subsystems of equal dimension): A 1-by-2 vector containing the dimensions of the two subsystems that X acts on.
  • STR (optional, default 2): An integer that determines how hard the script should work to determine separability before giving up (STR = -1 means that the script won't stop working until it finds an answer). Other valid values are 0, 1, 2, 3, .... In practice, if STR >= 4 then most computers will run out of memory and/or the sun will explode before computation completes.
  • VERBOSE (optional, default 1): A flag (either 0 or 1) that indicates that the script will or will not display a line of text describing how it proved that X is or is not separable.
  • TOL (optional, default eps^(3/8)): The numerical tolerance used throughout the script.

Examples

Determining entanglement of a bound entangled state

The following code constructs a two-qutrit bound entangled state based on the "Tiles" unextendible product basis. This state's entanglement can not be detected by the positive partial transpose criterion, but the following code shows that it is indeed entangled.

>> v = UPB('Tiles');
 
>> rho = eye(9);
>> for j = 1:5
       rho = rho - v(:,j)*v(:,j)';
   end
>> rho = rho/4; % we are now done constructing the bound entangled state
 
>> IsSeparable(rho)
Determined to be entangled via the realignment criterion. Reference:
K. Chen and L.-A. Wu. A matrix realignment method for recognizing entanglement.
Quantum Inf. Comput., 3:193-202, 2003.
 
ans =
 
     0

The following code performs the same computation, but has VERBOSE set to 0 so that the method of proving that rho is entangled is not displayed.

>> IsSeparable(rho,[3,3],2,0)
 
ans =
 
     0

Mixtures with the maximally-mixed state

Every state that is sufficiently close to the maximally-mixed state is separable, so for every entangled state $\rho \in M_A \otimes M_B$, the state $\sigma_p := p\rho + \tfrac{1}{d_A d_B}(1-p)I$ is separable as long as $p$ is small enough.

For the remainder of this example, $\rho$ is the same bound entangled state based on the "Tiles" UPB from the previous example. It was shown in [1] that the Filter Covariance Matrix Criterion detects the entanglement in $\sigma_p$ when $p = 0.8723$, which we can verify as follows:

>> p = 0.8723; sigma = p*rho + (1-p)*eye(9)/9;
>> IsSeparable(sigma)
Determined to be entangled via the Filter Covariance Matrix Criterion. Reference:
O. Gittsovich, O. Gühne, P. Hyllus, and J. Eisert. Unifying several separability
conditions using the covariance matrix criterion. Phys. Rev. A, 78:052319, 2008.
 
ans =
 
     0

However, if we decrease $p$ to $p = 0.8722$ then a stronger test of entanglement is needed to determine that $\sigma_p$ is entangled:

>> p = 0.8722; sigma = p*rho + (1-p)*eye(9)/9;
>> IsSeparable(sigma)
Determined to be entangled by not having a 2-copy symmetric extension. Reference:
A. C. Doherty, P. A. Parrilo, and F. M. Spedalieri. A complete family of separability
criteria. Phys. Rev. A, 69:022308, 2004.
 
ans =
 
     0

If we decrease $p$ further, then $\sigma_p$ becomes separable:

>> p = 0.2; sigma = p*rho + (1-p)*eye(9)/9;
>> IsSeparable(sigma)
Determined to be separable by closeness to the maximally mixed state. Reference:
L. Gurvits and H. Barnum. Largest separable balls around the maximally mixed bipartite
quantum state. Phys. Rev. A, 66:062311, 2002.
 
ans =
 
     1

If we increase $p$ back to $p = 0.4$ then $\sigma_p$ is still separable, but a stronger test is required to prove its separability:

>> p = 0.4; sigma = p*rho + (1-p)*eye(9)/9;
>> IsSeparable(sigma)
Determined to be separable via the semidefinite program based on 2-copy symmetric extensions from reference:
M. Navascues, M. Owari, and M. B. Plenio. A complete criterion for separability detection. Phys. Rev. Lett.,
103:160404, 2009.
 
ans =
 
     1

If we increase $p$ slightly more to $p = 0.45$ then $\sigma_p$ is still separable, but none of the default tests for separability are strong enough to determine this. To prove separability, we can increase the STR argument to 3:

>> p = 0.45; sigma = p*rho + (1-p)*eye(9)/9;
>> IsSeparable(sigma)
 
ans =
 
    -1
 
>> IsSeparable(sigma,[3,3],3)
Determined to be separable via the semidefinite program based on 3-copy symmetric extensions from reference:
M. Navascues, M. Owari, and M. B. Plenio. A complete criterion for separability detection. Phys. Rev. Lett., 103:160404, 2009.
 
ans =
 
     1

Isotropic states

It is well-known that isotropic states $(1-\alpha)I/d^2 + \alpha|\psi_+\rangle\langle\psi_+| \in M_d \otimes M_d$ are separable if and only if $\alpha \leq 1/(d+1)$,[2] which can be verified in the $d = 3$ case as follows:

>> d = 3;
>> IsSeparable(IsotropicState(d, 1/(d+1)))
Determined to be separable by being a small rank-1 perturbation of the maximally-mixed state. Reference:
G. Vidal and R. Tarrach. Robustness of entanglement. Phys. Rev. A, 59:141-155, 1999.
 
ans =
 
     1
 
>> IsSeparable(IsotropicState(d, 1/(d+1)+0.0001))
Determined to be entangled via the realignment criterion. Reference:
K. Chen and L.-A. Wu. A matrix realignment method for recognizing entanglement. Quantum Inf. Comput., 3:193-202, 2003.
 
ans =
 
     0

Werner states

It is well-known that Werner states $\frac{1}{d^2-d\alpha}(I\otimes I - \alpha S) \in M_d \otimes M_d$ are separable if and only if $\alpha \leq 1/d$,[3] which can be verified in the $d = 3$ case as follows:

>> d = 3;
>> IsSeparable(WernerState(d, 1/d))
Determined to be separable by closeness to the maximally mixed state. Reference:
L. Gurvits and H. Barnum. Largest separable balls around the maximally mixed bipartite quantum state. Phys. Rev. A, 66:062311, 2002.
 
ans =
 
     1
 
>> IsSeparable(WernerState(d, 1/d+0.0001))
Determined to be entangled via the PPT criterion. Reference:
A. Peres. Separability criterion for density matrices. Phys. Rev. Lett., 77:1413-1415, 1996.
 
ans =
 
     0

Notes

In general, we are much better at proving that an operator is entangled than we are at proving that it is separable. Thus, if this script returns a value of -1 (indicating that it was unable to prove that the operator is separable or entangled) then that is a pretty good indicator that the operator is probably separable.

Source code

Click on "expand" to the right to view the MATLAB source code for this function.

  1. %%  ISSEPARABLE    Determines whether or not a bipartite operator is separable
  2. %   This function has one required argument:
  3. %     X: a positive semidefinite matrix
  4. %
  5. %   SEP = IsSeparable(X) is either -1, 0, or 1. A value of 0 indicates that
  6. %   X is entangled, a value of 1 indicates that X is separable, and a value
  7. %   of -1 indicates that the separability of X could not be determined.
  8. %
  9. %   This function has four optional arguments:
  10. %     DIM (default has both subsystems of equal dimension)
  11. %     STR (default 2)
  12. %     VERBOSE (default 1)
  13. %     TOL (default eps^(3/8))
  14. %
  15. %   SEP = IsSeparable(X,DIM,STR,VERBOSE,TOL) is as above, where DIM is
  16. %   a 1-by-2 vector containing the dimensions of the subsystems on which X
  17. %   acts.
  18. %
  19. %   STR is an integer that determines how hard the script should work to
  20. %   determine separability before giving up (STR = -1 means that the script
  21. %   won't stop working until it finds an answer). Other valid values are 0,
  22. %   1, 2, 3, ... In practice, if STR >= 4 then most computers will run out
  23. %   of memory and/or the sun will explode before computation completes.
  24. %
  25. %   VERBOSE is a flag (either 0 or 1) that indicates that the script will
  26. %   or will not display a line of text describing how it proved that X is
  27. %   or is not separable.
  28. %
  29. %   TOL is the numerical tolerance used throughout the script.
  30. %
  31. %   URL: http://www.qetlab.com/IsSeparable
  32.  
  33. %   requires: ApplyMap.m, CVX (http://cvxr.com/cvx/), FilterNormalForm.m,
  34. %             iden.m, InSeparableBall.m, IsPPT.m, IsPSD.m, jacobi_poly.m,
  35. %             kpNorm.m, MaxEntangled.m, OperatorSchmidtDecomposition.m,
  36. %             OperatorSchmidtRank.m, opt_args.m, opt_disp.m, PartialMap.m,
  37. %             PartialTrace.m, PartialTranspose.m, PermutationOperator.m,
  38. %             PermuteSystems.m, Realignment.m, SchmidtDecomposition.m,
  39. %             SchmidtRank.m, sporth.m, Swap.m, SwapOperator.m,
  40. %             SymmetricExtension.m, SymmetricInnerExtension.m,
  41. %             SymmetricProjection.m, TraceNorm.m
  42. %             
  43. %   author: Nathaniel Johnston (nathaniel@njohnston.ca)
  44. %   package: QETLAB
  45. %   last updated: November 14, 2014
  46.  
  47. function sep = IsSeparable(X,varargin)
  48.  
  49.     X = full(X);
  50.     if(~IsPSD(X))
  51.         error('IsSeparable:NotPSD','X is not positive semidefinite, so the idea of it being separable does not make sense.');
  52.     end
  53.  
  54.     lX = length(X);
  55.     rX = rank(X);
  56.     X = X/trace(X);
  57.     sep = -1;
  58.  
  59.     % set optional argument defaults: dim=sqrt(length(X)), str=2, verbose=1, tol=eps^(1/4)
  60.     [dim,str,verbose,tol] = opt_args({ round(sqrt(lX)), 2, 1, eps^(3/8) },varargin{:});
  61.     if(str == -1)
  62.         str = 1/eps; % keep going forever!
  63.     end
  64.  
  65.     % allow the user to enter a single number for dim
  66.     if(length(dim) == 1)
  67.         dim = [dim,lX/dim];
  68.         if abs(dim(2) - round(dim(2))) >= 2*lX*eps
  69.             error('IsSeparable:InvalidDim','If DIM is a scalar, it must evenly divide length(X); please provide the DIM array containing the dimensions of the subsystems.');
  70.         end
  71.         dim(2) = round(dim(2));
  72.     end
  73.     nD = min(dim);
  74.     xD = max(dim);
  75.     pD = prod(dim);
  76.  
  77.     if(nD == 1)
  78.         sep = 1;
  79.         opt_disp(['Every positive semidefinite matrix is separable when one of the local dimensions is 1.\n'],verbose);
  80.         return        
  81.     end
  82.  
  83.     XA = PartialTrace(X,2,dim);
  84.     XB = PartialTrace(X,1,dim);
  85.  
  86.     % pre-load various references
  87.     refs = {'A. Peres. Separability criterion for density matrices. Phys. Rev. Lett., 77:1413-1415, 1996.', ... % refs[1]
  88.             'M. Horodecki, P. Horodecki, and R. Horodecki. Separability of mixed states: Necessary and sufficient conditions. Phys. Lett. A, 223:1-8, 1996.', ...
  89.             'P. Horodecki, M. Lewenstein, G. Vidal, and I. Cirac. Operational criterion and constructive checks for the separability of low-rank density matrices. Phys. Rev. A, 62:032310, 2000.', ...
  90.             'K. Chen and L.-A. Wu. A matrix realignment method for recognizing entanglement. Quantum Inf. Comput., 3:193-202, 2003.', ...
  91.             'F. Verstraete, J. Dehaene, and B. De Moor. Normal forms and entanglement measures for multipartite quantum states. Phys. Rev. A, 68:012103, 2003.', ...
  92.             'K.-C. Ha and S.-H. Kye. Entanglement witnesses arising from exposed positive linear maps. Open Systems & Information Dynamics, 18:323-337, 2011.', ... % refs[6]
  93.             'O. Gittsovich, O. Guehne, P. Hyllus, and J. Eisert. Unifying several separability conditions using the covariance matrix criterion. Phys. Rev. A, 78:052319, 2008.', ...
  94.             'L. Gurvits and H. Barnum. Largest separable balls around the maximally mixed bipartite quantum state. Phys. Rev. A, 66:062311, 2002.', ...
  95.             'H.-P. Breuer. Optimal entanglement criterion for mixed quantum states. Phys. Rev. Lett., 97:080501, 2006.', ...
  96.             'W. Hall. Constructions of indecomposable positive maps based on a new criterion for indecomposability. E-print: arXiv:quant-ph/0607035, 2006.', ...
  97.             'A. C. Doherty, P. A. Parrilo, and F. M. Spedalieri. A complete family of separability criteria. Phys. Rev. A, 69:022308, 2004.', ... % refs[11]
  98.             'M. Navascues, M. Owari, and M. B. Plenio. A complete criterion for separability detection. Phys. Rev. Lett., 103:160404, 2009.', ...
  99.             'N. Johnston. Separability from spectrum for qubit-qudit states. Phys. Rev. A, 88:062330, 2013.', ...
  100.             'C.-J. Zhang, Y.-S. Zhang, S. Zhang, and G.-C. Guo. Entanglement detection beyond the cross-norm or realignment criterion. Phys. Rev. A, 77:060301(R), 2008.', ...
  101.             'R. Hildebrand. Semidefinite descriptions of low-dimensional separable matrix cones. Linear Algebra Appl., 429:901-932, 2008,', ...
  102.             'R. Hildebrand. Comparison of the PPT cone and the separable cone for 2-by-n systems. http://www-ljk.imag.fr/membres/Roland.Hildebrand/coreMPseminar2005_slides.pdf', ... % refs[16]
  103.             'D. Cariello. Separability for weak irreducible matrices. E-print: arXiv:1311.7275 [quant-ph], 2013.', ...
  104.             'L. Chen and D. Z. Djokovic. Separability problem for multipartite states of rank at most four. J. Phys. A: Math. Theor., 46:275304, 2013.', ...
  105.             'G. Vidal and R. Tarrach. Robustness of entanglement. Phys. Rev. A, 59:141-155, 1999.'};
  106.  
  107.     % Start with the really easy separability checks (we always do these,
  108.     % regardless of str).
  109.  
  110.     % PPT is used in so many tests, just always compute it
  111.     ppt = IsPPT(X,2,dim,tol);
  112.  
  113.     % Check the PPT criterion
  114.     if(~ppt)
  115.         sep = 0;
  116.         opt_disp(['Determined to be entangled via the PPT criterion. Reference:\n',refs{1},'\n'],verbose);
  117.         return
  118.  
  119.     % Sometimes the PPT criterion is also sufficient for separability.
  120.     elseif(pD <= 6 || min(dim) <= 1)
  121.         sep = 1;
  122.         opt_disp(['Determined to be separable via sufficiency of the PPT criterion in small dimensions. Reference:\n',refs{2},'\n'],verbose);
  123.         return
  124.  
  125.     % Be careful with the next test! Checking that rX <= max(dim) is
  126.     % *incorrect* here, and is not what is meant by refs{3}.
  127.     elseif(rX <= 3 || rX <= rank(XB) || rX <= rank(XA))
  128.         sep = 1;
  129.         opt_disp(['Determined to be separable via sufficiency of the PPT criterion for low-rank operators. Reference:\n',refs{3},'\n'],verbose);
  130.         return
  131.     end
  132.  
  133.     % Realignment (aka computable cross norm) criterion. Only do it if
  134.     % VERBOSE = 1, since we are about to do another test that is strictly
  135.     % stronger.
  136.     if(verbose && TraceNorm(Realignment(X,dim)) > 1)
  137.         sep = 0;
  138.         opt_disp(['Determined to be entangled via the realignment criterion. Reference:\n',refs{4},'\n'],verbose);
  139.         return    
  140.     end
  141.  
  142.     % Another test that is strictly stronger than the realignment criterion
  143.     if(TraceNorm(Realignment(X - kron(XA,XB),dim)) > sqrt((1-trace(XA^2))*(1-trace(XB^2))))
  144.         sep = 0;
  145.         opt_disp(['Determined to be entangled by using Theorem 1 of reference:\n',refs{14},'\n'],verbose);
  146.         return    
  147.     end
  148.  
  149.     lam = sort(real(eig(X)),'descend'); % eigenvalues of X
  150.  
  151.     % There are some separability tests that work specifically in the
  152.     % qubit-qudit (i.e., 2 \otimes n) case. Do these now.
  153.     if(nD == 2)
  154.         % check if X is separable from spectrum
  155.         if((lam(1) - lam(2*xD-1))^2 <= 4*lam(2*xD-2)*lam(2*xD) + tol^2)
  156.             sep = 1;
  157.             opt_disp(['Determined to be separable by inspecting its eigenvalues. Reference:\n',refs{13},'\n'],verbose);
  158.             return
  159.         end
  160.  
  161.         % For the rest of the block matrix tests, we need the 2-dimensional
  162.         % subsystem to be the *first* subsystem, so swap accordingly.
  163.         if(dim(1) > 2)
  164.             Xt = Swap(X,[1,2],dim);
  165.         else
  166.             Xt = X;
  167.         end
  168.  
  169.         % Check if Lemma 1 of refs{13} applies to X. Also check the
  170.         % Hildebrand 2xn results.
  171.         A = Xt(1:xD,1:xD);
  172.         B = Xt(1:xD,xD+1:2*xD);
  173.         C = Xt(xD+1:2*xD,xD+1:2*xD);
  174.  
  175.         % This test is weaker than the next one, so only do it if VERBOSE =
  176.         % 1.
  177.         if(verbose && norm(B-B') < tol^2)
  178.             sep = 1;
  179.             opt_disp(['Determined to be separable by being a block Hankel matrix:\n',refs{15},'\n'],verbose);
  180.             return
  181.         end
  182.  
  183.         if(rank(B-B') <= 1 && ppt)
  184.             sep = 1;
  185.             opt_disp(['Determined to be separable by being a perturbed block Hankel matrix:\n',refs{16},'\n'],verbose);
  186.             return
  187.         end
  188.  
  189.         X_2n_ppt_check = [(5/6)*A-C/6, B; B', (5/6)*C-A/6];
  190.         if(IsPSD(X_2n_ppt_check) && IsPPT(X_2n_ppt_check,2,[2,xD]))
  191.             sep = 1;
  192.             opt_disp(['Determined to be separable via the homothetic images approach of:\n',refs{16},'\n'],verbose);
  193.             return
  194.         end
  195.  
  196.         if(norm(B)^2 <= min(real(eig(A)))*min(real(eig(C))) + tol^2)
  197.             sep = 1;
  198.             opt_disp(['Determined to be separable by using Lemma 1 of reference:\n',refs{13},'\n'],verbose);
  199.             return
  200.         end
  201.     end
  202.  
  203.     % There are conditions that are both necessary and sufficient when both
  204.     % dimensions are 3 and the rank is 4
  205.     if(rX == 4 && nD == 3 && xD == 3)
  206.         % this method computes more determinants than are actually
  207.         % necessary, but the speed loss isn't too great
  208.         ranX = orth(X);
  209.         for j = 6:-1:1
  210.             for k = 7:-1:j+1
  211.                 for l = 8:-1:k+1
  212.                     for m = 9:-1:l+1
  213.                         p(j,k,l,m) = det(ranX([j,k,l,m],:));
  214.                     end
  215.                 end
  216.             end
  217.         end
  218.  
  219.         F = det([p(1,2,4,5), p(1,3,4,6), p(2,3,5,6), p(1,2,4,6)+p(1,3,4,5), p(1,2,5,6)+p(2,3,4,5), p(1,3,5,6)+p(2,3,4,6);
  220.                  p(1,2,7,8), p(1,3,7,9), p(2,3,8,9), p(1,2,7,9)+p(1,3,7,8), p(1,2,8,9)+p(2,3,7,8), p(1,3,8,9)+p(2,3,7,9);
  221.                  p(4,5,7,8), p(4,6,7,9), p(5,6,8,9), p(4,5,7,9)+p(4,6,7,8), p(4,5,8,9)+p(5,6,7,8), p(4,6,8,9)+p(5,6,7,9);
  222.                  p(1,2,4,8)-p(1,2,5,7), p(1,3,4,9)-p(1,3,6,7), p(2,3,5,9)-p(2,3,6,8), p(1,2,4,9)-p(1,2,6,7)+p(1,3,4,8)-p(1,3,5,7), p(1,2,5,9)-p(1,2,6,8)+p(2,3,4,8)-p(2,3,5,7), p(1,3,5,9)-p(1,3,6,8)+p(2,3,4,9)-p(2,3,6,7);
  223.                  p(1,4,5,8)-p(2,4,5,7), p(1,4,6,9)-p(3,4,6,7), p(2,5,6,9)-p(3,5,6,8), p(1,4,5,9)-p(2,4,6,7)+p(1,4,6,8)-p(3,4,5,7), p(1,5,6,8)-p(2,5,6,7)+p(2,4,5,9)-p(3,4,5,8), p(1,5,6,9)-p(3,4,6,8)+p(2,4,6,9)-p(3,5,6,7);
  224.                  p(1,5,7,8)-p(2,4,7,8), p(1,6,7,9)-p(3,4,7,9), p(2,6,8,9)-p(3,5,8,9), p(1,5,7,9)-p(2,4,7,9)+p(1,6,7,8)-p(3,4,7,8), p(1,5,8,9)-p(2,4,8,9)+p(2,6,7,8)-p(3,5,7,8), p(1,6,8,9)-p(3,4,8,9)+p(2,6,7,9)-p(3,5,7,9)]);
  225.  
  226.         sep = (abs(F) < max(tol^2,eps^(3/4))); % X is separable iff F is zero (or suffiently close to it, for numerical reasons)
  227.  
  228.         if(sep == 1)
  229.             opt_disp(['Determined to be separable by checking the Chow form. Reference:\n',refs{18},'\n'],verbose);
  230.         else
  231.             opt_disp(['Determined to be entangled by checking the Chow form. Reference:\n',refs{18},'\n'],verbose);
  232.         end
  233.         return
  234.     end
  235.  
  236.     % Check proximity of X with the maximally mixed state.
  237.     if(InSeparableBall(X))
  238.         sep = 1;
  239.         opt_disp(['Determined to be separable by closeness to the maximally mixed state. Reference:\n',refs{8},'\n'],verbose);
  240.         return    
  241.     end
  242.  
  243.     % Check if X is a rank-1 perturbation of the identity, which is
  244.     % necessarily separable if it's PPT, which we have already checked.
  245.     if(lam(2) - lam(pD) < tol^2)
  246.         sep = 1;
  247.         opt_disp(['Determined to be separable by being a small rank-1 perturbation of the maximally-mixed state. Reference:\n',refs{19},'\n'],verbose);
  248.         return   
  249.     end
  250.  
  251.     % Check tensor rank equalling 2
  252.     if(OperatorSchmidtRank(X,dim) <= 2)
  253.         sep = 1;
  254.         opt_disp(['Determined to be separable by having operator Schmidt rank at most 2. Reference:\n',refs{17},'\n'],verbose);
  255.         return
  256.     end
  257.  
  258.     % There is a family of known optimal positive maps in the qutrit-qutrit
  259.     % case. Check for entanglement using these.
  260.     if(dim(1) == 3 && dim(2) == 3)
  261.         phi = MaxEntangled(3,0,0);
  262.         for t=0:0.1:0.9
  263.             for j=1:2
  264.                 if(t>0)
  265.                     t = 1/t; % this is a weird way of using both t and 1/t as indices for the maps Phi we generate
  266.                 elseif(j>1)
  267.                     break;
  268.                 end
  269.                 a = (1-t)^2/(1-t+t^2);
  270.                 b = t^2/(1-t+t^2);
  271.                 c = 1/(1-t+t^2);
  272.                 Phi = diag([a+1,c,b,b,a+1,c,c,b,a+1]) - phi*phi';
  273.  
  274.                 if(~IsPSD(PartialMap(X,Phi,2,dim)))
  275.                     sep = 0;
  276.                     opt_disp(['Determined to be entangled via the positive map Phi[a,b,c] with a=',num2str(a),', b=',num2str(b),', c=',num2str(c),'. Reference:\n',refs{6},'\n'],verbose);
  277.                     return
  278.                 end
  279.             end
  280.         end
  281.     end
  282.  
  283.     % Use the Breuer-Hall positive maps (in even dimensions only) based on
  284.     % antisymmetric unitary matrices.
  285.     for p = 1:2
  286.         if(mod(dim(p),2) == 0)
  287.             phi = MaxEntangled(dim(p),0,0);
  288.             U = kron(eye(dim(p)),fliplr(diag([ones(dim(p)/2,1);-ones(dim(p)/2,1)])));
  289.             Phi = diag(ones(dim(p)^2,1)) - phi*phi' - U*SwapOperator(dim(p))*U';
  290.  
  291.             if(~IsPSD(PartialMap(X,Phi,p,dim)))
  292.                 sep = 0;
  293.                 opt_disp(['Determined to be entangled via the Breuer-Hall positive maps based on antisymmetric unitary matrices. References:\n',refs{9},'\n',refs{10},'\n'],verbose);
  294.                 return
  295.             end
  296.         end
  297.     end
  298.  
  299.     % The next tests for entanglement and separability are slightly more
  300.     % time-intensive, so we only do them if str >= 1.
  301.     if(str >= 1)
  302.         % Use the filter covariance matrix criterion (Filter CMC) for entanglement.
  303.         % This strengthens the realignment criterion considerably, but we do the
  304.         % realignment first anyway, since it is quicker and simpler (and thus
  305.         % provides a more useful reference).
  306.         try
  307.             xi = FilterNormalForm(X,dim);
  308.  
  309.             if(sum(xi) > min( sqrt(nD*xD*(nD-1)*(xD-1)), nD*xD*(1 - 1/nD + (nD^2 - 1)/xD + min(0,(xD^2-nD^2)/xD - (xD-1)))/2 ) + 10*xD*eps)
  310.                 sep = 0;
  311.                 opt_disp(['Determined to be entangled via the Filter Covariance Matrix Criterion. Reference:\n',refs{7},'\n'],verbose);
  312.                 return    
  313.             end
  314.  
  315.         % FilterNormalForm often is slow and results in an error when X is
  316.         % low-rank, but the partial transpose criterion already took care of those
  317.         % situations anyway. Nonetheless, we err on the side of caution and ignore
  318.         % the cases when the state doesn't have a filter normal form.
  319.         catch err
  320.             if(~strcmpi(err.identifier,'FilterNormalForm:NoFNF'))
  321.                 rethrow(err);
  322.             end
  323.         end
  324.  
  325.         % Try the Guhne method of enhancing the separability within X.
  326.         Xsep = X;
  327.         XSeig = sort(real(eig(Xsep)));
  328.         e_norm = norm(XSeig);
  329.         it_ct = 0;
  330.         while it_ct < 10
  331.             lb = -1;
  332.             % Iteratively find a product state whose maximal overlap with Xsep
  333.             % is close to optimal. This has a non-zero chance of failing, so we
  334.             % make sure that the obtained overlap isn't terrible before moving
  335.             % on.
  336.             while lb < XSeig(sum(dim)-1) - pD*eps;
  337.                 [lb,v] = sk_iterate(Xsep,1,dim);
  338.             end
  339.  
  340.             % compute the minimum amount that we can subtract off v*v' and
  341.             % still be positive semidefinite
  342.             gen_eig = eig(Xsep,v*v');
  343.             gen_eig = min(real(gen_eig(~isinf(gen_eig))));
  344.  
  345.             % actually subtract half that amount
  346.             Xsep2 = Xsep - gen_eig*(v*v')/2;
  347.             Xsep2 = Xsep2/trace(Xsep2);
  348.  
  349.             XSeig2 = sort(real(eig(Xsep)));
  350.             e_norm2 = norm(XSeig);
  351.             if(e_norm - e_norm2 < 10^(-4)) % don't use abs() -- we want e_norm2 to be smaller, not larger
  352.                 it_ct = it_ct + 1;
  353.                 if(e_norm - e_norm2 > -10^(-4)) % if e_norm2 isn't larger, then do the update anyway
  354.                     XSeig = XSeig2;
  355.                     Xsep = Xsep2;
  356.                     e_norm = e_norm2;
  357.                 end
  358.             else
  359.                 it_ct = 0;
  360.                 XSeig = XSeig2;
  361.                 Xsep = Xsep2;
  362.                 e_norm = e_norm2;
  363.  
  364.                 if(InSeparableBall(Xsep))
  365.                     sep = 1;
  366.                     opt_disp('Determined to be separable by subtracting product states until the operator was close to the maximally-mixed state.',verbose);
  367.                     return
  368.                 end
  369.  
  370.             end
  371.         end
  372.     end
  373.  
  374.     % If str is high enough, also try the stronger separability tests based on
  375.     % symmetric extensions. These are much more time-intensive.
  376.     for k = 2:str
  377.         if(~SymmetricExtension(X,k,dim,1,1))
  378.             sep = 0;
  379.             opt_disp(['Determined to be entangled by not having a ',num2str(k),'-copy PPT symmetric extension. Reference:\n',refs{11},'\n'],verbose);
  380.             return        
  381.         end
  382.  
  383.         if(SymmetricInnerExtension(Xsep,k,dim,1)) % Use Xsep instead of X, since Xsep is in some sense "more" separable, so it is easier to detect
  384.             sep = 1;
  385.             opt_disp(['Determined to be separable via the semidefinite program based on ',num2str(k),'-copy PPT symmetric extensions from reference:\n',refs{12},'\n'],verbose);
  386.             return
  387.         end
  388.     end
  389. end

References

  1. O. Gittsovich, O. Gühne, P. Hyllus, and J. Eisert. Unifying several separability conditions using the covariance matrix criterion. Phys. Rev. A, 78:052319, 2008. E-print: arXiv:0803.0757 [quant-ph]
  2. M. Horodecki and P. Horodecki. Reduction criterion of separability and limits for a class of distillation protocols. Phys. Rev. A, 59:4206–4216, 1999. E-print: arXiv:quant-ph/9708015
  3. R. F. Werner. Quantum states with Einstein-Podolsky-Rosen correlations admitting a hidden-variable model. Phys. Rev. A, 40(8):4277–4281.