In [3]:
import numpy as np
import matplotlib as mpl
import matplotlib.pylab as plt
import qutip as qt
import seaborn as sns
%matplotlib inline
sns.set()
sns.set_context('poster')
sns.set_style('ticks')

In [385]:
def gauss_state(x,x0,sigma):
    return np.power((2*np.pi*sigma**2),-1/4)*np.exp(-1*np.power(x-x0,2)/(4*sigma**2))
def polProbs(sx,sy,sxp,syp,psi0,theta):
    result = 0
    basis  = [qt.Qobj([[1],[0]]).unit(), qt.Qobj([[0],[1]]).unit()]
    for v in basis:
        result  += ((v.dag()*theta)*(theta.dag()*sy)*(sy.dag()*sx)*(sx.dag()*psi0)\
                   *(psi0.dag()*sxp)*(sxp.dag()*syp)*(syp.dag()*v)).norm()
    return result
def two_gauss(x,x0,x1,sigma):
    return np.power((2*np.pi*sigma**2),-1/2)*np.exp(-1*np.power(x-x0,2)/(4*sigma**2)-1*np.power(x-x1,2)/(4*sigma**2))
def prob_n(n,psi0,eps,step,d0,xspace,yspace,sigma,theta):
    sigmas_x = [qt.Qobj([[1],[1]]).unit(), qt.Qobj([[1],[-1]]).unit()]
    eig_s = [1,-1]
    sigmas_y = [qt.Qobj([[1],[1j]]).unit(), qt.Qobj([[1],[-1j]]).unit()]
    z = 0
    x = 0
    y = 0
    for (sx,esx) in zip(sigmas_x,eig_s):
        for (sxp,esxp) in zip(sigmas_x,eig_s):
            for (sy,esyp) in zip(sigmas_y,eig_s):
                for (syp,esyp) in zip(sigmas_y,eig_s):
                    prob_amp = polProbs(sx,sy,sxp,syp,psi0,theta) 
                    yvals = prob_amp*two_gauss(xspace,n*(eps*esx)+step*d0,n*(eps*esxp)+step*d0,sigma)
                    xvals = two_gauss(yspace,n*(eps*esx),n*(eps*esxp),sigma)
                    #yvals = two_gauss(yspace,n*(eps*esy+d0),n*(eps*esyp+d0),sigma)
                    x += xvals
                    y += yvals
    return np.outer(x,y)

In [441]:
%matplotlib inline
f = plt.figure(figsize=(12,8))
sigma = 4
steps = 100*np.arange(6)
eps = 0.005*sigma
d0 = 12*sigma
xspace = np.arange(-5*sigma, (len(steps))*d0-5*sigma, 0.05)
yspace = np.arange(-7.5*sigma, 10*sigma, 0.05)
psi0 = qt.Qobj([[1],[1j]]).unit()
t = 0
theta = qt.Qobj([[np.cos(t)],[np.sin(t)]]).unit()
print(theta)

z=0
for s in range(len(steps)):
    z += prob_n(steps[s],psi0,eps,s,d0,xspace,yspace,sigma,theta)
plt.imshow(z,origin="lower",extent=[-5*sigma,(len(steps))*d0,-7.5*sigma,10*sigma],cmap=mpl.cm.afmhot,vmax=np.max(z)*0.1)
plt.xlabel(r"$x(arbitrary\ units)$")
plt.ylabel(r"$y(arbitrary\ units)$")
plt.xticks([],[])
plt.text(-7,25,r"$n=0$",color='white',fontsize=16)
plt.text(43,25,r"$n=100$",color='white',fontsize=16)
plt.text(92,25,r"$n=200$",color='white',fontsize=16)
plt.text(142,25,r"$n=300$",color='white',fontsize=16)
plt.text(192.5,25,r"$n=400$",color='white',fontsize=16)
plt.text(244,25,r"$n=500$",color='white',fontsize=16)
plt.title("Measurement of Non-Commuting Operators",y=1.09)
plt.savefig("non_commuting.png",dpi=f.dpi)


Quantum object: dims = [[2], [1]], shape = (2, 1), type = ket
Qobj data =
[[ 1.]
 [ 0.]]

In [ ]:
button=1, x=151, y=311, xdata=-0.221027, ydata=19.917085
button=1, x=264, y=311, xdata=51.756272, ydata=19.917085
button=1, x=380, y=316, xdata=105.113501, ydata=22.216965
button=1, x=485, y=314, xdata=153.410992, ydata=21.297013
button=1, x=597, y=316, xdata=204.928315, ydata=22.216965
button=1, x=710, y=319, xdata=256.905615, ydata=23.596894

In [183]:
sigmas_x = [qt.Qobj([[1],[1]]).unit(), qt.Qobj([[1],[-1]]).unit()]
sigmas_y = [qt.Qobj([[1],[1j]]).unit(), qt.Qobj([[1],[-1j]]).unit()]
eig_s = [1,-1]
psi0 = qt.Qobj([[0],[1]]).unit()
for (sx,esx) in zip(sigmas_x,eig_s):
        for (sxp,esxp) in zip(sigmas_x,eig_s):
            for (sy,esyp) in zip(sigmas_y,eig_s):
                for (syp,esyp) in zip(sigmas_y,eig_s):
                    print(polProbs(sx,sy,sxp,syp,psi0))


0.25
0
0
0.25
0.25
0
0
0.25
0.25
0
0
0.25
0.25
0
0
0.25

In [80]:
type(q)


Out[80]:
qutip.qobj.Qobj

In [ ]: