QuTiP example: Landau-Zener-Stuckelberg inteferometry

In [1]:

    

%pylab inline


    

    




Populating the interactive namespace from numpy and matplotlib

In [2]:

    

from qutip import *
from qutip.gui.progressbar import TextProgressBar as ProgressBar


    

In [3]:

    

import time


    

In [4]:

    

# set up the parameters and start calculation
delta  = 1.0  * 2 * pi  # qubit sigma_x coefficient
w      = 2.0  * 2 * pi  # driving frequency
T      = 2 * pi / w     # driving period 
gamma1 = 0.00001        # relaxation rate
gamma2 = 0.005          # dephasing  rate

eps_list = linspace(-10.0, 10.0, 31) * 2 * pi
A_list   = linspace(  0.0, 10.0, 31) * 2 * pi

# pre-calculate the necessary operators
sx = sigmax(); sz = sigmaz(); sm = destroy(2); sn = num(2)

# collapse operators
c_op_list = [sqrt(gamma1) * sm, sqrt(gamma2) * sz]  # relaxation and dephasing


    

In [5]:

    

# ODE settings (for list-str format)
options = Odeoptions()
options.rhs_reuse = True
options.atol = 1e-6 # reduce accuracy to speed
options.rtol = 1e-5 # up the calculation a bit


    

In [6]:

    

# for function-callback style time-dependence
def hamiltonian_t(t, args):
    """ evaluate the hamiltonian at time t. """
    H0 = args[0]
    H1 = args[1]
    w  = args[2]
    return H0 + H1 * sin(w * t)


    

In [8]:

    

# perform the calculation for each combination of eps and A, store the result
# in a matrix
def calculate():

    p_mat = zeros((len(eps_list),len(A_list)))

    pbar = ProgressBar(len(eps_list))
    
    for m, eps in enumerate(eps_list):
        H0 = - delta/2.0 * sx - eps/2.0 * sz

        pbar.update(m)

        for n, A in enumerate(A_list):
            H1 = (A/2) * sz

            # function callback format
            #args = (H0, H1, w); H_td = hamiltonian_t

            # list-str format
            #args = {'w': w}; H_td = [H0, [H1, 'sin(w * t)']]

            # list-function format
            args = w; H_td = [H0, [H1, lambda t, w: sin(w * t)]]

            U = propagator(H_td, T, c_op_list, args, options)
            rho_ss = propagator_steadystate(U)

            p_mat[m,n] = real(expect(sn, rho_ss))

    return p_mat


    

In [9]:

    

p_mat = calculate()


    

    




Completed:  0.0%. Elapsed time:   0.01s. Est. remaining time: 00:00:00:00.
Completed: 12.9%. Elapsed time: 133.62s. Est. remaining time: 00:00:15:01.
Completed: 22.6%. Elapsed time: 222.70s. Est. remaining time: 00:00:12:43.
Completed: 32.3%. Elapsed time: 298.59s. Est. remaining time: 00:00:10:27.
Completed: 41.9%. Elapsed time: 365.58s. Est. remaining time: 00:00:08:26.
Completed: 51.6%. Elapsed time: 435.31s. Est. remaining time: 00:00:06:48.
Completed: 61.3%. Elapsed time: 503.52s. Est. remaining time: 00:00:05:18.
Completed: 71.0%. Elapsed time: 570.32s. Est. remaining time: 00:00:03:53.
Completed: 80.6%. Elapsed time: 644.29s. Est. remaining time: 00:00:02:34.
Completed: 90.3%. Elapsed time: 730.05s. Est. remaining time: 00:00:01:18.






    




/usr/local/lib/python3.3/dist-packages/qutip/propagator.py:79: UserWarning: propagator is using previously defined rhs function (options.rhs_reuse = True)
  warnings.warn(msg)

In [10]:

    

figure(figsize=(8, 8))

A_mat, eps_mat = meshgrid(A_list/(2*pi), eps_list/(2*pi))
pcolor(eps_mat, A_mat, p_mat)
xlabel(r'Bias point $\epsilon$')
ylabel(r'Amplitude $A$')
title(r'Steadystate excitation probability\n' +
      r'$H = -\frac{1}{2}\Delta\sigma_x -\frac{1}{2}\epsilon\sigma_z - \frac{1}{2}A\sin(\omega t)$\n');


    

In [11]:

    

from qutip.ipynbtools import version_table

version_table()


    

    
Out[11]:




Software
Version
QuTiP
2.3.0.dev-eece79f
Cython
0.19
SciPy
0.14.0.dev-2a4ba40
OS
posix [linux]
IPython
2.0.0-dev
Python
3.3.1 (default, Apr 17 2013, 22:30:32) 
[GCC 4.7.3]
Numpy
1.7.1
matplotlib
1.4.x
Thu Oct 03 23:05:49 2013 JST