Unverified Commit 26cf40eb authored by dwuggh's avatar dwuggh
Browse files

fix bugs

parent 74ef7f27
......@@ -32,6 +32,14 @@ class DensityOperator(QOperator):
self.operator = op.operator
self.qubits = op.qubits
def scale(self, scalar):
self.operator = self.operator * scalar
def broadcast_with_self(self, qubits):
op = self.broadcast_with(qubits)
self.operator = op.operator
self.qubits = op.qubits
def evolution(self, unitary: QOperator, p_g = 0, on_qubits = None):
# ρ = U ρ U^\dagger
......@@ -60,32 +68,55 @@ class DensityOperator(QOperator):
self.qubits = result.qubits
return result
# bell measurement: only result 00 and 11 will be reserved
'''
In the paper, measurement error is modeled by perfect measurement preceded by
inversion of the state with probability $p_m$.
NOTE this function may cause porbability loss.
'''
def bell_measure(self, q1, q2, pauli, p_m = 0):
channel1 = measure_x(q1) if pauli == 'x' else measure_z(q2)
channel2 = measure_x(q1) if pauli == 'x' else measure_z(q2)
def bell_measure(self, q1, q2, pauli1, pauli2, p_m = 0):
channel1 = measure_x(q1) if pauli1 == 'x' else measure_z(q1)
channel2 = measure_x(q2) if pauli2 == 'x' else measure_z(q2)
# the 4 projection operator
P_00 = multiply(channel1.kraus_operators[0], channel2.kraus_operators[0])
P_01 = multiply(channel1.kraus_operators[0], channel2.kraus_operators[1])
P_10 = multiply(channel1.kraus_operators[1], channel2.kraus_operators[0])
P_11 = multiply(channel1.kraus_operators[1], channel2.kraus_operators[1])
# P_00.print()
# the 4 resulting density matrix, with weight of p_ij
ρ_00 = multiply(self, P_00).partial_trace([q1, q2])
ρ_01 = multiply(self, P_01).partial_trace([q1, q2])
ρ_10 = multiply(self, P_10).partial_trace([q1, q2])
ρ_11 = multiply(self, P_11).partial_trace([q1, q2])
p_00 = ρ_00.operator.trace()
p_01 = ρ_01.operator.trace()
p_10 = ρ_10.operator.trace()
p_11 = ρ_11.operator.trace()
p_sum = p_00 + p_11
p00 = (1 - p_m) ** 2 * p_00 + p_m ** 2 * p_11
p01 = p_m * (1 - p_m) * p_sum
p10 = p_m * (1 - p_m) * p_sum
p11 = (1 - p_m) ** 2 * p_11 + p_m ** 2 * p_00
# no rescaling: this is intentional
# p00 = p00 / p_sum
# p01 = p01 / p_sum
# p10 = p10 / p_sum
# p11 = p11 / p_sum
p1 = (1 - p_m) ** 2 + p_m ** 2
p2 = 2 * p_m * (1 - p_m)
print(p_00, p_01, p_10, p_11)
print(p00 , p01 , p10 , p11)
p1 = np.sqrt(p1)
p2 = np.sqrt(p2)
operator = ρ_00.operator * p00 + ρ_01.operator * p10 + ρ_10.operator * p10 +ρ_11.operator * p11
# this channel will cause probability loss because of post-selection
loss_channel = QChannel([P_00.scale_to(p1), P_01.scale_to(p2), P_10.scale_to(p2), P_11.scale_to(p1)])
self.channel(loss_channel)
operator = self.partial_trace([q1, q2])
self.operator = operator.operator
self.qubits = operator.qubits
self.operator = operator
self.qubits = ρ_00.qubits
......@@ -103,6 +134,10 @@ def measure_z(q = 0):
# generate bell pair in werner form
def bell_pair(p_n, qubits = [0, 1]):
Φ = np.array([1, 0, 0, 1], dtype=np.float64) / 2 ** 0.5
mat = (1 - 4 / 3 * p_n) * np.outer(Φ, Φ) + (p_n / 3) * np.identity(4)
ρ_0 = np.array([
[1, 0, 0, 1],
[0, 0, 0, 0],
[0, 0, 0, 0],
[1, 0, 0, 1]], dtype=np.float64) / 2
mat = (1 - 4 / 3 * p_n) * ρ_0 + (p_n / 3) * np.identity(4)
return DensityOperator(qubits, mat)
......@@ -27,7 +27,7 @@ class QOperator(object):
return QOperator(qubits, self.operator)
# scale this operator by a scalar.
def scale_to(self, scalar):
def scale(self, scalar):
operator = self.operator * scalar
return QOperator(self.qubits, operator)
......@@ -142,6 +142,64 @@ class QOperator(object):
mat_new[i_new, j_new] += self.operator[i, j]
return QOperator(qubits_reserved, mat_new)
'''
extract qubits at given basis('x' or 'z')
'''
def extract(self, qubit, basis, target):
# get index of qubit
index = -1
for (i, q) in enumerate(self.qubits):
if q == qubit:
index = i
break
qnum = self.qnum()
qubits_new = np.delete(self.qubits, index)
dim, _ = self.operator.shape
# get sub matrices of |0><0| and |1><1| (Z basis)
dim_new = 2 ** (qnum - 1)
mat_00 = np.zeros((dim_new, dim_new), dtype=np.float64)
mat_01 = np.zeros((dim_new, dim_new), dtype=np.float64)
mat_10 = np.zeros((dim_new, dim_new), dtype=np.float64)
mat_11 = np.zeros((dim_new, dim_new), dtype=np.float64)
for i in range(dim):
for j in range(dim):
di = np.flip(get_bin_digits(i, qnum))
dj = np.flip(get_bin_digits(j, qnum))
ql = di[index]
qr = dj[index]
di_new = np.delete(di, index)
dj_new = np.delete(dj, index)
i_new = from_bin_digits(np.flip(di_new))
j_new = from_bin_digits(np.flip(dj_new))
if ql == 0 and qr == 0:
mat_00[i_new, j_new] = self.operator[i, j]
if ql == 0 and qr == 1:
mat_01[i_new, j_new] = self.operator[i, j]
if ql == 1 and qr == 0:
mat_10[i_new, j_new] = self.operator[i, j]
if ql == 1 and qr == 1:
mat_11[i_new, j_new] = self.operator[i, j]
operator = None
if basis == 'z':
if target == 0:
operator = mat_00
else:
operator = mat_11
else:
if target == 0:
mat_pp = (mat_00 + mat_01 + mat_10 + mat_11) / 2
operator = mat_pp
else:
mat_nn = (mat_00 - mat_01 - mat_10 + mat_11) / 2
operator = mat_nn
return QOperator(qubits_new, operator)
# multiply 2 QOperators, often needs to broadcast them
def multiply(lhs: QOperator, rhs: QOperator):
......
......@@ -132,3 +132,116 @@ def stablizer_measurement(err_model: ErrorModel, stringent = True, stringent_plu
ρ1 = make_bell(err_model, stringent, stringent_plus)
ρ2 = make_bell(err_model, stringent, stringent_plus)
ρ2 = ρ2.alter_qubits([6, 9])
'''
|0>x|0>z + |1>x|1>z
'''
def BBPSSW_bell_pair(qubits, p_n = 0):
ρ_0 = np.array([
[1, 1, 1, -1],
[1, 1, 1, -1],
[1, 1, 1, -1],
[-1, -1, -1, 1]], dtype=np.float64) / 4
mat = (1 - 4 / 3 * p_n) * ρ_0 + (p_n / 3) * np.identity(4)
return DensityOperator(qubits, mat)
'''
BBPSSW protocol in 1996, Bennett et al.
q2 ------∎--X--M(z)
q1 ---∎--|--⊥------
| |
q4 ---|--∎--⊤--M(x)
q3 ---∎-----X------
'''
def BBPSSW(err_model: ErrorModel, detail = False):
# ρ1 = bell_pair(err_model.p_n, [1, 3])
# ρ2 = bell_pair(err_model.p_n, [2, 4])
ρ1 = BBPSSW_bell_pair([1, 3], err_model.p_n)
ρ2 = BBPSSW_bell_pair([2, 4], err_model.p_n)
if detail:
ρ1.print()
ρ1.merge(ρ2)
ρ1.broadcast_with_self([1, 3, 2, 4])
if detail:
ρ1.print()
c1 = cnot([1, 2])
# c1 = cnot([2, 1])
c2 = cnot([4, 3])
c2 = multiply(c1, c2)
# c2.broadcast_with([1, 3, 2, 4]).print()
# ρ1.evolution(c1, err_model.p_g)
ρ1.evolution(c2, err_model.p_g)
ρ1.broadcast_with_self([1, 3, 2, 4])
if detail:
ρ1.print()
# ρ1.bell_measure(2, 4, 'z')
probs = ρ1.operator.trace()
if detail:
print("probs", probs)
ρ1.bell_measure(2, 4, 'z', 'x')
# ρ1 = ρ1.partial_trace([2, 4])
probs = ρ1.operator.trace()
ρ1.scale(1 / probs)
if detail:
print("probs", probs)
ρ1.print()
return ρ1
'''
this circuit does nothing.
q2 ------∎--⊤--M(z)
q1 ---∎--|--X------
| |
q4 ---|--∎--⊤--M(z)
q3 ---∎-----X------
'''
def BBPSSW_2(err_model: ErrorModel, detail = False):
ρ1 = bell_pair(err_model.p_n, [1, 3])
ρ2 = bell_pair(err_model.p_n, [2, 4])
if detail:
ρ1.print()
ρ1.merge(ρ2)
ρ1.broadcast_with_self([1, 3, 2, 4])
if detail:
ρ1.print()
# c1 = cnot([1, 2])
c1 = cnot([2, 1])
c2 = cnot([4, 3])
c2 = multiply(c1, c2)
c2.broadcast_with([1, 3, 2, 4]).print()
# ρ1.evolution(c1, err_model.p_g)
ρ1.evolution(c2, err_model.p_g)
ρ1.broadcast_with_self([1, 3, 2, 4])
if detail:
ρ1.print()
# ρ1.bell_measure(2, 4, 'z')
ρ1.bell_measure(2, 4, 'z', 'z')
# ρ1 = ρ1.partial_trace([2, 4])
probs = ρ1.operator.trace()
ρ1.scale(1 / probs)
if detail:
print("probs", probs)
ρ1.print()
return ρ1
......@@ -8,7 +8,7 @@ if __name__ == "__main__":
formatter=dict(float=lambda x: "%6.4g" % x)
)
perfect_bell = cc.bell_pair(0, [0, 1])
err_model = cc.ErrorModel(0.1, 0.005, 0.005)
err_model = cc.ErrorModel(0.05, 0.005, 0.005)
noise_bell = cc.bell_pair(err_model.p_n, [0, 1])
print("no purification: ", cc.entanglement_fidelity(perfect_bell, noise_bell))
......
from calc_channel.QOperator import multiply
import calc_channel as cc
import numpy as np
'''
A: 0, 1, 2
B: 3, 4, 5
C: 6, 7, 8
D: 9, 10, 11
'''
def test_1():
err_model = cc.ErrorModel()
ρ1 = cc.bell_pair(err_model.p_n, [1, 3])
ρ2 = cc.bell_pair(err_model.p_n, [2, 4])
# ρ1.print()
# print(ρ1.operator.trace())
ρ1.merge(ρ2)
# ρ1.print()
print("trace of ρ1", ρ1.operator.trace())
# ρ1.print()
c1 =cc.cphase([2, 1])
c2 =cc.cphase([4, 3])
ρ1.evolution(c1)
ρ1.evolution(c2)
print("trace of ρ1", ρ1.operator.trace())
ρ1.bell_measure(2, 4, 'z', 'z')
ρ1.print()
return ρ1
def test_2():
x = cc.pauli_x(1)
z = cc.pauli_z(0)
a = x.broadcast_with([1, 0])
# x.print()
# z.print()
a = cc.multiply(x, z)
a.print()
a = a.broadcast_with([0, 1])
a.print()
def test_3():
ρ = cc.bell_pair(0)
ρ.bell_measure(0, 1, 'x', 'x')
ρ.print()
def test_BBPSSW_2(p_n = 0.1):
err_model = cc.ErrorModel(p_n)
# perfect_bell = cc.bell_pair([0, 1])
perfect_bell = cc.bell_pair(0, [0, 1])
noise_bell = cc.bell_pair(err_model.p_n, [0, 1])
ρ = cc.BBPSSW_2(err_model, True)
# f1 is just 1 - p_n
f1 = cc.entanglement_fidelity(noise_bell, perfect_bell)
f2 = cc.entanglement_fidelity(ρ, perfect_bell)
p_theory = (f1 ** 2 + 2 * f1 * (1 - f1) / 3 + 5 * ((1 - f1) / 3) ** 2)
f_theory = (f1 ** 2 + ((1 - f1) / 3) ** 2) / p_theory
print("no purification: ", f1)
print("after purification: ", f2)
print("fidelity in theory: ", f_theory)
print("success probability in theory: ", p_theory)
def test_BBPSSW(p_n = 0.1):
np.set_printoptions(edgeitems=16, linewidth=200,
formatter=dict(float=lambda x: "%8.4g" % x)
)
err_model = cc.ErrorModel(p_n)
perfect_bell = cc.BBPSSW_bell_pair([0, 1])
noise_bell = cc.BBPSSW_bell_pair([0, 1], err_model.p_n)
ρ = cc.BBPSSW(err_model, True)
# f1 is just 1 - p_n
f1 = cc.entanglement_fidelity(noise_bell, perfect_bell)
f2 = cc.entanglement_fidelity(ρ, perfect_bell)
p_theory = (f1 ** 2 + 2 * f1 * (1 - f1) / 3 + 5 * ((1 - f1) / 3) ** 2)
f_theory = (f1 ** 2 + ((1 - f1) / 3) ** 2) / p_theory
print("no purification: ", f1)
print("after purification: ", f2)
print("fidelity in theory: ", f_theory)
print("success probability in theory: ", p_theory)
def test_4():
a = cc.bell_pair(0, [1, 3])
# a.broadcast_with_self([1, 3, 2, 4])
# a.print()
c1 = cc.cnot([2, 1])
c2 = cc.cnot([4, 3])
c = multiply(c1, c2).broadcast_with([1, 3, 2, 4])
# c.print()
a.evolution(c)
a.print()
# 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15
array = np.array(
[[1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0], # 0
[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], # 1
[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], # 2
[0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1], # 3
[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], # 4
[0, 0, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0], # 5
[0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 0, 0], # 6
[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], # 7
[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], # 8
[0, 0, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0], # 9
[0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 0, 0], # 10
[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], # 11
[1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0], # 12
[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], # 13
[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], # 14
[0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1], # 15
]
) / 2
print("\n--------------\n")
print(array)
if __name__ == '__main__':
np.set_printoptions(edgeitems=16, linewidth=200,
formatter=dict(float=lambda x: "%8.4g" % x)
)
# test_BBPSSW(0.1)
a = cc.bell_pair(0)
test_BBPSSW(0.1)
# test_BBPSSW_2(0.1)
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment