Unverified Commit cf69f674 authored by dwuggh's avatar dwuggh
Browse files

update

parent 09360a35
from qutip import *
import numpy as np
import math
def photon_loss_error2(op,kappa):
tlists=[0.0,1.0]
H=tensor(qeye(5),qeye(5))
result=mesolve(H,op,tlists,[np.sqrt(kappa)*tensor(destroy(5),qeye(5))],[])
op1=result.states[-1]
return op1
def photon_loss_error1(op,kappa):
tlists=[0.0,1.0]
H=qeye(5)
result=mesolve(H,op,tlists,[np.sqrt(kappa)*destroy(5)],[])
op1=result.states[-1]
return op1
def kraus(kappa,l):
factor=math.pow(1-math.exp(-kappa),l)/math.factorial(l)
op=math.sqrt(factor)*(-kappa/2*num(5)).expm()*destroy(5)**l
return op
def get_bell_pair(p, kappa=0.1):
# p=0.005
roi1=tensor((basis(5,0)+basis(5,4)).unit(),(basis(5,0)+basis(5,4)).unit())+tensor(basis(5,2),basis(5,2))
roi1=roi1.unit()
logical0=(basis(5,0)+basis(5,4)).unit()
logical1=basis(5,2)
project=(basis(5,0)+basis(5,4)).unit()*((basis(5,0)+basis(5,4)).unit()).dag()+basis(5,2)*basis(5,2).dag()
projectA=project
project=tensor(project,qeye(5))
epsi=photon_loss_error1(projectA,kappa)
poi=photon_loss_error2(roi1,kappa)
E0=(-kappa/2*num(5)).expm()
E0=tensor(E0,qeye(5))
optm = Qobj(dims=epsi.dims)
for i in range(len(epsi.eigenenergies())):
if np.real(epsi.eigenenergies()[i]) > 1e-6 and np.isnan(epsi.eigenenergies()[i]) == False:
optm += 1 / np.sqrt(epsi.eigenenergies()[i]) * ket2dm(epsi.eigenstates()[1][i])
optm = tensor(optm, qeye(5))
poi1=project*E0.dag()*optm*poi*optm*E0*project
l=1
while l<4:
El=kraus(kappa,l)
El=tensor(El,qeye(5))
poi1=poi1+project*El.dag()*optm*poi*optm*El*project
l=l+1
X1=tensor(sigmax(),qeye(2))
Y1=tensor(sigmay(),qeye(2))
Z1=tensor(sigmaz(),qeye(2))
X2=tensor(qeye(2),sigmax())
Y2=tensor(qeye(2),sigmay())
Z2=tensor(qeye(2),sigmaz())
A=tensor(basis(2,0)*((basis(5,0)+basis(5,4)).unit()).dag()+basis(2,1)*basis(5,2).dag(),basis(2,0)*((basis(5,0)+basis(5,4)).unit()).dag()+basis(2,1)*basis(5,2).dag())
poi2=A*poi1*A.dag()
poi3=(1-3*p)*(1-3*p)*poi2+p*(1-3*p)*(X1*poi2*X1+Y1*poi2*Y1+Z1*poi2*Z1)+p*(1-3*p)*(X2*poi2*X2+Y2*poi2*Y2+Z2*poi2*Z2)+p*p*(X1*X2*poi2*X1*X2+Y1*X2*poi2*Y1*X2+Z1*X2*poi2*Z1*X2+X1*Y2*poi2*X1*Y2+Y1*Y2*poi2*Y1*Y2+Z1*Y2*poi2*Z1*Y2+X1*Z2*poi2*X1*Z2+Y1*Z2*poi2*Y1*Z2+Z1*Z2*poi2*Z1*Z2)
phi=tensor(basis(2,0),basis(2,0))+tensor(basis(2,1),basis(2,1))
phi=phi.unit()
# print(poi3/poi3.tr(),phi.dag()*(poi3/poi3.tr())*phi,poi1.tr())
ρ = poi3 / poi3.tr()
print(ρ)
return np.array(ρ)
if __name__ == '__main__':
ρ1 = get_bell_pair(0.01)
print(ρ1)
import calc_channel as cc
import numpy as np
import sys
from bosonic_code import get_bell_pair
def write_coeffs(name, correct, error):
with open(name + ".txt", 'w') as f:
for coeff in correct:
f.write("%.8f " % coeff)
f.write('\n')
for coeff in error:
f.write("%.8f " % coeff)
def ph_expedient(p_ph, p_local, name):
bell_operator = get_bell_pair(p_ph)
err_model = cc.ErrorModel(0.1, p_local, p_local)
ρ = cc.make_bell_with_initial(bell_operator, err_model, False, False)
ρ = cc.make_GHZ_with(ρ, err_model, False, False)
result_correct, result_error = cc.get_result(ρ.operator, err_model)
r1 = result_correct.reshape((1, -1))[0]
r2 = result_error.reshape((1, -1))[0]
write_coeffs(name, r1, r2)
return ρ
if __name__ == '__main__':
# for p_ph in np.arange(0.004, 0.0085, 0.0005):
p_local = 0.003
for p_ph in np.arange(0.000, 0.004, 0.0005):
# for p_local in [0, 0.003, 0.0045]:
name = "bosonic_channels/" + str(p_ph) + '_' + str(p_local)
print("running " + name + " ...")
ph_expedient(p_ph, p_local, name)
print("finished " + name)
# for p_ph in [0.005, 0.006]:
# for p_local in np.arange(0.003, 0.006, 0.0005):
# name = "bosonic_channels/{:.3f}_{:.4f}".format(p_ph, p_local)
# print("running {} ...".format(name))
# ph_expedient(p_ph, p_local, name)
# print("finished {}".format(name))
......@@ -140,6 +140,10 @@ qubit indexing:
def make_bell(err_model: ErrorModel, stringent = True, stringent_plus = True, step = False):
return make_bell_with(bell_pair(err_model.p_n, [0, 3]), err_model, stringent, stringent_plus)
def make_bell_with_initial(ρ: np.ndarray, err_model: ErrorModel, stringent = True, stringent_plus = True, step = False):
ρ = DensityOperator([0, 3], ρ)
return make_bell_with(ρ, err_model, stringent, stringent_plus)
def make_bell_with(ρ: DensityOperator, err_model: ErrorModel, stringent = True, stringent_plus = True, step = False):
ρ.alter_qubits([0, 3])
bell_purify_1(ρ, err_model, 0, 3, 1, 2, 4, 5, 'z')
......
import calc_channel as cc
import numpy as np
import sys
from bosonic_code import get_bell_pair
def write_coeffs(name, correct, error):
with open(name + ".txt", 'w') as f:
for coeff in correct:
f.write("%.6f " % coeff)
f.write('\n')
for coeff in error:
f.write("%.6f " % coeff)
def custom(p, stringent):
err_model = cc.ErrorModel(0.1, p, p)
......@@ -34,6 +42,19 @@ def expedient():
# ρ1.print()
ρ = cc.make_GHZ_with(ρ1, err_model, False, False)
return ρ
def ph_expedient(p_ph, p_local):
bell_operator = get_bell_pair(p_ph)
err_model = cc.ErrorModel(0.1, p_local, p_local)
ρ = cc.make_bell_with_initial(bell_operator, err_model, False, False)
ρ = cc.make_GHZ_with(ρ, err_model, False, False)
result_correct, result_error = cc.get_result(ρ.operator, err_model)
r1 = result_correct.reshape((1, -1))[0]
r2 = result_error.reshape((1, -1))[0]
name = str(p_ph) + str(p_local)
write_coeffs(name, r1, r2)
return ρ
def stringent():
err_model = cc.ErrorModel(0.1, 0.0075, 0.0075)
......@@ -58,19 +79,29 @@ if __name__ == "__main__":
print("expedient")
# print(argv[2])
p = float(argv[2]) if len(argv) > 2 else 0.006
name = 'expedient' + str(p)
ρ = custom(p, False)
elif argv[1] == 'stringent' or argv[1] == 's':
print("stringent")
# ρ = stringent()
p = float(argv[2]) if len(argv) > 2 else 0.0075
name = 'stringent' + str(p)
ρ = custom(p, True)
elif argv[1] == 'boson-expedient' or argv[1] == 'be':
p_ph = float(argv[2]) if len(argv) > 2 else 0.005
p_local = float(argv[3]) if len(argv) > 2 else 0.006
print("boson(expedient) p_ph: " + str(p_ph) + " p_local: " + str(p_local))
name = str(p_ph) + str(p_local)
ρ = ph_expedient(p_ph, p_local)
else:
print("custom")
name = 'custom'
ρ = custom(0.003, False)
ρ.print()
GHZ_perfect = np.array([1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1])
GHZ_perfect = np.outer(GHZ_perfect, GHZ_perfect) / 2
GHZ_perfect = cc.DensityOperator([0, 3, 6, 9], GHZ_perfect)
# GHZ_perfect = np.array([1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1])
# GHZ_perfect = np.outer(GHZ_perfect, GHZ_perfect) / 2
# GHZ_perfect = cc.DensityOperator([0, 3, 6, 9], GHZ_perfect)
# print("fidelity of GHZ: ", cc.entanglement_fidelity(ρ, GHZ_perfect))
......
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