← Back to branch summary

~jrjohansson/qutip/master

« back to all changes in this revision

Viewing changes to xsqueeze.py

  • Committer: Paul Nation
  • Date: 2011-04-21 04:46:56 UTC
  • Revision ID: git-v1:dd4c966b490aa468dfbd28cef66694df4bf235c8

Show diffs side-by-side

added added

removed removed

Lines of Context:
xsqueeze.py
 
1
#This file is part of QuTIP.
 
2
#
 
3
#    QuTIP is free software: you can redistribute it and/or modify
 
4
#    it under the terms of the GNU General Public License as published by
 
5
#    the Free Software Foundation, either version 3 of the License, or
 
6
#   (at your option) any later version.
 
7
#
 
8
#    QuTIP is distributed in the hope that it will be useful,
 
9
#    but WITHOUT ANY WARRANTY; without even the implied warranty of
 
10
#    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 
11
#    GNU General Public License for more details.
 
12
#
 
13
#    You should have received a copy of the GNU General Public License
 
14
#    along with QuTIP.  If not, see <http://www.gnu.org/licenses/>.
 
15
###########################################################################
 
16
 
 
17
 
 
18
# XSQUEEZE illustrates the operator exponential and its use in making a squeezed state
 
19
#
 
20
# Derived from xsqueeze.m from the Quantum Optics toolbox by Sze M. Tan
 
21
 
 
22
 
 
23
from qutip import *
 
24
from mpl_toolkits.mplot3d import Axes3D
 
25
from matplotlib import cm
 
26
from pylab import *
 
27
 
 
28
#-----------------------------------------------------------------------------
 
29
# XSQUEEZE illustrates the operator exponential and its use in making a squeezed state
 
30
#-----------------------------------------------------------------------------
 
31
N = 20;
 
32
alpha = -1.0;     # Coherent amplitude of field
 
33
epsilon = 0.5j;   # Squeezing parameter 
 
34
a = destroy(N);
 
35
#-----------------------------------------------------------------------------
 
36
# Define displacement and squeeze operators
 
37
#-----------------------------------------------------------------------------
 
38
D = qexpm(alpha*trans(a)-conj(alpha)*a);                    # Displacement
 
39
S = qexpm(0.5*conj(epsilon)*a*a-0.5*epsilon*trans(a)*trans(a));  # Squeezing
 
40
psi = D*S*basis(N,0); # Apply to vacuum state
 
41
g = 2;
 
42
print "psi = ", psi
 
43
 
 
44
#-----------------------------------------------------------------------------
 
45
#pause # Press [Enter] to calculate Wigner function
 
46
#-----------------------------------------------------------------------------
 
47
xvec = arange(-40.,40.)*5./40
 
48
X,Y = meshgrid(xvec, xvec)
 
49
 
 
50
W=wigner(psi,xvec,xvec)
 
51
 
 
52
print "W = ", W
 
53
 
 
54
fig1 = plt.figure()
 
55
ax = Axes3D(fig1)
 
56
ax.plot_surface(X, Y, W, rstride=2, cstride=2, cmap=cm.jet, alpha=0.7)
 
57
ax.contour(X, Y, W, 15,zdir='x', offset=-6)
 
58
ax.contour(X, Y, W, 15,zdir='y', offset=6)
 
59
ax.contour(X, Y, W, 15,zdir='z', offset=-0.3)
 
60
ax.set_xlim3d(-6,6)
 
61
ax.set_xlim3d(-6,6)
 
62
ax.set_zlim3d(-0.3,0.4)
 
63
#plt.show()
 
64
 
 
65
#shading interp; 
 
66
#title('Wigner function of squeezed state');
 
67
#-----------------------------------------------------------------------------
 
68
#pause # Press [Enter] to calculate Q function
 
69
#-----------------------------------------------------------------------------
 
70
Q = qfunc(psi,xvec,xvec,g);
 
71
 
 
72
print "Q = ", Q
 
73
 
 
74
fig2 = plt.figure()
 
75
ax = Axes3D(fig2)
 
76
ax.plot_surface(X, Y, Q, rstride=2, cstride=2, cmap=cm.jet, alpha=0.7)
 
77
ax.contour(X, Y, Q,zdir='x', offset=-6)
 
78
ax.contour(X, Y, Q,zdir='y', offset=6)
 
79
ax.contour(X, Y, Q, 15,zdir='z', offset=-0.4)
 
80
ax.set_xlim3d(-6,6)
 
81
ax.set_xlim3d(-6,6)
 
82
ax.set_zlim3d(-0.3,0.4)
 
83
 
 
84
 
 
85
 
 
86
#f2 = figure(2); pcolor(xvec,yvec,real(Q));
 
87
#shading interp;
 
88
#title('Q function of squeezed state');
 
89
#-----------------------------------------------------------------------------
 
90
#pause # Press [Enter] to end demonstration
 
91
#-----------------------------------------------------------------------------
 
92
#delete(f1); delete(f2);
 
93
plt.show()
 
94
print "done..."
 
95