# QuTiP example: Wigner function animation for the dynamics of the Jaynes-Cumming model¶

J.R. Johansson and P.D. Nation

In [1]:
%pylab inline

Populating the interactive namespace from numpy and matplotlib

In [2]:
from mpl_toolkits.mplot3d import Axes3D
from matplotlib import cm
import time

In [3]:
from qutip import *
from qutip.ipynbtools import plot_animation

In [4]:
def jc_integrate(N, wc, wa, g, kappa, gamma, psi0, use_rwa, tlist):

# Hamiltonian
idc = qeye(N)
ida = qeye(2)

a  = tensor(destroy(N), ida)
sm = tensor(idc, destroy(2))

if use_rwa:
# use the rotating wave approxiation
H = wc * a.dag() * a + wa * sm.dag() * sm + g * (a.dag() * sm + a * sm.dag())
else:
H = wc * a.dag() * a + wa * sm.dag() * sm + g * (a.dag() + a) * (sm + sm.dag())

# collapse operators
c_op_list = []

n_th_a = 0.0 # zero temperature

rate = kappa * (1 + n_th_a)
if rate > 0.0:
c_op_list.append(sqrt(rate) * a)

rate = kappa * n_th_a
if rate > 0.0:
c_op_list.append(sqrt(rate) * a.dag())

rate = gamma
if rate > 0.0:
c_op_list.append(sqrt(rate) * sm)

# evolve and calculate return state vectors
result = mesolve(H, psi0, tlist, c_op_list, [])

return result

In [5]:
# parameters
wc = 1.0 * 2 * pi   # cavity frequency
wa = 1.0 * 2 * pi   # atom frequency
g  = 0.05 * 2 * pi  # coupling strength
kappa = 0.05        # cavity dissipation rate
gamma = 0.15        # atom dissipation rate
N = 10              # number of cavity fock states

use_rwa = True

# intial state
psi0 = tensor(basis(N,0), basis(2,1))    # start with an excited atom
#psi0 = tensor(coherent(N,1.5), basis(2,0))   # or a coherent state the in cavity
#psi0 = tensor((coherent(N,2.0)+coherent(N,-2.0)).unit(), basis(2,0))
# or a superposition of coherent states

tlist = linspace(0, 30, 150)

In [6]:
result = jc_integrate(N, wc, wa, g, kappa, gamma, psi0, use_rwa, tlist)

In [7]:
xvec = linspace(-5.,5.,100)
X,Y = meshgrid(xvec, xvec)

In [8]:
def plot_setup(result):

fig = plt.figure(figsize=(12, 6))
ax = Axes3D(fig, azim=-107, elev=49)

return fig, ax

In [9]:
cb = None

def plot_result(result, n, fig=None, axes=None):

global cb

if fig is None or axes is None:
fig, ax = plot_setup(result)

axes.cla()

# trace out the atom
rho_cavity = ptrace(result.states[n], 0)

W = wigner(rho_cavity, xvec, xvec)

surf = axes.plot_surface(X, Y, W, rstride=1, cstride=1, cmap=cm.jet,
alpha=1.0, linewidth=0.05, vmax=0.25, vmin=-0.25)
axes.set_xlim3d(-5, 5)
axes.set_ylim3d(-5, 5)
axes.set_zlim3d(-0.25, 0.25)

if not cb:
cb = plt.colorbar(surf, shrink=0.65, aspect=20)

#check_update

return fig, axes

In [10]:
plot_animation(plot_setup, plot_result, result)

---------------------------------------------------------------------------
UnicodeDecodeError                        Traceback (most recent call last)
<ipython-input-10-78c818104ba7> in <module>()
----> 1 plot_animation(plot_setup, plot_result, result)

/usr/local/lib/python3.4/dist-packages/qutip/ipynbtools.py in plot_animation(plot_setup_func, plot_func, result, name, verbose)
282     anim = animation.FuncAnimation(fig, update, frames=len(result.times), blit=True)
283
--> 284     anim.save(name + '.mp4', fps=10, writer="avconv", codec="libx264")
285
286     plt.close(fig)

/usr/lib/python3/dist-packages/matplotlib/animation.py in save(self, filename, writer, fps, dpi, codec, bitrate, extra_args, metadata, extra_anim, savefig_kwargs)
717                     #TODO: Need to see if turning off blit is really necessary
718                     anim._draw_next_frame(d, blit=False)
--> 719                 writer.grab_frame(**savefig_kwargs)
720
721         # Reconnect signal for first draw if necessary

/usr/lib/python3/dist-packages/matplotlib/animation.py in grab_frame(self, **savefig_kwargs)
203             # frame format and dpi.
204             self.fig.savefig(self._frame_sink(), format=self.frame_format,
--> 205                              dpi=self.dpi, **savefig_kwargs)
206         except RuntimeError:
207             out, err = self._proc.communicate()

/usr/lib/python3/dist-packages/matplotlib/figure.py in savefig(self, *args, **kwargs)
1420             self.set_frameon(frameon)
1421
-> 1422         self.canvas.print_figure(*args, **kwargs)
1423
1424         if frameon:

/usr/lib/python3/dist-packages/matplotlib/backend_bases.py in print_figure(self, filename, dpi, facecolor, edgecolor, orientation, format, **kwargs)
2223                 restore_bbox()
2224
-> 2225             self.figure.dpi = origDPI
2226             self.figure.set_facecolor(origfacecolor)
2227             self.figure.set_edgecolor(origedgecolor)

/usr/lib/python3/dist-packages/matplotlib/figure.py in _set_dpi(self, dpi)
383     def _set_dpi(self, dpi):
384         self._dpi = dpi
--> 385         self.dpi_scale_trans.clear().scale(dpi, dpi)
386         self.callbacks.process('dpi_changed', self)
387     dpi = property(_get_dpi, _set_dpi)

/usr/lib/python3/dist-packages/matplotlib/transforms.py in clear(self)
1784         Reset the underlying matrix to the identity transform.
1785         """
-> 1786         self._mtx = np.identity(3)
1787         self.invalidate()
1788         return self

/usr/lib/python3/dist-packages/numpy/core/numeric.py in identity(n, dtype)
2051
2052     """
-> 2053     from numpy import eye
2054     return eye(n, dtype=dtype)
2055

/usr/lib/python3.4/importlib/_bootstrap.py in _handle_fromlist(module, fromlist, import_)

UnicodeDecodeError: 'utf-8' codec can't decode byte 0xff in position 0: invalid start byte

## Versions¶

In []:
from qutip.ipynbtools import version_table

version_table()