# IsBlockPositive

(diff) ← Older revision | Latest revision (diff) | Newer revision → (diff)
 Other toolboxes required IsBlockPositive Determines whether or not an operator is block positive CVX IsPPTIsPSDIsSeparable Entanglement and separability

IsBlockPositive is a function that determines whether or not a bipartite operator is block positive (roughly speaking, it determines whether the operator is an entanglement witness or not). A value of 0 indicates that the operator is not block positive (and hence is not an entanglement witness), a value of 1 indicates that the operator is block positive (and hence is either an entanglement witness or is positive semidefinite), and a value of -1 indicates that the script was unable to determine whether or not the operator is block positive.

## Syntax

• IBP = IsBlockPositive(X)
• IBP = IsBlockPositive(X,K)
• IBP = IsBlockPositive(X,K,DIM)
• IBP = IsBlockPositive(X,K,DIM,STR)
• IBP = IsBlockPositive(X,K,DIM,STR,TOL)
• [IBP,WIT] = IsBlockPositive(X,K,DIM,STR,TOL)

## Argument descriptions

• X: A bipartite Hermitian operator.
• K (optional, default 1): A positive integer indicating that the script should determine whether or not X is K-block positive (i.e., whether or not it remains nonnegative under left and right multiplication by vectors with Schmidt rank <= K).
• 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 block positivity 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.
• TOL (optional, default eps^(3/8)): The numerical tolerance used throughout the script.

## Examples

### The swap operator is block positive

The swap operator always block positive, since it is the Choi matrix of the transpose map (and it is easy to see that the transpose map is positive). We can verify this in small dimensions as follows:

>> for d = 2:10
res(d) = IsBlockPositive(SwapOperator(d));
end
>> res(2:10)

ans =

1     1     1     1     1     1     1     1     1

However, the swap operator is not 2-block positive, which we can verify as follows:

>> IsBlockPositive(SwapOperator(5),2)

ans =

0

### The Choi map is positive

The Choi map is a well-known example of a positive linear map that acts on 3-by-3 matrices. The following code verifies that this map is indeed positive, but not 2-positive:

>> IsBlockPositive(ChoiMap())

ans =

1

>> IsBlockPositive(ChoiMap(),2)

ans =

0

### A complicated operator on 4 local dimensions

A somewhat messy positive linear map was introduced in Lemma 3 of . The following code verifies that this map is indeed positive in the $t = 1$ case.

>> b0 = Bell(0);
>> b3 = Bell(3);
>> v0 = kron(b0,b0);
>> Y = (kron(eye(4),b0*b0')/2 + kron(b3*b3',PartialTranspose(b3*b3',1)))/3 - v0*v0'/4;
>> Y = Swap(Y,[2,3],[2,2,2,2]); % this is the Choi matrix of the linear map given in Lemma 3
>> IsBlockPositive(Y)

ans =

1

## Source code

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

1. %%  ISBLOCKPOSITIVE    Determines whether or not an operator is block positive
2. %   This function has one required input argument:
3. %     X: a square matrix
4. %
5. %   IBP = IsBlockPositive(X) is either -1, 0, or 1. A value of 1 indicates
6. %   that X is block positive, a value of 0 indicates that X is not block
7. %   positive, and a value of -1 indicates that the block positivity of X
8. %   could not be determined.
9. %
10. %   This function has four optional input arguments:
11. %     K (default 1)
12. %     DIM (default has both subsystems of equal dimension)
13. %     STR (default 2)
14. %     TOL (default eps^(3/8))
15. %
16. %   [IBP,WIT] = IsBlockPositive(X,K,DIM,STR,TOL) is as above, where DIM is
17. %   a 1-by-2 vector containing the dimensions of the subsystems on which X
18. %   acts.
19. %
20. %   K is the "level" of positivity -- the script checks whether or not X is
21. %   K-block positive. That is, it checks whether or not X has positive
22. %   expectation with all states of Schmidt rank <= K.
23. %
24. %   STR is an integer that determines how hard the script should work to
25. %   determine block positivity before giving up (STR = -1 means that the
26. %   script won't stop working until it finds an answer). Other valid values
27. %   are 0, 1, 2, 3, ... In practice, if STR >= 4 then most computers will
28. %   run out of memory and/or the sun will explode before computation
29. %   completes.
30. %
31. %   TOL is the numerical tolerance used throughout the script.
32. %
33. %   WIT is a witness that verifies that X is (or is not) block positive.
34. %
35. %   URL: http://www.qetlab.com/IsBlockPositive
36. 
37. %   requires: cvx (http://cvxr.com/cvx/), iden.m, IsPSD.m, kpNorm.m,
38. %             MaxEntangled.m, normalize_cols.m, opt_args.m, PartialMap.m,
39. %             PartialTrace.m, PartialTranspose.m, PermuteSystems.m,
40. %             Realignment.m, SchmidtDecomposition.m, SchmidtRank.m,
41. %             sk_iterate.m, SkOperatorNorm.m, SkVectorNorm.m, sporth.m,
42. %             Swap.m
43. %
44. %   author: Nathaniel Johnston (nathaniel@njohnston.ca)
45. %   package: QETLAB
46. %   last updated: September 23, 2014
47. 
48. function [ibp,wit] = IsBlockPositive(X,varargin)
49. 
50.     wit = 0;
51.     X = full(X);
52.     dX = size(X);
53. 
54.     % set optional argument defaults: k=1, dim=sqrt(length(X)), str=2,
55.     % tol=eps^(3/8)
56.     [k,dim,str,tol] = opt_args({ 1, round(sqrt(dX(1))), 2, eps^(3/8) },varargin{:});
57.     if(str == -1)
58.         str = 1/eps; % keep going forever!
59.     end
60. 
61.     % Make sure that X is Hermitian.
62.     if(dX(1) ~= dX(2) || (max(max(abs(X-X'))) > tol))
63.         ibp = 0;
64.         return
65.     end
66. 
67.     % allow the user to enter a single number for dim
68.     if(length(dim) == 1)
69.         dim = [dim,dX(1)/dim];
70.         if abs(dim(2) - round(dim(2))) >= 2*dX(1)*eps
71.             error('IsBlockPositive:InvalidDim','If DIM is a scalar, it must evenly divide length(X); please provide the DIM array containing the dimensions of the subsystems.');
72.         end
73.         dim(2) = round(dim(2));
74.     end
75. 
76.     % When a local dimension is small, block positivity is trivial.
77.     if(min(dim) <= k)
78.         ibp = IsPSD(X);
79.         wit = zeros(dX);
80.         return
81.     end
82. 
83.     op_norm = norm(X);
84.     Y = op_norm*speye(dX(1)) - X; % we compute the S(k)-norm of this operator
85. 
86.     [lb,lwit,ub,uwit] = SkOperatorNorm(Y,k,dim,str,op_norm,tol); % compute the norm
87. 
88.     if(ub <= op_norm*(1 + tol)) % block positive
89.         ibp = 1;
90.         wit = uwit;
91.     elseif(lb >= op_norm*(1 - tol)) % not block positive
92.         ibp = 0;
93.         wit = lwit;
94.     else % not sure :(
95.         ibp = -1;
96.     end
97. end