[1]:
%matplotlib inline
%config InlineBackend.figure_format = 'svg'

import scqubits as scq
import scqubits.utils.sweep_plotting as splot
from scqubits import HilbertSpace, InteractionTerm, ParameterSweep

import numpy as np

Composite Hilbert Spaces, QuTiP Interface

The HilbertSpace class provides data structures and methods for handling composite Hilbert spaces which may consist of multiple qubits or qubits and oscillators coupled to each other. To harness the power of QuTiP, a toolbox for studying stationary and dynamical properties of closed and open quantum systems (and much more), HilbertSpace provides a convenient interface: it generates qutip.qobj objects which are then directly handled by QuTiP.

Example: two transmons coupled to a harmonic mode

Transmon qubits can be capacitively coupled to a common harmonic mode, realized by an LC oscillator or a transmission-line resonator. The Hamiltonian describing such a composite system is given by:

\begin{equation} H=H_\text{tmon,1} + H_\text{tmon,2} + \omega_r a^\dagger a + \sum_{j=1,2}g_j n_j(a+a^\dagger), \end{equation}

where \(j=1,2\) enumerates the two transmon qubits, \(\omega_r\) is the (angular) frequency of the resonator. Furthermore, \(n_j\) is the charge number operator for qubit \(j\), and \(g_j\) is the coupling strength between qubit \(j\) and the resonator.

Create Hilbert space components

The first step consists of creating the objects describing the individual building blocks of the full Hilbert space. Here, these will be the two transmons and one oscillator:

[2]:
tmon1 = scq.Transmon(
    EJ=40.0,
    EC=0.2,
    ng=0.3,
    ncut=40,
    truncated_dim=4     # after diagonalization, we will keep 3 levels
)

tmon2 = scq.Transmon(
    EJ=15.0,
    EC=0.15,
    ng=0.0,
    ncut=30,
    truncated_dim=4
)

resonator = scq.Oscillator(
    E_osc=4.5,
    truncated_dim=4  # up to 3 photons (0,1,2,3)
)

The system objects are next grouped into a Python list, and in this form used for the initialization of a HilbertSpace object. Once created, a print call to this object outputs a summary of the composite Hilbert space.

[3]:
hilbertspace = scq.HilbertSpace([tmon1, tmon2, resonator])
print(hilbertspace)
====== HilbertSpace object ======

TRANSMON
 ———— PARAMETERS ————
EJ      : 40.0
EC      : 0.2
ng      : 0.3
ncut    : 40
truncated_dim   : 4
Hilbert space dimension : 81

TRANSMON
 ———— PARAMETERS ————
EJ      : 15.0
EC      : 0.15
ng      : 0.0
ncut    : 30
truncated_dim   : 4
Hilbert space dimension : 61

OSCILLATOR
 ———— PARAMETERS ————
E_osc   : 4.5
truncated_dim   : 4
Hilbert space dimension : 4

One useful method of the HilbertSpace class is .bare_hamiltonian(). This yields the bare Hamiltonian of the non-interacting subsystems, expressed as a qutip.Qobj:

[4]:
bare_hamiltonian = hilbertspace.bare_hamiltonian()
bare_hamiltonian
[4]:
Quantum object: dims = [[4, 4, 4], [4, 4, 4]], shape = (64, 64), type = oper, isherm = True\begin{equation*}\left(\begin{array}{*{11}c}-48.968 & 0.0 & 0.0 & 0.0 & 0.0 & \cdots & 0.0 & 0.0 & 0.0 & 0.0 & 0.0\\0.0 & -44.468 & 0.0 & 0.0 & 0.0 & \cdots & 0.0 & 0.0 & 0.0 & 0.0 & 0.0\\0.0 & 0.0 & -39.968 & 0.0 & 0.0 & \cdots & 0.0 & 0.0 & 0.0 & 0.0 & 0.0\\0.0 & 0.0 & 0.0 & -35.468 & 0.0 & \cdots & 0.0 & 0.0 & 0.0 & 0.0 & 0.0\\0.0 & 0.0 & 0.0 & 0.0 & -44.881 & \cdots & 0.0 & 0.0 & 0.0 & 0.0 & 0.0\\\vdots & \vdots & \vdots & \vdots & \vdots & \ddots & \vdots & \vdots & \vdots & \vdots & \vdots\\0.0 & 0.0 & 0.0 & 0.0 & 0.0 & \cdots & -4.723 & 0.0 & 0.0 & 0.0 & 0.0\\0.0 & 0.0 & 0.0 & 0.0 & 0.0 & \cdots & 0.0 & -14.477 & 0.0 & 0.0 & 0.0\\0.0 & 0.0 & 0.0 & 0.0 & 0.0 & \cdots & 0.0 & 0.0 & -9.977 & 0.0 & 0.0\\0.0 & 0.0 & 0.0 & 0.0 & 0.0 & \cdots & 0.0 & 0.0 & 0.0 & -5.477 & 0.0\\0.0 & 0.0 & 0.0 & 0.0 & 0.0 & \cdots & 0.0 & 0.0 & 0.0 & 0.0 & -0.977\\\end{array}\right)\end{equation*}

Set up the interaction between subsystems

The pairwise interactions between subsystems are assumed to have the general form

\(V=\sum_{i\not= j} g_{ij} A_i B_j\),

where \(g_{ij}\) parametrizes the interaction strength between subsystems \(i\) and \(j\). The operator content of the coupling is given by the two coupling operators \(A_i\), \(B_j\), which are operators in the two respective subsystems. This structure is captured by setting up an InteractionTerm object:

[5]:
g1 = 0.1  # coupling resonator-CPB1 (without charge matrix elements)
g2 = 0.2  # coupling resonator-CPB2 (without charge matrix elements)

interaction1 = InteractionTerm(
    hilbertspace = hilbertspace,
    g_strength = g1,
    op1 = tmon1.n_operator(),
    subsys1 = tmon1,
    op2 = resonator.creation_operator() + resonator.annihilation_operator(),
    subsys2 =resonator
)

interaction2 = InteractionTerm(
    hilbertspace = hilbertspace,
    g_strength = g2,
    op1 = tmon2.n_operator(),
    subsys1 = tmon2,
    op2 = resonator.creation_operator() + resonator.annihilation_operator(),
    subsys2 = resonator
)

Each InteractionTerm object is initialized by specifying 1. the Hilbert space object to which it will belong 2. the interaction strength coefficient \(g_{ij}\) 3. op1, op2: the subsystem operators \(A_i\), \(B_j\) (these should be operators within the subsystems’ respective Hilbert spaces only) 4. subsys1: the subsystem objects to which op1 and op2 belong

Note: interaction Hamiltonians of the alternative form \(V=g_{ij}A_i B_j^\dagger + g_{ij}^* A_i^\dagger B_J\) (a typical form when performing rotating-wave approximation) can be specified by setting op1 to \(A_i\) and op2 to \(B_j^\dagger\), and providing the additional keyword parameter add_hc = True.

Now, collect all interaction terms in a list, and insert into the HilbertSpace object.

[6]:
interaction_list = [interaction1, interaction2]
hilbertspace.interaction_list = interaction_list

With the interactions specified, the full Hamiltonian of the coupled system can be obtained via the method .hamiltonian(). Again, this conveniently results in a qubit.Qobj operator:

[7]:
dressed_hamiltonian = hilbertspace.hamiltonian()
dressed_hamiltonian
[7]:
Quantum object: dims = [[4, 4, 4], [4, 4, 4]], shape = (64, 64), type = oper, isherm = True\begin{equation*}\left(\begin{array}{*{11}c}-48.968 & 0.030 & 0.0 & 0.0 & 0.0 & \cdots & 0.0 & 0.0 & 0.0 & 0.0 & 0.0\\0.030 & -44.468 & 0.042 & 0.0 & 0.261 & \cdots & 0.0 & 0.0 & 0.0 & 0.0 & 0.0\\0.0 & 0.042 & -39.968 & 0.052 & 0.0 & \cdots & 0.0 & 0.0 & 0.0 & 0.0 & 0.0\\0.0 & 0.0 & 0.052 & -35.468 & 0.0 & \cdots & 0.0 & 0.0 & 0.0 & 0.0 & 0.0\\0.0 & 0.261 & 0.0 & 0.0 & -44.881 & \cdots & 0.0 & 0.0 & 0.0 & 0.0 & 0.0\\\vdots & \vdots & \vdots & \vdots & \vdots & \ddots & \vdots & \vdots & \vdots & \vdots & \vdots\\0.0 & 0.0 & 0.0 & 0.0 & 0.0 & \cdots & -4.723 & 0.0 & 0.0 & -0.749 & 0.0\\0.0 & 0.0 & 0.0 & 0.0 & 0.0 & \cdots & 0.0 & -14.477 & 0.030 & 0.0 & 0.0\\0.0 & 0.0 & 0.0 & 0.0 & 0.0 & \cdots & 0.0 & 0.030 & -9.977 & 0.042 & 0.0\\0.0 & 0.0 & 0.0 & 0.0 & 0.0 & \cdots & -0.749 & 0.0 & 0.042 & -5.477 & 0.052\\0.0 & 0.0 & 0.0 & 0.0 & 0.0 & \cdots & 0.0 & 0.0 & 0.0 & 0.052 & -0.977\\\end{array}\right)\end{equation*}

Obtaining the eigenspectrum via QuTiP

Since the Hamiltonian obtained this way is a proper qutip.qobj, all QuTiP routines are now available. In the first case, we are still making use of the scqubit HilbertSpace.eigensys() method. In the second, case, we use QuTiP’s method .eigenenergies():

[11]:
evals, evecs = hilbertspace.eigensys(evals_count=4)
print(evals)
[-48.97770317 -45.02707241 -44.36656205 -41.18438832]
[17]:
dressed_hamiltonian = hilbertspace.hamiltonian()
dressed_hamiltonian.eigenenergies()
[17]:
array([-48.97770317, -45.02707241, -44.36656205, -41.18438832,
       -41.1776098 , -40.46448065, -39.76202396, -37.45533167,
       -37.22705488, -36.65156212, -36.56694139, -35.88617024,
       -35.14255357, -33.58960343, -33.38439485, -33.12061816,
       -32.66494366, -32.04065582, -31.96284558, -30.93536847,
       -29.65538466, -29.63912664, -28.97946785, -28.85287223,
       -28.64748897, -28.09427075, -27.35745603, -26.97461415,
       -26.21586637, -25.79648462, -25.3216198 , -25.07747135,
       -24.37567532, -24.24760531, -23.66618527, -23.15199748,
       -22.2580349 , -22.06752989, -21.57303944, -21.26626732,
       -20.85288717, -20.51511596, -19.78521899, -19.1938092 ,
       -18.41228666, -17.73463904, -17.66080329, -16.93618597,
       -16.66699841, -15.89010361, -15.58181323, -14.68288545,
       -13.84725844, -13.27037775, -13.08388811, -12.34366698,
       -11.62655644, -10.308913  ,  -9.23168892,  -8.32791939,
        -8.13883089,  -5.82301882,  -4.18180351,  -0.88957289])