-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathIsAQAssemblage.m
164 lines (139 loc) · 6.28 KB
/
IsAQAssemblage.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
function [output,Fr,Fx,Fy,F,beta,Gamma] = IsAQAssemblage(Sr,Sx,Sy,S,ma,mb)
% given an assemblage in CG notation, this function determines whether it belongs to the
% Almost-Quantum outer approximation to the set of quantum assemblages.
% This function returns out=0 if the assemblage is outside and out=1 if inside.
%
% Steering scenario with two uncharacterised parties (Alice and Bob), and
% one characterised one (Charlie).
%
% mb = number of measurements for Alice
% mc = number of measurements for Bob
%
% this code assumes that all measurements are dichotomic
d=size(Sr,1);
dimG=(ma+1)*(mb+1)*d;
cvx_begin sdp quiet
variable Gamma(dimG,dimG) hermitian % moment matrix
variable pi_r(d,d) hermitian
variable pix(d,d,ma) hermitian
variable piy(d,d,mb) hermitian
variable pixy(d,d,ma*mb) hermitian semidefinite
dual variable Fr
dual variables Fx{ma}
dual variables Fy{mb}
dual variables F{ma*mb}
minimise trace(pi_r)
Gamma == hermitian_semidefinite(dimG);
for x = 1:ma
for y = 1:mb
pix(:,:,x)-pixy(:,:,x+ma*(y-1)) == hermitian_semidefinite(d);
piy(:,:,y)-pixy(:,:,x+ma*(y-1)) == hermitian_semidefinite(d);
pi_r-pix(:,:,x)-piy(:,:,y)+pixy(:,:,x+ma*(y-1)) == hermitian_semidefinite(d);
end
end
% First row
Fr : Gamma(1:d,1:d)==Sr+pi_r;
% Alice's marginal assemblage
for k=1:ma
Fx{k} : Gamma(1:d,k*d+1:(k+1)*d)==Sx(:,:,k)+pix(:,:,k);
end
% Bob's marginal assemblage
for k=1:mb
Fy{k} : Gamma(1:d,k*d*(ma+1)+1:k*d*(ma+1)+d)==Sy(:,:,k) + piy(:,:,k);
end
% the assemblage for all-but-one outcome
for k=1:ma
for j=1:mb
F{ma*(j-1)+k} : Gamma(1:d,j*d*(ma+1)+k*d+1:j*d*(ma+1)+k*d+d)...
== S(:,:,ma*(j-1)+k) + pixy(:,:,ma*(j-1)+k);
end
end
% Conditions on Gamma to be AQ
for y=0:ma
for z=0:mb
for yp=0:ma
for zp=0:mb
if y==0 && z==0 % reductions of type AB, where it can be A=1 and/or B=1
r1=1;
c1=zp*d*(ma+1)+yp*d+1;
r2=yp*d+1;
c2=zp*d*(ma+1)+1;
r3=yp*d+1;
c3=zp*d*(ma+1)+yp*d+1;
r4=zp*d*(ma+1)+1;
c4=zp*d*(ma+1)+yp*d+1;
r5=zp*d*(ma+1)+yp*d+1;
c5=zp*d*(ma+1)+1;
r6=zp*d*(ma+1)+yp*d+1;
c6=zp*d*(ma+1)+yp*d+1;
r7=zp*d*(ma+1)+1;
c7=yp*d+1;
r8=zp*d*(ma+1)+yp*d+1;
c8=1;
r9=zp*d*(ma+1)+yp*d+1;
c9=yp*d+1;
Gamma(r1:r1+d-1,c1:c1+d-1)==Gamma(r2:r2+d-1,c2:c2+d-1);
Gamma(r1:r1+d-1,c1:c1+d-1)==Gamma(r3:r3+d-1,c3:c3+d-1);
Gamma(r1:r1+d-1,c1:c1+d-1)==Gamma(r4:r4+d-1,c4:c4+d-1);
Gamma(r1:r1+d-1,c1:c1+d-1)==Gamma(r5:r5+d-1,c5:c5+d-1);
Gamma(r1:r1+d-1,c1:c1+d-1)==Gamma(r6:r6+d-1,c6:c6+d-1);
Gamma(r1:r1+d-1,c1:c1+d-1)==Gamma(r7:r7+d-1,c7:c7+d-1);
Gamma(r1:r1+d-1,c1:c1+d-1)==Gamma(r8:r8+d-1,c8:c8+d-1);
Gamma(r1:r1+d-1,c1:c1+d-1)==Gamma(r9:r9+d-1,c9:c9+d-1);
end
if y==0 && z~=0 && yp~=0 && zp~=0 % Reduction of type ABB', none = 1.
r1=z*d*(ma+1)+1;
c1=zp*d*(ma+1)+yp*d+1;
r2=z*d*(ma+1)+yp*d+1;
c2=zp*d*(ma+1)+1;
r3=z*d*(ma+1)+yp*d+1;
c3=zp*d*(ma+1)+yp*d+1;
Gamma(r1:r1+d-1,c1:c1+d-1)==Gamma(r2:r2+d-1,c2:c2+d-1);
Gamma(r1:r1+d-1,c1:c1+d-1)==Gamma(r3:r3+d-1,c3:c3+d-1);
end
if z~=0 && yp==0 && zp~=0 % (AB,A')=(A'B,A)', B can be =1.
if z~=zp
r1=z*d*(ma+1)+y*d+1;
c1=zp*d*(ma+1)+1;
r2=zp*d*(ma+1)+y*d+1;
c2=z*d*(ma+1)+1;
Gamma(r1:r1+d-1,c1:c1+d-1)==Gamma(r2:r2+d-1,c2:c2+d-1)';
end
end
if z==0 && y~=0 && yp~=0 && zp~=0 % Reduction of type AA'B, none = 1.
r1=y*d+1;
c1=zp*d*(ma+1)+yp*d+1;
r2=zp*d*(ma+1)+y*d+1;
c2=yp*d+1;
r3=zp*d*(ma+1)+y*d+1;
c3=zp*d*(ma+1)+yp*d+1;
Gamma(r1:r1+d-1,c1:c1+d-1)==Gamma(r2:r2+d-1,c2:c2+d-1);
Gamma(r1:r1+d-1,c1:c1+d-1)==Gamma(r3:r3+d-1,c3:c3+d-1);
end
if y~=0 && yp~=0 && zp==0 % (AB,B')=(AB',B)'
if y~=yp
r1=z*d*(ma+1)+y*d+1;
c1=yp*d+1;
r2=z*d*(ma+1)+yp*d+1;
c2=y*d+1;
Gamma(r1:r1+d-1,c1:c1+d-1)==Gamma(r2:r2+d-1,c2:c2+d-1)';
end
end
if y~=0 && z~=0 && yp~=0 && zp==0 % (AB, A'B') = (A'B',AB)'
if y~=yp && z~=zp
r1=z*d*(ma+1)+y*d+1;
c1=zp*d*(ma+1)+yp*d+1;
r2=zp*d*(ma+1)+yp*d+1;
c2=z*d*(ma+1)+y*d+1;
Gamma(r1:r1+d-1,c1:c1+d-1)==Gamma(r2:r2+d-1,c2:c2+d-1)';
end
end
end
end
end
end
cvx_end
threshold = 1e-8; % threshold for being 'inside' AQ (i.e. anything smaller will be considered inside, anything larger outside)
output=cvx_optval <= threshold;
beta = cvx_optval;
end