# Distinguishability

 Other toolboxes required Distinguishability Computes the maximum probability of distinguishing quantum states CVX ChannelDistinguishabilityLocalDistinguishability Distinguishing objects

Distinguishability is a function that computes the maximum probability of distinguishing two or more quantum states. That is, this function computes the maximum probability of winning the following game: You are given a complete description of a set of $k$ quantum states $\rho_1, \ldots, \rho_k$, and then are given one of those $k$ states, and asked to determine (via quantum measurement) which state was given to you.

## Syntax

• DIST = Distinguishability(X)
• DIST = Distinguishability(X,P)
• [DIST,MEAS] = Distinguishability(X,P)

## Argument descriptions

### Input arguments

• X: The quantum states to be distinguished. X can either be a cell containing 2 or more density matrices, or X can be a matrix whose columns are pure vector states.
• P (optional, default [1/k, 1/k, ..., 1/k], where k is the number of quantum states): A vector whose j-th entry is the probability that the state $\rho_j$ is given to you in the game described above. All entries must be non-negative, and the entries of this vector must sum to 1.

### Output arguments

• DIST: The maximum probability of distinguishing the states specified by X.
• MEAS (optional): A cell containing optimal measurement operators that distinguish the states specified by X with probability DIST.

## Examples

### Orthogonal states can be perfectly distinguished

Any number of quantum states can be perfectly distinguished (i.e., distinguished with probability 1) if they are mutually orthogonal. The following code generates a random $6\times 6$ unitary matrix (i.e., a matrix with orthogonal pure states as columns) and verifies that those pure states are perfectly distinguishable:

>> Distinguishability(RandomUnitary(6))

ans =

1

### Two states

The maximum probability of distinguishing two quantum states $\rho$ and $\sigma$ is exactly $\frac{1}{2} + \frac{1}{4}\|\rho - \sigma\|_1$[1], where $\|\cdot\|_1$ is the trace norm. We can verify this in a special case as follows:

>> rho = RandomDensityMatrix(4);
>> sigma = RandomDensityMatrix(4);
>> Distinguishability({rho, sigma})

ans =

0.7762

>> 1/2 + TraceNorm(rho - sigma)/4

ans =

0.7762

### Three or more states

We can also compute the maximum probability of distinguishing three or more states, but no simple formula is known in this case.

>> for j = 1:6
rho{j} = RandomDensityMatrix(4);
end
>> Distinguishability(rho)

ans =

0.4156

## Source code

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

1. %%  DISTINGUISHABILITY    Computes the maximum probability of distinguishing quantum states
2. %   This function has one required input argument:
3. %     X: a cell containing density matrices or a single matrix containing
4. %        pure states as its column vectors
5. %
6. %   DIST = Distinguishability(X) is the maximum probability of
7. %   distinguishing the quantum states specified by X. X can either be a
8. %   cell containing 2 or more density matrices to be distinguished, or X
9. %   can be a matrix whose columns are pure vector states to be
10. %   distinguished.
11. %
12. %   This function has one optional input argument:
13. %     P (default [1/k, ..., 1/k], where k is the number of states)
14. %
15. %   [DIST,MEAS] = Distinguishability(X,P) is the maximum probability of
16. %   distinguishing the quantum states specified by X, where the vector P
17. %   contains the probability that each state is chosen (by default, the
18. %   states are chosen uniformly at random). MEAS is a cell containing the
19. %   optimal measurement operators.
20. %
21. %   URL: http://www.qetlab.com/Distinguishability
22. 
23. %   requires: cvx (http://cvxr.com/cvx/), kpNorm.m, opt_args.m, TraceNorm.m
24. %
25. %   author: Nathaniel Johnston (nathaniel@njohnston.ca)
26. %   package: QETLAB
27. %   last updated: October 6, 2014
28. 
29. function [dist,meas] = Distinguishability(X,varargin)
30. 
31. if(iscell(X))
32.     num_ops = length(X);
33.     dim = length(X{1});
34. else
35.     [dim,num_ops] = size(X);
36. end
37. 
38. % set optional argument defaults: p = [1/k,1/k,...,1/k], where k = # of states
39. [p] = opt_args({ ones(1,num_ops)/num_ops },varargin{:});
40. 
41. if(abs(sum(p) - 1) > num_ops^2*eps || length(p) ~= num_ops)
42.     error('Distinguishability:InvalidP','The vector P must be a probability distribution of the same length as the number of states: its elements must be non-negative and they must sum to 1.');
43. end
44. 
45. if(num_ops == 1 || max(p) >= 1) % of course we can distinguish 1 object
46.     dist = 1;
47.     if(nargout > 1)
48.         meas = eye(dim); % optimal measurements is trivial in this case
49.     end
50.     return
51. end
52. 
53. % X is a cell of density matrices
54. if(iscell(X))
55.     for j = 1:num_ops % make sure that the density operators are scaled
56.         X{j} = X{j}/trace(X{j});
57.     end
58. 
59.     % There is a closed-form expression for the distinguishability of two
60.     % density matrices.
61.     if(num_ops == 2)
62.         dist = 1/2 + TraceNorm(p(1)*X{1} - p(2)*X{2})/2;
63.         if(nargout > 1) % construct optimal measurements, if requested
64.             [v,d] = eig(p(1)*X{1} - p(2)*X{2});
65.             d = diag(d);
66.             pind = find(d >= 0);
67.             meas{1} = v(:,pind)*v(:,pind)';
68.             meas{2} = eye(dim) - meas{1};
69.         end
70.         return
71.     end
72. 
73.     % Check to see if the states are mutually orthogonal, and return 1 if
74.     % they are.
75.     if(num_ops <= dim)
76.         mut_orth = 1;
77.         for j = 1:num_ops
78.             for k = j+1:num_ops
79.                 if(max(max(abs(X{j}*X{k}))) > eps*dim^2)
80.                     mut_orth = 0;
81.                     break;
82.                 end
83.             end
84.         end
85.         if(mut_orth == 1) % states are mutually orthogonal
86.             dist = 1;
87.             if(nargout > 1) % construct optimal measurements, if requested
88.                 meas_sum = zeros(dim);
89.                 for j = num_ops:-1:2 % pre-allocate for speed
90.                     oX = orth(X{j});
91.                     meas{j} = oX*oX';
92.                     meas_sum = meas_sum + meas{j};
93.                 end
94.                 meas{1} = eye(dim) - meas_sum;
95.             end
96.             return;
97.         end
98.     end
99. 
100. % X is a matrix whose columns are pure states
101. else
102.     X = normalize_cols(X); % make sure that the columns have unit length
103. 
104.     % There is a closed-form expression for the distinguishability of two
105.     % pure states.
106.     if(num_ops == 2)
107.         dist = 1/2 + sqrt(2*(p(1)^2 + p(2)^2) - 4*p(1)*p(2)*abs(X(:,1)'*X(:,2))^2)/2;
108.         if(nargout > 1) % construct optimal measurements, if requested
109.             [v,d] = eig(p(1)*X(:,1)*X(:,1)' - p(2)*X(:,2)*X(:,2)');
110.             d = diag(d);
111.             pind = find(d >= 0);
112.             meas{1} = v(:,pind)*v(:,pind)';
113.             meas{2} = eye(dim) - meas{1};
114.         end
115.         return
116.     end
117. 
118.     % Check to see if the states are mutually orthogonal, and return 1 if
119.     % they are.
120.     if(num_ops <= dim)
121.         X2 = X'*X;
122.         if(max(max(abs(X2 - diag(diag(X2))))) < eps*dim^2)
123.             dist = 1; % the states are orthogonal, so perfectly distinguishable
124.             if(nargout > 1) % construct optimal measurements, if requested
125.                 meas_sum = zeros(dim);
126.                 for j = num_ops:-1:2 % pre-allocate for speed
127.                     meas{j} = X(:,j)*X(:,j)';
128.                     meas_sum = meas_sum + meas{j};
129.                 end
130.                 meas{1} = eye(dim) - meas_sum;
131.             end
132.             return
133.         end
134.     end
135. 
136.     % Turn pure states into density operator form.
137.     Y = X; X = {};
138.     for j = num_ops:-1:1 % pre-allocate for speed
139.         X{j} = Y(:,j)*Y(:,j)';
140.     end
141. end
142. 
143. % For 3 or more density matrices, we have to use semidefinite programming.
144. cvx_begin sdp quiet
145.     cvx_precision default;
146.     variable P(dim,dim,num_ops) hermitian
147. 
148.     P_tr = 0;
149.     P_sum = zeros(dim);
150.     for j = 1:num_ops
151.         P_tr = P_tr + p(j)*trace(P(:,:,j)*X{j});
152.         P_sum = P_sum + P(:,:,j);
153.     end
154.     P_tr = P_tr + P_tr';
155. 
156.     maximize P_tr
157.     subject to
158.         P_sum == eye(dim);
159.         P >= 0;
160. cvx_end
161. dist = real(cvx_optval)/2;
162. 
163. % Also return the optimal measurements, if requested.
164. if(nargout > 1)
165.     meas = mat2cell(reshape(P,dim,dim*num_ops),dim,dim*ones(1,num_ops));
166. end

## References

1. John Watrous. Theory of Quantum Information lecture notes, Fall 2011.