In [88]:
%pylab inline
import sys
from sympy import MatrixSymbol, Matrix, Symbol
from sympy import S, simplify, count_ops, oo
from sympy.physics.quantum import TensorProduct
import sympy as sy
sys.path.append('../')
import Icarus
In [14]:
bench = Icarus.OpticalBench()
qd = Icarus.QuantumDot()
pcm = Icarus.PhotonCountingModule()
pcm.register_detector('D1',
pcm.Detector(
delay = 0,
efficiency = 0,
sigma = 0,
matrix = bench.jxh
)
)
pcm.register_detector('D3',
pcm.Detector(
delay = 0,
efficiency = 0,
sigma = 0,
matrix = bench.ixxh
)
)
pcm.register_detector('D4',
pcm.Detector(
delay = 0,
efficiency = 0,
sigma = 0,
matrix = bench.jxxv
)
)
pcm.register_channel('D1D3',
pcm.Channel(
1.0,
pcm.detector('D3'),
pcm.detector('D1'),
'D1D3'
)
)
pcm.register_channel('D1D4',
pcm.Channel(
1.0,
pcm.detector('D4'),
pcm.detector('D1'),
'D1D4'
)
)
In [15]:
def calculate_HWP(t1, t2):
"""
Calculates the HWP matrix, with symbolic angle.
"""
return Matrix([
[sy.cos(t1)**2 - sy.sin(t1)**2, 2*sy.cos(t1)*sy.sin(t1), 0, 0, 0, 0, 0, 0],
[2*sy.cos(t1)*sy.sin(t1), -sy.cos(t1)**2 + sy.sin(t1)**2, 0, 0, 0, 0, 0, 0],
[0, 0, sy.cos(t1)**2 - sy.sin(t1)**2, 2*sy.cos(t1)*sy.sin(t1), 0, 0, 0, 0],
[0, 0, 2*sy.cos(t1)*sy.sin(t1), -sy.cos(t1)**2 + sy.sin(t1)**2, 0, 0, 0, 0],
[0, 0, 0, 0, sy.cos(t2)**2 - sy.sin(t2)**2, 2*sy.cos(t2)*sy.sin(t2), 0, 0],
[0, 0, 0, 0, 2*sy.cos(t2)*sy.sin(t2), -sy.cos(t2)**2 + sy.sin(t2)**2, 0, 0],
[0, 0, 0, 0, 0, 0, sy.cos(t1)**2 - sy.sin(t1)**2, 2*sy.cos(t1)*sy.sin(t1)],
[0, 0, 0, 0, 0, 0, 2*sy.cos(t1)*sy.sin(t1), -sy.cos(t1)**2 + sy.sin(t1)**2]
])
def calculate_QWP(t1, t2):
"""
Calculates the QWP matrix, with symbolic angle.
"""
return Matrix([
[sy.cos(t1)**2 + 1j*sy.sin(t1)**2, (1-1j)*sy.cos(t1)*sy.sin(t1), 0, 0, 0, 0, 0, 0],
[(1-1j)*sy.cos(t1)*sy.sin(t1), 1j*sy.cos(t1)**2 + sy.sin(t1)**2, 0, 0, 0, 0, 0, 0],
[0, 0, sy.cos(t1)**2 + 1j*sy.sin(t1)**2, (1-1j)*sy.cos(t1)*sy.sin(t1), 0, 0, 0, 0],
[0, 0, (1-1j)*sy.cos(t1)*sy.sin(t1), 1j*sy.cos(t1)**2 + sy.sin(t1)**2, 0, 0, 0, 0],
[0, 0, 0, 0, sy.cos(t2)**2 + 1j*sy.sin(t2)**2, (1-1j)*sy.cos(t2)*sy.sin(t2), 0, 0],
[0, 0, 0, 0, (1-1j)*sy.cos(t2)*sy.sin(t2), 1j*sy.cos(t2)**2 + sy.sin(t2)**2, 0, 0],
[0, 0, 0, 0, 0, 0, sy.cos(t2)**2 + 1j*sy.sin(t2)**2, (1-1j)*sy.cos(t2)*sy.sin(t2)],
[0, 0, 0, 0, 0, 0, (1-1j)*sy.cos(t2)*sy.sin(t2), 1j*sy.cos(t2)**2 + sy.sin(t2)**2],
])
In [16]:
phi = Symbol('phi')
state = (1.0/np.sqrt(2.0))*(Matrix(np.kron(bench.ixh, bench.ixxh)) + phi*Matrix(np.kron(bench.ixv, bench.ixxv)))
In [74]:
alpha = Symbol('alpha')
bench.HWP = calculate_HWP(alpha, alpha).evalf(subs={alpha:np.pi})
bench.QWP = calculate_QWP(alpha, alpha).evalf(subs={alpha:np.pi})
bench.HWPHWP = TensorProduct(bench.HWP, bench.HWP)
bench.QWPQWP = TensorProduct(bench.QWP, bench.QWP)
bench.SS = Matrix(bench.SS)
bench.NBSNBS = Matrix(bench.NBSNBS)
bench.PBSPBS = Matrix(bench.PBSPBS)
In [75]:
bench.setLabMatrix('NBSNBS HWPHWP SS PBSPBS')
bench.matrix_diag = Matrix(bench.matrix)
bench.setLabMatrix('NBSNBS QWPQWP SS PBSPBS')
bench.matrix_circ = Matrix(bench.matrix)
In [76]:
D1D3 = pcm.channel('D1D3')
D1D4 = pcm.channel('D1D4')
D1D3.matrix = Matrix(D1D3.matrix)
D1D4.matrix = Matrix(D1D4.matrix)
In [77]:
p_state_diag = (bench.matrix_diag)*state
p_state_circ = (bench.matrix_circ)*state
In [78]:
diag_TT = 4*(D1D3.matrix.T*p_state_diag)[0]**2
diag_TR = 4*(D1D4.matrix.T*p_state_diag)[0]**2
circ_TT = 4*(D1D3.matrix.T*p_state_circ)[0]**2
circ_TR = 4*(D1D4.matrix.T*p_state_circ)[0]**2
In [79]:
def plot_prob_diag(p, a=0):
return np.abs(diag_TT.evalf(subs={phi:p, alpha:a})), np.abs(diag_TR.evalf(subs={phi:p, alpha:a}))
def plot_prob_circ(p, a=0):
return np.abs(circ_TT.evalf(subs={phi:p, alpha:a})), np.abs(circ_TR.evalf(subs={phi:p, alpha:a}))
def plot_prob_rect(p, a=0):
return 0.5, 0
In [80]:
t = np.linspace(0, 10, 100)
def make_phase(fss, crosstau = 1e10):
qd.FSS = 0.
qd.xlifetime = fss
return qd.generate_phase()
phases = np.array([make_phase(tt) for tt in t])
diag_probs = np.array([plot_prob_diag(p) for p in phases])
circ_probs = np.array([plot_prob_circ(p) for p in phases])
rect_probs = np.array([plot_prob_rect(p) for p in phases])
In [81]:
gdiag = (diag_probs[:,0] - diag_probs[:,1])/(diag_probs[:,0] + diag_probs[:,1])
gcirc = (circ_probs[:,0] - circ_probs[:,1])/(circ_probs[:,0] + circ_probs[:,1])
grect = (rect_probs[:,0] - rect_probs[:,1])/(rect_probs[:,0] + rect_probs[:,1])
plt.plot(t, diag_probs[:,0])
plt.plot(t, circ_probs[:,0])
#plt.plot(t, gcirc/grect)
plt.ylim([-1.05, 1.05])
Out[81]:
In [87]:
(bench.QWPQWP*state).evalf(subs={phi:1, alpha:0}) - (bench.QWPQWP*state).evalf(subs={phi:0.2, alpha:np.pi/5})
Out[87]:
In [ ]: