<jupyter_start><jupyter_text>Studying avoided crossing for a 1 cavity-2 qubit system, with and without thermal losses 1. **Introduction**2. **Problem parameters**3. **Setting up operators, Hamiltonian's, and the initial state**4. **Demonstrating avoided crossing**    * Plotting the ramp pulse generated    * Solving the Master equation and plotting the results (without thermal losses)5. **Studying the effect of various ramp times on avoided crossing**    * { Case I } No thermal losses     * { Case II } Thermal losses    * Plotting the results6. Calculating the Fidelity and Concurrence**Author** : <NAME> (<EMAIL>)**Supervisor** : <NAME> (<EMAIL>) **Date**: 9th August 2019This script was coded as part of the Helen Edwards Summer Internship program at Fermilab. The code studies the effect of avoided crossing for loading a photon from a qubit into the cavity. This is done by generating pulses with varying ramp times, and raising one of the qubit's frequecy above the cavity.   1. IntroductionThe Jaynes-Cumming model is used to explain light-matter interaction in a system with a qubit and a single cavity mode. The Hamiltonian $H$ can be extended to describe a 2-qubit and cavity system as,$$ H = \hbar \omega_c a^{\dagger}a+ \sum_{i=1}^2\frac{1}{2}\hbar \omega_{qi}\ \sigma_i^z + \sum_{i=1}^2\frac{1}{2} \hbar g(a^{\dagger} + a)(\sigma_i^-+\sigma_i^+)$$which simplifies under the rotating-wave approximation as$$ H_{\rm RWA} = \hbar \omega_c a^{\dagger}a+ \sum_{i=1}^2\frac{1}{2}\hbar \omega_a \sigma_i^z + \sum_{i=1}^2\frac{1}{2} \ \hbar g\ (a^{\dagger}\sigma_i^- + a\  \sigma_i^+)$$where $\omega_c$ and $\omega_{qi}$ are the cavity and qubit frequencies, while $a$ and $\sigma_i^-$ are the annihalation operators for the cavity and qubit respectively. Note that $i=1,2$ represents the 2 qubits.<jupyter_code>%matplotlib inline
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d.axes3d import Axes3D
from matplotlib.ticker import LinearLocator, FormatStrFormatter
from matplotlib import cm
plt.rcParams.update({'font.size': 16})

import numpy as np
from numpy import ones,vstack
from numpy.linalg import lstsq
from math import pi
from scipy.signal import find_peaks

from time import sleep
import sys

from qutip import *<jupyter_output><empty_output><jupyter_text>2. Problem parametersHere we use $\hbar=1$; the coupling terms are redefined with a multiple of $2\pi$ before them for convinience.<jupyter_code>def generateTimePulse(tr, th):
    """
    Function that generates the pulse based on the input parameters
    @param tr :: ramp up/down time for the pulse
    @param th :: hold time for the pulse
    @return t_pulse :: np array with 4 times that define the pulse
    """
    t_pulse = [0, 0, 0, 0]
    t_pulse[0] = 0
    t_pulse[1] = tr + t_pulse[0]
    t_pulse[2] = t_pulse[1] + th
    t_pulse[3] = t_pulse[2] + tr
    print("The time pulse is: ", t_pulse)
    return t_pulse
"""------------- FREQUENCIES -----------------"""
w_q1 = 2*pi*6.5;     # Qubit 1 frequency
w_q2 = 2*pi*6.8;     # Qubit 2 frequency: range from 1-9 GHz
w_f = 2*pi*7.1;      # Resonator/ Filter frequency

"""------------- COUPLING --------------------"""
g_q1f = 2*pi*0.135   # qubit 1-fitler coupling
#g_q2f = 2*pi*0.415   # qubit 2-fitler coupling

numF = 1             # number of filters
N = 2                # number of fock states 
times = np.linspace(0,200,1500)

"""------------- DISSIPATION PARAMETERS -----"""
kappa = 5*10**-3        # cavity dissipation rate
n_th_a = 3*10**-3       # avg. no. of thermal bath excitation
r1 = 5*10**-6           # qubit relaxation rate
r2 = 1*10**-5           # qubit dephasing rate

"""------------- PULSE CONTROL PARAMETERS -----"""
tr = 0              # ramp up and ramp down times
th = 110              # hold time
t_pulse = generateTimePulse(tr, th)

# amplitude to raise pulse above cavity frequency (optional)
d = 0.25
w_top = w_f + 2*pi*d
no_ramps = 800        # number of ramps pulses sent into the Hamiltonian<jupyter_output>The time pulse is:  [0, 0, 110, 110]<jupyter_text>3. Setting up the operators, Hamiltonian's, and Initial stateFor every qubit:  **sm** $\ \rightarrow \ \hat{\sigma}^{+(-)}$ is the raising and lowering operator of the *qubit* **sz** $\ \ \rightarrow \ \sigma_z $ is the Pauli-z matrix of the *qubit* **n** $\ \ \ \rightarrow \ n$ is the number operator<jupyter_code>def numOp(m):
  """
  Computes the number operator 
  @param loweringMat :: lowering matrix operator for a system
  """
  return m.dag()*m

def rwaCoupling(m1, m2):
    return m1.dag()*m2 + m2.dag()*m1

def setXYlabel(ax, x, y, req_title, title_):
    """
    Generic function to set labels for plots
    """
    ax.set_xlabel(x)
    ax.set_ylabel(y)
    if req_title == True:
        ax.set_title(title_)
    return<jupyter_output><empty_output><jupyter_text>3.1 Operators<jupyter_code># cavity
a = tensor(destroy(N), qeye(2), qeye(2))
nc = numOp(a)

# qubit 1
sm1 = tensor(qeye(N), sigmam(), qeye(2))
sz1 = tensor(qeye(N), sigmaz(), qeye(2))
n1 = numOp(sm1)

# qubit 2
sm2 = tensor(qeye(N), qeye(2), sigmam())
sz2 = tensor(qeye(N), qeye(2), sigmaz())
n2 = numOp(sm2)
# collapse operators
c_ops = []

# cavity relaxation
rate = kappa * (1 + n_th_a)
c_ops.append(np.sqrt(rate) * a)
    
# cavity excitation

# qubit 1 relaxation
c_ops.append(np.sqrt(r1 * (1+n_th_a)) * sm1)
c_ops.append(np.sqrt(r1 * n_th_a) * sm1.dag())
c_ops.append(np.sqrt(r2) * sz1)

# qubit 2 relaxation
c_ops.append(np.sqrt(r1 * (1+n_th_a)) * sm2)
c_ops.append(np.sqrt(r1 * n_th_a) * sm2.dag())
c_ops.append(np.sqrt(r2) * sz2)<jupyter_output><empty_output><jupyter_text>3.2 Hamiltonian's and initial state<jupyter_code># Qubit Hamiltonians (Hq1+Hq2)
Hq1 =  0.5*sz1       
Hq2 =  0.5*sz2

# Filter Hamiltonians (refer formula in the Introduction)
Hf = numOp(a)

# Qubit-Filter Hamiltonian
Hqf = g_q1f*(rwaCoupling(a, sm1) + rwaCoupling(a, sm2))

# time-independent Hamiltonian (see later)
H0 = w_f*Hf + w_q2*Hq2 + Hqf
 
H = H0 + w_q1*Hq1  # Resultant Hamiltonian<jupyter_output><empty_output><jupyter_text>3.3 Initial State<jupyter_code># initial state of the system. Qubit 1: excited, Qubit 2: ground st.
psi0 = tensor(basis(N,0), basis(2,0), basis(2,1))<jupyter_output><empty_output><jupyter_text>4. Demonstrating avoided crossingIn this section the qubit frequency is raised above the cavity frequency by applying a linearly varying ramp time $t$ (ns). The signal is held for a time $T-2t$ before it is ramped down again. Tranformations on closed quantum states can be modelled by unitary operators. The combined time-dependent Hamiltonian for a system undergoing a tranformation that can be representated as,$$ H(t) = H_0 + \sum_{i=0}^n c_i(t)H_i$$where $H_0$ is called the time-independent drift Hamiltonian and $H_i$ are the control Hamiltonians with a time varying amplitude $c_i(t)$.Here we write the Hamiltonian in a function-based time dependent way. See other ways [here](http://qutip.org/docs/latest/guide/dynamics/dynamics-time.html). Here the time-dependent coefficients, $f_n(t)$ of the Hamiltonian (e.g. `wf_t, w1_t,w2_t`) are expressed using Python functions  4.1 Functions<jupyter_code>"""----------------------------------------
                  PULSE FUNCTIONS
------------------------------------------"""
def fitLine(t_pulse, i, j, w1, w2, t):
    """
    Function generates a best fit line between [x1, y1] ->[x2, y2]
    Input:
        @param t_pulse :: np array containing the 4 points parameterizing the pulse
        @param i,j :: indicies of t_pulse determining the start-stop times
        @param w1, w2 :: lower and higher frequencies of the ramp pulse
        @param t :: interable time variable
    Returns:
        @polynomial(t) :: best-fit y value at t
    """
    # compute coefficients
    coefficients = np.polyfit([t_pulse[i], t_pulse[j]], [w1, w2], 1) 
    # generate best-fit polynmial
    polynomial = np.poly1d(coefficients)    
    return polynomial(t)

def rampUp(t_pulse, w1, w2, t):
    """
    Generating a ramp up pulse
    Input:
        @param t_pulse :: np array containing the 4 points parameterizing the pulse
        @param w1, w2 :: lower and higher frequencies of the ramp pulse
        @param t :: interable time variable
    Returns:
        @param w :: int giving the y-value based on t
    """
    t0 = t_pulse[0]
    t1 = t_pulse[1]
    if t0 != t1:
        if t < t1:
            return w1 + fitLine(t_pulse, 0, 1, 0, (w2-w1), t)*(t>t0) 
        if t > t1:
            return w1 + (w2-w1)*(t>t1)
    else:
        return w1 + (w2 - w1)*(t > t1)

def rampDown(t_pulse, w1, w2, t):
    """
    Generating a ramp Down pulse
    Same as the ramp Up pulse given above only with the 
    """
    t2 = t_pulse[2]
    t3 = t_pulse[3]
    if t2 != t3:
        if t > t2:
            return w1 + fitLine(t_pulse, 2, 3, (w2-w1), 0, t)*(t>t2 and t<t3) 
        if t < t2:
            return w1 + (w2-w1)*(t<t2)
    else:
        return w1 + (w2-w1)*(t<t2)

def wq1_t(t, args=None):
    """
    Function defines the time depended co-efficent of the qubit 1
    w_q1(t) is a pulse wave going from 0 to height (w_f-w_q1) at T0_1 
    """
    return (rampUp(t_pulse, w_q1, w_top, t) + rampDown(t_pulse, w_q1, w_top, t)-w_top)

def wq1_tdown(t, args=None):
    """
    Function defines the time depended co-efficent of the qubit 1
    w_q1(t) is a pulse wave going from 0 to height (w_f-w_q1) at T0_1 
    """
    return rampDown(t_pulse, w_q1, w_top, t)

def wf_t(t, args=None):
    """
    Function defines the time depended co-efficent of the filters
    (Although, there is no frequency change of the filters with time)
    so w_f(t) = constant
    """
    return w_f


def wq2_t(t, args=None):
    """
    Function defines the time depended co-efficent of qubit 2
    (Although, there is no frequency change of the quibit 2 with time)
    so w_q2(t) = constant
    """
    return w_q2

"""---------------------------------------------
                  HAMILTONIAN FUNCTIONS
---------------------------------------------"""
def plotPulse(ax, times, t_pulse, w_q1, w_top, colorCode, label_, ramp):
    """
    Plots the required pulse
    """
    if ramp == True:
        plotting = ax.plot(times, [rampUp(t_pulse, w_q1, w_top, t)/(2*pi) for t in times], colorCode, label=label_)
    elif ramp == False:
        plotting = ax.plot(times, [rampDown(t_pulse, w_q1, w_top, t)/(2*pi) for t in times], colorCode, label=label_)
    if ramp == 'Custom':
        plotting = ax.plot(times, [(rampUp(t_pulse, w_q1, w_top, t) + rampDown(t_pulse, w_q1, w_top, t)-w_top)/(2*pi) for t in times], colorCode, label=r"$\Delta$t = %.1f ns"%(t_pulse[1]-t_pulse[0]))
    return plotting

def labelTimes(t_r, t_H):
    return r"$\Delta t = %.2f {\ \rm ns}, t_{\rm H} = %.2f {\ \rm ns}$"%(t_r, t_H)    

def plotFrequencies(ax, times, wf_t, Colour, labels_, linestyle_):
    """
    Function plots the frequencies as a function of times
    """
    ax.plot(times, np.array(list(map(wf_t, times)))/(2*pi), Colour, linewidth=2, label=labels_, linestyle=linestyle_)
    ax.legend(loc = 'center left', bbox_to_anchor = (1.0, 0.5))
    return 
    
def setLabels(ax, tr, th, plot_no):
    """
    Function sets the labels of the x-y axis in the plot below
    """
    if plot_no == 0:
        ax.set_ylabel("Frequency (GHz)", fontsize=16)
        ax.set_title(labelTimes(tr, th))
    else:
        ax.set_xlabel("Time (ns)")
        ax.set_ylabel("Occupation \n probability")
    return 

def plotProb(ax, times, component, res, Colour, labels_, linestyle_):
    """
    Function plots the occupation probabilities of the components after running mesolve
    """
    ax.plot(times, np.real(expect(component, res.states)), Colour, linewidth=1.5, label=labels_, linestyle=linestyle_)
    ax.legend(loc = 'center left', bbox_to_anchor = (1.0, 0.5))
    return<jupyter_output><empty_output><jupyter_text>4.1 Plotting the ramp pulse generatedThe figure below demonstrated how the combination of ramping up and down forms the required pulse.<jupyter_code>fig, ax = plt.subplots(1, 1, figsize=(7,5))

t_pulse1 = [t_pulse[0], t_pulse[1]+2.5, t_pulse[2]-2.5, t_pulse[3]]
t_pulse2 = [t_pulse[0], (t_pulse[3]-t_pulse[0])/2+t_pulse[0], (t_pulse[3]-t_pulse[0])/2+t_pulse[0], t_pulse[3]]

# plotting the pulses
plotPulse(ax, times, t_pulse, w_q1, w_top, 'g--', r"$\Delta$t = Ramp up", True)
plotPulse(ax, times, t_pulse, w_q1, w_top, 'b--', r"$\Delta$t = Ramp down", False)
plotPulse(ax, times, t_pulse, w_q1, w_top, 'r', ' ', 'Custom')
plotPulse(ax, times, t_pulse1, w_q1, w_top, '#03fcba', ' ', 'Custom')
plotPulse(ax, times, t_pulse2, w_q1, w_top, '#c4f2f1', ' ', 'Custom')

# guide lines
ax.axvline(x=t_pulse[0], color='#f2d4c4', linestyle='--')
ax.axvline(x=t_pulse[3], color='#f2d4c4', linestyle='--')
ax.axvline(x=t_pulse2[2], color='#f2d4c4', linestyle='--')

setXYlabel(ax, 'Time (ns)', 'Frequency (Hz)', False, '')
ax.legend(loc="upper right")
fig.tight_layout()<jupyter_output><empty_output><jupyter_text>4.2 Solving the Master equation and plotting the results (without thermal losses)<jupyter_code>opts = Options(nsteps = 50000, atol = 1e-30)

# time dependent Hamiltonian
H_t = [H0, [Hq1, wq1_t]]

# Evolving the system
res1 = mesolve(H_t, psi0, times, [], [])
fig, axes = plt.subplots(2, 1, sharex=True, figsize=(12,7))

labels_ = ["cavity", "qubit 1", "qubit 2"]
w_list = [wf_t, wq1_t, wq2_t]
colors_ = ['#b4bfbc', 'b', '#b0ed3e']
linestyle_ = ['--', '-', '-']
components_ = [nc, n1, n2]

for i in [0, 1, 2]:
    plotFrequencies(axes[0], times, w_list[i], colors_[i], labels_[i], linestyle_[i])
setLabels(axes[0], tr, th, 0)

for i in [0, 1, 2]:
    plotProb(axes[1], times, components_[i], res1, colors_[i], labels_[i], linestyle_[i])
setLabels(axes[1], tr, th, 1)

fig.tight_layout()<jupyter_output><empty_output><jupyter_text>5. Studying the effect of various ramp times on avoided crossing<jupyter_code>def showProgress(idx, n):
    """
    Function prints the progress bar for a running function
    @param idx :: iterating index
    @param n :: total number of iterating variables/ total length
    """
    j = (idx+1)/n
    sys.stdout.write('\r')
    sys.stdout.write("[%-20s] %d%%" % ('='*int(20*j), 100*j))
    sys.stdout.flush()
    sleep(0.25)
    return 

def findIndex(times, t4):
    """
    Function finds the index in the times array at required point t4
    @param times :: np array contains the times at which H is evaluated
    @param t4 :: the point at which the pulse ends
    @returns param idx_array[0] :: the index of t4 in the times array
    """
    idx_array = []
    for i, t in enumerate(times):
        if t >= t4 and t < t4+1:
            idx_array.append(i)
    return idx_array[0]

def genTimePulses(rampList):
    """
    Generates pulses with variable ramp times
    @param rampList :: List with 
    """
    ramp_vals = np.empty((0, 4))
    for dt in rampList:
        t_new = [t_pulse[0], t_pulse[1]+dt, t_pulse[2]-dt, t_pulse[3]]
        ramp_vals = np.append(ramp_vals, [t_new], axis=0)  
    return ramp_vals  

def printShape(ramp_dt_array):
    print("\nDimensions of the resultant 2D array:",   np.shape(ramp_dt_array))
    return 

# get the point after the ramp down excitation
t_idx = findIndex(times, t_pulse[3])

# generating a range of pulse with varying ramp times
rampList = np.linspace(t_pulse[1], (t_pulse[3]-t_pulse[0])/2+t_pulse[0], no_ramps)-t_pulse[1]

# generates the pulses
ramp_vals = genTimePulses(rampList)<jupyter_output><empty_output><jupyter_text>{ Case I } No thermal losses: Evaluating the excited state population at all times of the pulse. The excited state population is studied for a range of different ramp pulses.<jupyter_code>#ramp_dt_array2D = evaluateHam2D(ramp_vals, True, no_ramps, H0, Hq1, wq1_t)
no_loss = True
exp_vals = []
ramp_exp_arr = []
ramp_dt_array2D = np.empty((0, len(times)))
for i in range(no_ramps):
    t_pulse = ramp_vals[i][:]        

    # time dependent Hamiltonian
    H_t = [H0, [Hq1, wq1_t]]

    # Evolving the system with/without thermal losses
    if no_loss == True:         
        output = mesolve(H_t, psi0, times, [], [])
    else:
        output = mesolve(H_t, psi0, times, c_ops, [])

    exp_vals = np.real(expect(n1, output.states))
    exp_val = np.mean(exp_vals[t_idx:-1])
    ramp_dt_array2D = np.append(ramp_dt_array2D, [exp_vals], axis=0)
    ramp_exp_arr.append(exp_val)

    # progress bar
    showProgress(i, no_ramps)
printShape(ramp_dt_array2D)<jupyter_output>[====================] 100%
Dimensions of the resultant 2D array: (800, 1500)<jupyter_text>{ Case II } Thermal losses: Evaluating the excited state population at the end of ramp down of the pulse. The excited state population is studied for a range of different ramp pulses.<jupyter_code>no_loss = False
exp_valsi = []
ramp_exp_arri = []
ramp_dt_array2Di = np.empty((0, len(times)))
for i in range(no_ramps):
    t_pulse = ramp_vals[i][:]        

    # time dependent Hamiltonian
    H_t = [H0, [Hq1, wq1_t]]

    # Evolving the system with/without thermal losses
    if no_loss == True:         
        output = mesolve(H_t, psi0, times, [], [])
    else:
        output = mesolve(H_t, psi0, times, c_ops, [])

    exp_valsi = np.real(expect(n1, output.states))
    exp_vali = np.mean(exp_valsi[t_idx:-1])
    ramp_dt_array2Di = np.append(ramp_dt_array2Di, [exp_valsi], axis=0)
    ramp_exp_arri.append(exp_vali)

    # progress bar
    showProgress(i, no_ramps)
printShape(ramp_dt_array2Di)<jupyter_output>[====================] 100%
Dimensions of the resultant 2D array: (100, 1500)<jupyter_text>5.1 Plotting the result obtained for different ramp times without thermal losses<jupyter_code>def plotForVariousRamps(rampList, times, ramp_exp_arr, t_eval):
    """
    Plots the variation in the excitation probability as a function of times and ramp up/down times
    @param rampList :: array of times by which the ramp time is increased
    @param times :: array of times at which H is evaluated
    @praem ramp_dt_array2D :: 2D array of occupation probabilities resulting for evaluating at various ramp times
    """
    fig, ax = plt.subplots(1, 2, figsize=(11,4))
    ax[0].plot(rampList, ramp_exp_arr, 'k.-', markerfacecolor='r', markeredgecolor='r', markersize=8)
    setXYlabel(ax[0], r'Ramp times $t$ (ns)', 'Excited population', True, '%d cavity'%(numF) )

    Colors_ = ['r', 'b', 'g', '#ffd500']
    for i,j in enumerate([0, findIndex(rampList, 18.3), findIndex(rampList, 36.7), findIndex(rampList, 55)]):
        ax[1].hlines(ramp_exp_arr[j], times[t_eval], times[-1], color=Colors_[i], linewidth=2.5, label=r'$\Delta t =$ %.2f'%rampList[j])
        ax[1].legend()
    setXYlabel(ax[1], 'Times (ns)', 'Final occupation probabilty', False, 'Occupation probabilty vs times for various ramps\n' )

    fig.tight_layout()
    return 

def plot3Dramps(rampList, times, ramp_dt_array2D):
    """
    3D plot of the variation in the excitation probability as a function of times and ramp up/down times
    @param rampList :: array of times by which the ramp time is increased
    @param times :: array of times at which H is evaluated
    @praem ramp_dt_array2D :: 2D array of occupation probabilities resulting for evaluating at various ramp times
    """
    fig = plt.figure(figsize=(12,7))
    ax = fig.gca(projection='3d')

    X, Y = np.meshgrid(rampList, times)

    surf = ax.plot_surface(X, Y, np.transpose(ramp_dt_array2D), rstride=1, cstride=1, cmap=cm.gist_heat, linewidth=1, antialiased=False)

    #surf2 = ax.plot_wireframe(X, Y, np.transpose(ramp_dt_array2D), rstride=40, cstride=40, color='k', linewidth=0.5)

    # Add a color bar, axis properties
    fig.colorbar(surf, shrink=0.5, aspect=10)
    ax.set_xlabel('\nRamp times' + r'$\ \Delta t$ (ns)')
    ax.set_ylabel('\nTime (ns)')
    ax.set_zlabel('\nOccupation Probabilities');
    ax.set_title(labelTimes(tr, th))
    ax.view_init(16, 25)
    plt.show()
    return


def FourierTransformOf(rampList, ramp_exp_arr):
    """
    Function calculates the Fourier Transform of the input x-y data
    @param rampLish :: x-values e.g. array of times
    @param ramp_exp_arr :: real valued array whose FFT is calculated
    
    @returns freq_arr :: x-vales in the freqeuncy domain
             power :: Fourier transformed values of input ramp_exp_arr
    """
    # fft of ram_exp_arr
    ramp_FFT = np.fft.rfft(ramp_exp_arr)
    power = np.real(ramp_FFT)*np.real(ramp_FFT)+np.imag(ramp_FFT)*np.imag(ramp_FFT)
    
    # generating the FFT frequency array
    start_pt = 1/rampList[-1]
    freq_arr = np.linspace(start_pt, start_pt*len(power), len(power))
    return freq_arr, power

def plotFFT(ax, rampList, ramp_exp_arr):
    """
    Function finds the peaks in the FFT spectrum and plots the results
    @param rampList :: x-vales e.g. array of times
    @param ramp_exp_arr :: real valued array whose FFT is calculated
    """
    rampList_FFT, ramp_exp_arr_FFT = FourierTransformOf(rampList, ramp_exp_arr)     
    # find peaks
    peak, _ = find_peaks(ramp_exp_arr_FFT, distance=100)
    
    # plot
    ax.plot(rampList_FFT[1:], ramp_exp_arr_FFT[1:],  color='#d97829', linestyle=':', marker= '.', markersize=8)
    ax.plot(rampList_FFT[peak], ramp_exp_arr_FFT[peak], 'ro')
    setXYlabel(ax, 'Frequency (GHz)', r'$\mathcal{F}\ [n_1]:$ 1 cavity', True, '(x, y) = (%.1f, %.2f)'%(ramp_exp_arr_FFT[peak], rampList_FFT[peak]))
    fig.tight_layout()
    return ramp_exp_arr_FFT[peak], rampList_FFT[peak]

def printResults(y, x):
    print(' Power value:     ', y)
    print(' Frequency value: ', x)
    return
plotForVariousRamps(rampList, times, ramp_exp_arr, t_idx)<jupyter_output><empty_output><jupyter_text>Plotting the Fourier Transform of the above plot showing Excited population as a function of Ramp times (ns). The plot below helps to summarize the shift between slow and fast modes.<jupyter_code>fig, ax = plt.subplots(1, 2, figsize=(8,4))
br_pt = 20
yf_peak, xf_peak = plotFFT(ax[0], rampList[:findIndex(rampList, br_pt)], ramp_exp_arr[:findIndex(rampList, br_pt)])
ax[0].set_xlim(0.01, 1.5)

yf_peak1, xf_peak1 = plotFFT(ax[1], rampList[findIndex(rampList, br_pt+5):], ramp_exp_arr[findIndex(rampList, br_pt+5):])
ax[1].set_xlim(0, 0.5)

print('Small ramp times (t<%.2f):'%br_pt)
printResults(yf_peak, xf_peak)
print('\nLarge ramp tines (t>%.2f):'%(br_pt+5)) 
printResults(yf_peak1, xf_peak1)<jupyter_output>Small ramp times (t<20.00):
 Power value:      [2131.80923626]
 Frequency value:  [0.45084639]

Large ramp tines (t>25.00):
 Power value:      [43.3649052]
 Frequency value:  [0.2]<jupyter_text>3D plot summing up the above two plots.<jupyter_code>plot3Dramps(rampList, times, ramp_dt_array2D)<jupyter_output><empty_output><jupyter_text>5.2 Plotting the result obtained for different ramp times with thermal losses<jupyter_code>plotForVariousRamps(rampList, times, ramp_exp_arri, t_idx)
plot3Dramps(rampList, times, ramp_dt_array2Di)<jupyter_output><empty_output><jupyter_text>6. Calculating the Fidelity and Concurrence<jupyter_code># extract the final state from the result of the simulation
rho_final = res1.states[-1]

# trace out the resonator mode and print the two-qubit density matrix
rho_qubits = ptrace(rho_final, [1, 2])
rho_qubits
# compare to the ideal result of the sqrtiswap gate (plus phase correction) for the current initial state
rho_qubits_ideal = ket2dm(tensor(phasegate(0), phasegate(-pi/2)) * sqrtiswap() * tensor(basis(2,0), basis(2,1)))
rho_qubits_ideal
print('Fidelity = ', fidelity(rho_qubits, rho_qubits_ideal))
print('Concurrence = ', concurrence(rho_qubits))<jupyter_output>Fidelity =  0.9411114916058715
Concurrence =  0.8895488058045645
{'fraction_non_alphanumeric': 0.141912792914954, 'fraction_numerical': 0.030043645917199883, 'mean_word_length': 3.699322978892871, 'pattern_counts': {'":': 0, '<': 76, '<?xml version=': 0, '>': 78, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '17', 'dir': 'jupyter-structured-clean-dedup', 'id': '13012', 'license': 'MIT', 'n_tokens_mistral': 9023, 'n_tokens_neox': 8479, 'n_words': 2419, 'path': 'Noteboooks/09. Studying avoided crossing for 2-qubit and 1-cavity system.ipynb', 'repo_name': 'SoumyaShreeram/Qubit', 'size': '23638'}