SymmetricInnerExtension

From QETLAB
Jump to: navigation, search
SymmetricInnerExtension
Determines whether or not an operator has a symmetric inner extension

Other toolboxes required cvx
Related functions IsPPT
IsSeparable
SymmetricExtension
SymmetricProjection
Function category Entanglement and separability

SymmetricInnerExtension is a function that determines whether or not a given positive semidefinite operator is in the cone defined in reference [1] that approximates the set of separable operators from the inside, based on operators with symmetric extensions. Importantly, if this script returns a value of 1 then we know that the given operator is separable (if the script returns a value of 0, we know nothing about whether the operator is separable or entangled).

Syntax

  • EX = SymmetricInnerExtension(X)
  • EX = SymmetricInnerExtension(X,K)
  • EX = SymmetricInnerExtension(X,K,DIM)
  • EX = SymmetricInnerExtension(X,K,DIM,PPT)
  • EX = SymmetricInnerExtension(X,K,DIM,PPT,TOL)
  • [EX,WIT] = SymmetricInnerExtension(X,K,DIM,PPT,TOL)

Argument descriptions

Input arguments

  • X: A positive semidefinite operator.
  • K (optional, default 2): The number of copies of the second subsystem in the desired symmetric extension.
  • 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.
  • PPT (optional, default 0): A flag (either 1 or 0) that indicates whether or not the desired symmetric extension must have positive partial transpose.
  • TOL (optional, default eps^(1/4)): The numerical tolerance used throughout this script. It is recommended that this is left at the default value unless numerical problems arise and the script has difficulty determining whether or not X has a symmetric inner extension.

Output arguments

  • EX: A flag (either 1 or 0) indicating that X does or does not have a symmetric inner extension of the desired type. A value of 1 means that X is separable.
  • WIT (optional): A witness that verifies that the answer provided by EX is correct. If EX = 1 (i.e., X is in the "inner" symmetric extension cone) then WIT is a symmetric extension of the operator \(\sigma_{AB}\) from [1], and thus acts as a witness that verifies that EX = 1 is correct. If X is not in this cone (i.e., EX = 0) then WIT is an operator with trace(WIT*X) = -1 but trace(WIT*Y) >= 0 for all operators in the described cone. Note that WIT may not be an entanglement witness!

Examples

The following code takes an unnormalized quantum state and shows that it is separable, since it has a symmetric inner extension.

>> rho = [11     4    -1    -1    -1    -1    -1    -1    -1
           4    11    -1    -1    -1    -1    -1    -1    -1
          -1    -1    11    -1    -1     4    -1    -1    -1
          -1    -1    -1    11    -1    -1     4    -1    -1
          -1    -1    -1    -1    16    -1    -1    -1    -1
          -1    -1     4    -1    -1    11    -1    -1    -1
          -1    -1    -1     4    -1    -1    11    -1    -1
          -1    -1    -1    -1    -1    -1    -1    11     4
          -1    -1    -1    -1    -1    -1    -1     4    11];
>> SymmetricInnerExtension(rho,2,[3,3],1) % we set PPT=1, since this makes the separability test stronger
 
ans =
 
     1

If we subtract the identity from rho, then intuitively it becomes "less separable". We can verify this as follows, where we see that setting K = 2 is no longer sufficient to prove that the state is separable, and instead we must set K = 3. In general, the higher K is, the stronger the test for separability is, but the longer and more memory-intensive the computation is.

>> SymmetricInnerExtension(rho-eye(9),2,[3,3],1)
 
ans =
 
     0
 
>> SymmetricInnerExtension(rho-eye(9),3,[3,3],1)
 
ans =
 
     1

Notes

If your goal is to prove that an operator is separable, it is probably in your best interest to see the optional argument PPT to 1. It only defaults to 0 to be consistent with the SymmetricExtension function.

Source code

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

  1. %%  SYMMETRICINNEREXTENSION    Determines whether or not an operator has a symmetric inner extension
  2. %   This function has one required argument:
  3. %     X: a positive semidefinite matrix
  4. %
  5. %   EX = SymmetricInnerExtension(X) is either 1 or 0, indicating that X is
  6. %   or is not in the cone defined in reference [1] that approximates the
  7. %   set of separable operators from the inside, based on operators with
  8. %   2-copy positive partial transpose Bosonic symmetric extensions.
  9. %
  10. %   This function has four optional arguments:
  11. %     K (default 2)
  12. %     DIM (default has both subsystems of equal dimension)
  13. %     PPT (default 0)
  14. %     TOL (default eps^(1/4))
  15. %
  16. %   [EX,WIT] = SymmetricInnerExtension(X,K,DIM,PPT,TOL) determines whether
  17. %   or not X is in the cone defined in [1] that approximates the set of
  18. %   separable operators from the insides, based on operators with K-copy
  19. %   Bosonic symmetric extensions, which is positive partial transpose if
  20. %   PPT = 1 and does not have to be positive partial transpose if PPT = 0.
  21. %   If X is in this cone (i.e., EX = 1) then WIT is a symmetric extension
  22. %   of the operator sigma_{AB} from [1], and thus acts as a witness that
  23. %   verifies that EX = 1 is correct. If X is not in this cone (i.e., EX =
  24. %   0) then WIT is an operator with trace(WIT*X) = -1 but trace(WIT*Y) >= 0
  25. %   for all operators in the described cone. Note that WIT may not be an
  26. %   entanglement witness!
  27. %
  28. %   K is the desired number of copies of the second subsystem. DIM is a
  29. %   1-by-2 vector containing the dimensions of the subsystems on which X
  30. %   acts. PPT is a flag (either 0 or 1) indicating whether the desired
  31. %   symmetric extension must have positive partial transpose. TOL is the
  32. %   numerical tolerance used when performing the semidefinite program
  33. %   feasibility check.
  34. %
  35. %   URL: http://www.qetlab.com/SymmetricInnerExtension
  36. %
  37. %   References:
  38. %   [1] M. Navascues, M. Owari, and M. B. Plenio. Complete Criterion for
  39. %       Separability Detection. Physical Review Letters, 103:160404, 2009.
  40.  
  41. %   requires: cvx (http://cvxr.com/cvx/), IsPPT.m, IsPSD.m, jacobi_poly.m,
  42. %             opt_args.m, PartialTrace.m, PartialTranspose.m,
  43. %             PermutationOperator.m, PermuteSystems.m, sporth.m,
  44. %             SymmetricProjection.m
  45. %             
  46. %   author: Nathaniel Johnston (nathaniel@njohnston.ca)
  47. %   package: QETLAB
  48. %   last updated: December 16, 2014
  49.  
  50. function [ex,wit] = SymmetricInnerExtension(X,varargin)
  51.  
  52. lX = length(X);
  53.  
  54. % set optional argument defaults: k=2, dim=sqrt(length(X)), ppt=1, tol=eps^(1/4)
  55. [k,dim,ppt,tol] = opt_args({ 2, round(sqrt(lX)), 0, eps^(1/4) },varargin{:});
  56.  
  57. % allow the user to enter a single number for dim
  58. if(length(dim) == 1)
  59.     dim = [dim,lX/dim];
  60.     if abs(dim(2) - round(dim(2))) >= 2*lX*eps
  61.         error('SymmetricInnerExtension:InvalidDim','If DIM is a scalar, it must evenly divide length(X); please provide the DIM array containing the dimensions of the subsystems.');
  62.     end
  63.     dim(2) = round(dim(2));
  64. end
  65.  
  66. if(k >= 2)
  67.     if(ppt)
  68.         en = min(1-roots(jacobi_poly(dim(2)-2,mod(k,2),floor(k/2)+1)))*dim(2)/(2*(dim(2)-1));
  69.     end
  70.  
  71.     sdp_dim = [dim(1),dim(2)*ones(1,k)];
  72.     sdp_prod_dim = dim(1)*nchoosek(k+dim(2)-1, dim(2)-1);
  73.     P = kron(speye(dim(1)),SymmetricProjection(dim(2),k,1));
  74.  
  75.     cvx_begin sdp quiet
  76.         cvx_precision(tol);
  77.         variable rho(sdp_prod_dim,sdp_prod_dim) hermitian
  78.         if(nargout > 1) % don't want to compute the dual solution in general (especially not if this is called within CVX)
  79.             dual variable W
  80.         end
  81.         subject to
  82.             rho >= 0;
  83.             if(nargout > 1) % the code here gets a bit repetitive, but I'm not sure of a better way to do it
  84.                 if(ppt)
  85.                     W : (1-en)*PartialTrace(P*rho*P',3:k+1,sdp_dim) + en*kron(PartialTrace(P*rho*P',2:k+1,sdp_dim),eye(dim(2)))/dim(2) == X;
  86.                 else
  87.                     W : PartialTrace(P*rho*P',3:k+1,sdp_dim)*k/(k+dim(2)) + kron(PartialTrace(P*rho*P',2:k+1,sdp_dim),eye(dim(2)))/(dim(2)+k) == X;                
  88.                 end
  89.             else
  90.                 if(ppt)
  91.                     (1-en)*PartialTrace(P*rho*P',3:k+1,sdp_dim) + en*kron(PartialTrace(P*rho*P',2:k+1,sdp_dim),eye(dim(2)))/dim(2) == X;
  92.                 else
  93.                     PartialTrace(P*rho*P',3:k+1,sdp_dim)*k/(k+dim(2)) + kron(PartialTrace(P*rho*P',2:k+1,sdp_dim),eye(dim(2)))/(dim(2)+k) == X;                
  94.                 end
  95.             end
  96.             if(ppt)
  97.                 PartialTranspose(P*rho*P',1:ceil(k/2)+1,sdp_dim) >= 0;
  98.             end
  99.     cvx_end
  100.  
  101.     ex = 1-min(cvx_optval,1); % CVX-safe way to map (0,Inf) to (1,0)
  102.     if(~isa(X,'cvx')) % make the output prettier if it's not a CVX input
  103.         ex = round(ex);
  104.  
  105.         % Deal with error messages and witnesses.
  106.         if(nargout > 1 && strcmpi(cvx_status,'Solved'))
  107.             wit = P*rho*P';
  108.         elseif(strcmpi(cvx_status,'Inaccurate/Solved'))
  109.             wit = P*rho*P';
  110.             warning('SymmetricInnerExtension:NumericalProblems','Numerical problems encountered by CVX. Consider adjusting the tolerance level TOL and re-running the script.');
  111.         elseif(nargout > 1 && strcmpi(cvx_status,'Infeasible'))
  112.             wit = -W;
  113.         elseif(strcmpi(cvx_status,'Inaccurate/Infeasible'))
  114.             if(nargout > 1)
  115.                 wit = -W;
  116.             end
  117.             warning('SymmetricInnerExtension:NumericalProblems','Numerical problems encountered by CVX. Consider adjusting the tolerance level TOL and re-running the script.');
  118.         elseif(strcmpi(cvx_status,'Unbounded') || strcmpi(cvx_status,'Inaccurate/Unbounded') || strcmpi(cvx_status,'Failed'))
  119.             error('SymmetricInnerExtension:NumericalProblems',strcat('Numerical problems encountered (CVX status: ',cvx_status,'). Please try adjusting the tolerance level TOL.'));
  120.         end
  121.     end
  122. else
  123.     error('SymmetricInnerExtension:InvalidK','K must be an integer >= 2.');
  124. end

References

  1. 1.0 1.1 M. Navascués, M. Owari, and M. B. Plenio. Complete Criterion for Separability Detection. Physical Review Letters, 103:160404, 2009. E-print: arXiv:0906.2735 [quant-ph]