Homec4science
Diffusion BosonicCodes (master) Browse

R10106/example_calculation.py
d508aeea0430master

Last Change
Raw File

example_calculation.py
View Options

#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Sat Mar 21 11:19:05 2020
@author: ddavidsgl
"""
import numpy as np
import qutip
from scipy.optimize import minimize, shgo
""""
This is an example calculation for using some of the functionality of the
code provided in this repository.
To illustrate how the code and qutip in general is working, we try to construct
a Hamiltonian as in 10.1103/PhysRevLett.115.180501 to reproduce the results
for their parameter set
"""
#import hamiltonian and construction of operators:
from hm import *
""""setting parameters"""
# The index 0 represents the storage cavity and index 1 the cooling cavity
# Delta_s = 0 # Detuning for the storage cavity
# The detuning for the cooling cavity was obtained by optimization
A_s, A_c = 30.77, 2.31 # Anharmonicity of storage and cooling cavity
U = [-A_s, -A_c] # Anharmonicity of storage and cooling cavity
#cross-Kerr: By the construction of the Hamiltonian, only the upper triangle
#of the interaction matrix is needed.
xK = 4*np.sqrt(A_s*A_c) # cross-Kerr amplitude, should be xK approxeq 4*np.sqrt(A_s*A_c) 39.85
K = np.array([[ 0. , -xK],[ 0. , 0. ]]) # cross-Kerr interaction matrix
#This matrix is trivial for 2 modes, but can be more complex for other topologies
J = 0 # hopping
kappa_s = 1 # Dissipation rate of the storage cavity
# ALL values are expressed in units of kappa_s
# kappa_c = 26.15 # Dissipation rate of the cooling cavity
# Omega_s = kappa_c # Drive amplitudes of storage and cooling cavity
Nmodes = 2 # number of modes
# The drive amplitude of the cooling cavity was obtained through optimization
Nmax = [5, 19] #Max Fock-number for the storage and cooling cavity
"""getting list collapse operators for the dissipation process"""
# This only need to be computed once.
aa = a_tensor_list(Nmodes, Nmax) #contains a list of two qutip QObj
# the first is the annihilation op for the storage cavity the second for the
# cooling cavity
def pol(rho_s):
dims = rho_s.dims[0][0]
f1 = qutip.basis(dims, 1)
f0 = qutip.basis(dims, 0)
p1 = np.abs(qutip.expect(rho_s, f1))**2
p0 = np.abs(qutip.expect(rho_s, f0))**2
print("p1: ", p1, " p0:", p0)
return (p0-p1)/(p0 + p1)
def fidelity(rho_s):
""" calculate fidelity between single-photon Fock state and
steady state of the storage cavity"""
rho_opt = qutip.basis(Nmax[0]+1, 1)
rho_opt = qutip.ket2dm(rho_opt.unit())
return qutip.fidelity(rho_s, rho_opt)
def optimize_hamiltonian(paramlist):
""""constructing hamiltonian"""
[Delta_s, Delta_c, Omega_s, Omega_c, kappa_c] = paramlist
Delta = [Delta_s, Delta_c]
F = [Omega_s, Omega_c]
h = hamiltonian(Delta, U, K, J, F, Nmodes, Nmax, verbose=False)
""""solving for the steady-state"""
# Including the dissipation rates into the loss operators:
aa_scaled = [np.sqrt(kappa_s) * aa[0], np.sqrt(kappa_c) * aa[1]]
rho_ss = qutip.steadystate(h, aa_scaled)
rho_s = rho_ss.ptrace(0)
# or use pol(rho_s) as a minimization function
return 1 - fidelity(rho_s)
# starting optimization procedure to find optimal cooling drive amplitude
# and detuning of the cooling cavity
# here x0 is the initial guess
print("starting optimization routine...")
optimizer_args = {'options':{'disp': True}}
bounds = [(0,1*A_s), (1*xK, 3*xK), (1, 1*xK), (1, 2*xK), (10, 2*xK)]
# opt = minimize(optimize_hamiltonian, x0 = [2*A_s, 2*xK, kappa_c, kappa_c],
# callback=print, options = {'disp': True})
opt = shgo(optimize_hamiltonian, bounds, iters = 4,
minimizer_kwargs = optimizer_args,
options = {'disp' : True})
# get the optimal parameters from the minimization
[Delta_s, Delta_c, Omega_s, Omega_c, kappa_c] = opt.x
Delta = [Delta_s, Delta_c]
F = [Omega_s, Omega_c]
h_opt = hamiltonian([Delta_s, Delta_c], U, K, J, F,
Nmodes, Nmax, verbose=False)
# Including the dissipation rates into the loss operators:
aa_scaled = [np.sqrt(kappa_s) * aa[0], np.sqrt(kappa_c) * aa[1]]
rho_ss = qutip.steadystate(h_opt, aa_scaled)
#plot results
fig1, axs1 = qutip.plot_wigner_fock_distribution(rho_ss.ptrace(0), alpha_max=3,
colorbar=True)
fig2, axs2 = qutip.plot_wigner_fock_distribution(rho_ss.ptrace(1), alpha_max=5,
colorbar=True)
fig1.suptitle('storage cavity')
fig2.suptitle('cooling cavity')
fig1.tight_layout()
fig2.tight_layout()
c4science · Help