Skip to content

Commit 508bb8a

Browse files
committed
first version of rabi
1 parent c09ed7a commit 508bb8a

File tree

2 files changed

+276
-0
lines changed

2 files changed

+276
-0
lines changed
Lines changed: 275 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,275 @@
1+
---
2+
title: "Vacuum Rabi oscillation"
3+
author: Li-Xun Cai
4+
date: 2025-01-14 # last update (keep this comment as a reminder)
5+
6+
engine: julia
7+
---
8+
9+
Inspirations taken from this [QuTiP tutorial](https://nbviewer.org/urls/qutip.org/qutip-tutorials/tutorials-v5/time-evolution/004_rabi-oscillations.ipynb) by J.R. Johansson, P.D. Nation, and C. Staufenbiel
10+
11+
In this notebook, the usage of `QuantumToolbox.sesolve()` and `QuantumToolbox.mesolve()` will be demonstrated with Jaynes-Cummings model to observe Rabi Oscillation in the isolated case and the dissipative case. In dissipative case, a bosonic environment would interact with JC model's cavity and two-level atom.
12+
13+
## `using` Packages
14+
15+
```{julia}
16+
using QuantumToolbox
17+
using CairoMakie # For high-quality plotting
18+
```
19+
20+
## Introduction to Jaynes-Cumming model
21+
Jaynes-Cummings (JC) model, simplest quantum mechanical model for light-matter interaction, describes an atom interacting with electromagnetic field. To simplify the interaction at the time of the model proposal, JC model considered a two-level atom interacting with a single bosonic mode (or you can consider it a single-mode cavity), i.e., a two-dimensional Hilbert space interacting with a Fock space. <!-- (with finite truncation $N$) -->
22+
23+
The Hamiltonian of JC model is given by
24+
25+
$$
26+
\hat{H}_\text{tot} = \hat{H}_{\text{a}} + \hat{H}_{\text{c}} + \hat{H}_{\text{int}}
27+
$$
28+
29+
where
30+
31+
- $\hat{H}_{\text{a}} = \frac{\omega_\text{a}}{2} \hat{\sigma}_z$: Hamiltonian for atom alone
32+
- $\hat{H}_{\text{c}} = \omega_\text{c} \hat{a}^\dagger \hat{a}$: Hamiltonian for cavity alone
33+
- $\hat{H}_{\text{int}} = \Omega \left( \hat{\sigma}^\dagger + \hat{\sigma} \right) \cdot \left( \hat{a}^\dagger + \hat{a} \right)$: Interaction Hamiltonian for coherent interaction
34+
35+
with
36+
37+
- $\omega_\text{a}$: Frequency gap of the two-level atom
38+
- $\omega_\text{c}$: Frequency of the cavity's EM mode (This is not specified whether to be in resonance with atom or not.)
39+
- $\Omega$ : Coupling strength between the atom and the cavity
40+
- $\hat{\sigma}_z$ : Pauli's $z$ matrix. Equivalent to $|e\rangle\langle e| - |g\rangle\langle g|$
41+
- $\hat{a}$ : Annihilation operator in Fock space <!-- $N$-truncated -->
42+
- $\hat{\sigma}$ : Lowering operator of atom. Equivalent to $|g\rangle\langle e|$
43+
44+
By applying rotating wave approximation (RWA), the counter rotating terms ($\hat{\sigma} \cdot \hat{a}$ and its hermition conjugate) which rotate considerably faster than the others in the interaction picture is ignored, yeilding
45+
46+
$$
47+
\hat{H}_\text{tot} \approx \hat{H}_{\text{a}} + \hat{H}_{\text{c}} + \Omega \left( \hat{\sigma} \cdot \hat{a}^\dagger + \hat{\sigma}^\dagger \cdot \hat{a} \right)
48+
$$
49+
50+
### Usage of `QuantumToolbox` for JC model in general
51+
52+
Since the packages are already imported with `using` at the beginning, we hereby start by defining some parameters and system Hamiltonian.
53+
54+
```{julia}
55+
N = 15 # Fock space truncated dimension
56+
57+
ωa = 1
58+
ωc = 1 * ωa # this frequency considers cavity and atom are IN resonance
59+
ωc_ = 0.8 * ωa # this frequency considers cavity and atom are OFF resonance
60+
σz = sigmaz() ⊗ qeye(N) # order of tensor product should be consistent throughout
61+
a = qeye(2) ⊗ destroy(N)
62+
Ω = 0.05
63+
σ = sigmam() ⊗ qeye(N)
64+
65+
Ha = ωa / 2 * σz
66+
Hc = ωc * a' * a # the symbol `'` after a `QuantumObject` act as †
67+
Hc_ = ωc_ * a' * a
68+
Hint = Ω * (σ * a' + σ' * a)
69+
70+
Htot = Ha + Hc + Hint
71+
Htot_ = Ha + Hc_ + Hint
72+
73+
print(Htot)
74+
println("\n####################\n")
75+
print(Htot_)
76+
```
77+
78+
## Isolated case
79+
For the case of JC system being isolated, i.e., with no interaction with the surrounding environment, the system's time-evolution is governed solely by Schrödinger equation $\hat{H}|\psi\rangle = \partial_t|\psi(t)\rangle$. Using `QuantumToolbox.sesolve()` is ideal for pure state evolution.
80+
81+
82+
```{julia}
83+
# define initial state ψ0
84+
ψ0 = basis(2,0) ⊗ fock(N, 1)
85+
ψ0 = ψ0 / norm(ψ0)
86+
87+
tlist = 0:2:200 # a list of time points of interest
88+
89+
# define a list of operators whose expectation value dynamics exhibit Rabi oscillation
90+
eop_ls = [
91+
a' * a, # number operator of cavity Fock space
92+
σ' * σ # excited state population in atom
93+
]
94+
95+
sol = sesolve(Htot , ψ0, tlist; e_ops = eop_ls)
96+
sol_ = sesolve(Htot_, ψ0, tlist; e_ops = eop_ls)
97+
print(sol)
98+
println("\n####################\n")
99+
print(sol_)
100+
```
101+
102+
```{julia}
103+
# Create the first figure for expectation value of photon number
104+
n = real.(sol.expect[1, :])
105+
n_ = real.(sol_.expect[1, :])
106+
fign = Figure(size = (600, 350))
107+
axn = Axis(
108+
fign[1, 1],
109+
title = "⟨n⟩",
110+
xlabel = "time [1/ωa]",
111+
ylabel = "expectation value",
112+
xlabelsize = 15,
113+
ylabelsize = 15,
114+
width = 400,
115+
height = 220
116+
)
117+
lines!(axn, tlist, n, label = "resonant")
118+
lines!(axn, tlist, n_, label = "off resonant")
119+
axislegend(axn; position = :rt, labelsize = 12)
120+
display(fign);
121+
```
122+
123+
```{julia}
124+
# Create the second figure for excited state population
125+
e = real.(sol.expect[2, :])
126+
e_ = real.(sol_.expect[2, :])
127+
fige = Figure(size = (600, 350))
128+
axe = Axis(
129+
fige[1, 1],
130+
title = "⟨e⟩",
131+
xlabel = "time [1/ωa]",
132+
ylabel = "probability",
133+
xlabelsize = 15,
134+
ylabelsize = 15,
135+
width = 400,
136+
height = 220
137+
)
138+
ylims!(axe, 0, 1)
139+
lines!(axe, tlist, e, label = "resonant")
140+
lines!(axe, tlist, e_, label = "off resonant")
141+
axislegend(axe; position = :rb, labelsize = 12)
142+
display(fige);
143+
```
144+
145+
## Dissipative case
146+
147+
In contrast to isolated evolution, a factual system interacts with its surrounding environments, resulting in energy/particle exchange. We are currently interested in observing JC model's Rabi oscillation with the addition of interaction with the external EM field.
148+
149+
We start by reviewing the interaction Hamiltonians between the EM environment and atom/cavity
150+
151+
- Atom: $$
152+
\hat{H}_{\text{a}}^\text{int} = \sum_l \alpha_l \left( \hat{b}_l + \hat{b}_l^\dagger \right) \left( \hat{\sigma} + \hat{\sigma}^\dagger \right)
153+
$$
154+
- Cavity: $$
155+
\hat{H}_{\text{c}}^\text{int} = \sum_l \beta_l \left( \hat{b}_l + \hat{b}_l^\dagger \right) \left( \hat{a} + \hat{a}^\dagger \right)
156+
$$
157+
158+
where for the $l$-th mode
159+
160+
- $\alpha_l$ is the coupling strength with the atom
161+
- $\beta_l$ is the coupling strength with the cavity
162+
- $\hat{b}_l$ is the annihilation operator
163+
164+
165+
Following the RWA approach previously mentioned and the standard procedure of Born-Markovian master equation, we obtain $\kappa$, the cavity dissipation rate, and $\gamma$, the atom dissipation rate.
166+
167+
Therefore, the time evolution of the dissipative JC model can be described by the master equation
168+
169+
$$
170+
\dot{\hat{\rho}} = -\frac{i}{\hbar} [\hat{H}, \hat{\rho}] + \sum_{i = 1}^4 \mathcal{D}[\sqrt{\Gamma_i} \hat{S}_i]\left(\hat{\rho}\right)
171+
$$
172+
173+
where the $(\Gamma_i, \hat{S}_i)$ pairs are
174+
175+
1. $\kappa \cdot n(\omega_c, T), \hat{a}^\dagger$
176+
2. $\kappa \cdot (1 + n(\omega_c, T)), \hat{a}$
177+
3. $\gamma \cdot n(\omega_a, T), \hat{\sigma}^\dagger$
178+
4. $\gamma \cdot (1 + n(\omega_a, T)), \hat{\sigma}$
179+
180+
with $n(\omega, T)$ being the Bose-Einstein distribution for the EM environment and $\mathcal{D}[\hat{\mathcal{O}}]\left(\cdot\right) = \hat{\mathcal{O}} \left(\cdot\right) \hat{\mathcal{O}}^\dagger - \frac{1}{2} \{ \hat{\mathcal{O}}^\dagger \hat{\mathcal{O}}, \cdot \}$ being the Lindblad dissipator.
181+
182+
### Solve for evolutions in dissipative case
183+
184+
We can now define variables in `julia` and solve the evolution of dissipative JC model
185+
```{julia}
186+
# Bose-Einstein distribution
187+
function nb(_ω, _KT)
188+
if _KT < 0
189+
println("Temperatur should not go below absolute zero, returned 0")
190+
_n = 0.0
191+
else
192+
_n = (_KT == 0) ? 0 : (1 / (exp(_ω / _KT) - 1))
193+
end
194+
return _n
195+
end
196+
197+
# Collapse operators for interaction with the environment with variable
198+
# dissipation rates, cavity frequency, and thermal energy of the environment
199+
cop_ls(_γ, _κ, _ωc, _KT) = (
200+
√(_κ * nb(_ωc, _KT)) * a',
201+
√(_κ * (1 + nb(_ωc, _KT))) * a,
202+
√(_γ * nb(ωa, _KT)) * σ',
203+
√(_γ * (1 + nb(ωa, _KT))) * σ,
204+
)
205+
```
206+
207+
```{julia}
208+
# use the same ψ0, tlist, and eop_ls from isolated case
209+
γ = 4e-3
210+
κ = 0.005
211+
KT = 0 # for theoretical vacuum EM field
212+
213+
# `mesolve()` only has one additional keyword argument `c_ops` from `sesolve()`
214+
sol_me = mesolve(Htot, ψ0, tlist, cop_ls(γ, κ, ωc, KT), e_ops = eop_ls)
215+
sol_me_ = mesolve(Htot_, ψ0, tlist, cop_ls(γ, κ, ωc_, KT), e_ops = eop_ls)
216+
217+
print(sol_me)
218+
println("\n####################\n")
219+
print(sol_me_)
220+
```
221+
222+
```{julia}
223+
# Create the first figure for expectation value of photon number
224+
n_me = real.(sol_me.expect[1, :])
225+
n_me_ = real.(sol_me_.expect[1, :])
226+
fign_me = Figure(size = (600, 350))
227+
axn_me = Axis(
228+
fign_me[1, 1],
229+
title = "⟨n⟩",
230+
xlabel = "time [1/ωa]",
231+
ylabel = "expectation value",
232+
xlabelsize = 15,
233+
ylabelsize = 15,
234+
width = 400,
235+
height = 220
236+
)
237+
lines!(axn_me, tlist, n_me, label = "resonant")
238+
lines!(axn_me, tlist, n_me_, label = "off resonant")
239+
axislegend(axn_me; position = :rt, labelsize = 12)
240+
display(fign_me);
241+
```
242+
243+
244+
```{julia}
245+
# Create the second figure for atom excited state population
246+
e_me = real.(sol_me.expect[2, :])
247+
e_me_ = real.(sol_me_.expect[2, :])
248+
fige_me = Figure(size = (600, 350))
249+
axe_me = Axis(
250+
fige_me[1, 1],
251+
title = "⟨e⟩",
252+
xlabel = "time [1/ωa]",
253+
ylabel = "probability",
254+
xlabelsize = 15,
255+
ylabelsize = 15,
256+
width = 400,
257+
height = 220
258+
)
259+
ylims!(axe_me, 0, 1)
260+
lines!(axe_me, tlist, e_me, label = "resonant")
261+
lines!(axe_me, tlist, e_me_, label = "off resonant")
262+
axislegend(axe_me; position = :rb, labelsize = 12)
263+
display(fige_me);
264+
```
265+
266+
### Note for RWA
267+
268+
Without applying RWA, things are more interesting in the \emph{Ultrastrong Coupling} (USC) regime, where $\Omega \sim \omega_c$. In that regime, the \emph{localized} master equation (the one we used for this tutorial) fails. As that phenomenon is worth itself an independent tutorial, we will pass that to [this tutorial for `liouvillian_generalized()`]().
269+
270+
271+
272+
## Version Information
273+
```{julia}
274+
QuantumToolbox.versioninfo()
275+
```

QuantumToolbox.jl/toc.md

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -11,6 +11,7 @@ listing:
1111
date: "Last Update"
1212
contents:
1313
- "time_evolution/*.qmd"
14+
- "rabi/*.qmd"
1415
---
1516

1617
The following tutorials demonstrate and introduce specific functionality of `QuantumToolbox.jl`.

0 commit comments

Comments
 (0)