Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

New dipole-dipole example #27

Open
wants to merge 1 commit into
base: master
Choose a base branch
from
Open
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
265 changes: 265 additions & 0 deletions notebooks/atomic_dipole_arrays.ipynb
Original file line number Diff line number Diff line change
@@ -0,0 +1,265 @@
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Excitation transfer between dipole-coupled nanorings"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"In this example we look at the radiation pattern and transfer of collective states in rings of two-level emitters with transition frequency $\\omega_0$ following closely the work of M. Moreno-Cardoner et al, 10.1103/PhysRevA.100.023806.\n",
"\n",
"By considering only a single excitation for the system the Hamiltonian can be cast in non-Hermitian form and the dynamics follow the Schroedinger equation $\\dot{|\\Psi\\rangle} = \\frac{1}{i\\hbar} H_\\mathrm{nh} {|\\Psi\\rangle}$ with\n",
"\n",
"$$ H_\\mathrm{nh} = \\sum_{i\\neq j} (\\Omega_{ij}-\\frac{i}{2}\\Gamma_{ij}) \\sigma^-_i \\sigma^+_j$$\n",
"\n",
"For a symmetric ring, where the dipole orientations preserve the rotational symmetry (e.g. if the dipoles are\n",
"oriented perpendicularly to the plane of the ring, or tangentially to the ring), the collective modes in the\n",
"single-excitation manifold are perfect spin waves given by $|\\Psi_m\\rangle = S^+_m |g\\rangle$, with\n",
"\n",
"$$S^+_m = \\frac{1}{\\sqrt{N}} \\sum_{j=1}^N e^{i m \\phi_j} \\sigma_j^+.$$\n",
"\n",
"Here $\\phi_j = 2\\pi (j-1)/N$ is the angle associated with position $j \\ (j=1,...,N)$, and $m = 0,\\pm 1,\\pm 2,...,\\lceil \\pm (N-1)/2 \\rceil$ corresponds to the angular momentum of the mode.\n",
"\n",
"First we load the necessary libraries."
]
},
{
"cell_type": "code",
"execution_count": 98,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"true"
]
},
"execution_count": 98,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"using QuantumOptics, PyPlot, CollectiveSpins;\n",
"pygui(true)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Next we define the parameters, operators, basis and the geometry for two rings with 10 emitters each. Here we choose a tangential polarization for the dipole transitions."
]
},
{
"cell_type": "code",
"execution_count": 99,
"metadata": {},
"outputs": [],
"source": [
"N1 = 10 # number of emitters in first ring\n",
"Γ0 = 1.0 # decay rate\n",
"d = 0.1 # dipole distance\n",
"k0 = 2pi # corresponds to λ=1\n",
"\n",
"# Ring Geometry\n",
"ϕ(j) = 2π/N1*(j-1) # Angle separating atoms\n",
"Rot(φ) = [cos(φ) -sin(φ) 0;sin(φ) cos(φ) 0; 0 0 1] # Tangential polarization\n",
"R = 0.5d/sin((ϕ(2)/2))\n",
"pos = [[R*cos(ϕ(j)), R*sin(ϕ(j)), 0.0] for j=1:N1];\n",
"dips = [normalize(Rot(pi/2)*pos[i]) for i=1:N1];\n",
"\n",
"# add second ring of same size\n",
"for i = 1:N1\n",
" push!(pos,pos[i]+[0.0,2R+2d,0.0])\n",
" push!(dips,normalize(Rot(-pi/2)*pos[i]))\n",
"end\n",
"\n",
"# create N-Level basis for the whole system\n",
"N2 = length(pos)\n",
"b = NLevelBasis(N2)\n",
"sm(i) = transition(b, N2, i)\n",
"sp(i) = transition(b, i, N2)\n",
"σ = sm.([1:N2;])\n",
"σp = sp.([1:N2;]);"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The collective dipole and decay matrix elements are related to the Green's tensor $G(\\vec{r}-\\vec{r}_i,\\omega_0)$ in free space and can be generated via the CollectiveSpins.jl package. In this example we start with the most subradiant state of the first ring corresponding to $m = \\lceil \\pm (N-1)/2 \\rceil$ and let the state evolve in time.\n",
"\n",
"We also calculate the output field $E^+(\\vec{r}) = \\frac{|\\vec{\\mu}|k_0^2}{\\epsilon_0}\\sum_i G(\\vec{r}-\\vec{r}_i,\\omega_0)\\cdot \\vec{\\mu}_i\\sigma^-_i$ of the system and then the intensity $I(\\vec{r}) = \\langle E^+ E^- \\rangle$ at a point $\\vec{r}$, which is generated by all the dipoles."
]
},
{
"cell_type": "code",
"execution_count": 100,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"I_exp (generic function with 1 method)"
]
},
"execution_count": 100,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"G(i,j) = transpose(dips[i]) * GreenTensor(pos[i] - pos[j], k0) * dips[j]\n",
"\n",
"function Ω(i,j)\n",
" if i==j\n",
" return 0.0\n",
" else\n",
" return -3pi/k0*0.5*real(G(i,j))\n",
" end\n",
"end\n",
"\n",
"function Γ(i,j)\n",
" if i==j\n",
" return 1.0\n",
" else\n",
" return 6pi/k0*0.5*imag(G(i,j))\n",
" end\n",
"end\n",
"\n",
"# Eigenstates for the 1. ring\n",
"H_nh = sum((Ω(i,j)-1im*Γ(i,j)*0.5)*sp(i)*sm(j) for i=1:N1,j=1:N1) \n",
"λ = eigenstates(dense(H_nh);warning = false)\n",
"psi0 = λ[2][end] # initial state, most subradiant\n",
"\n",
"\n",
"# Schroedinger time evolution for both rings\n",
"tspan = [0.0:0.2:100;]\n",
"Heff = sum((Ω(i,j)-1im*Γ(i,j)*0.5)*sp(i)*sm(j) for i=1:N2,j=1:N2)\n",
"tout,psit = timeevolution.schroedinger(tspan,psi0,H_nh)\n",
"\n",
"function Ep(r,N) # electric field\n",
" out = [0*transition(NLevelBasis(N), N, N) for i=1:3]\n",
" for i=1:N\n",
" G_ = GreenTensor(r - pos[i], k0)\n",
" Gpi = (G_ * dips[i])\n",
" @inbounds for j=1:3\n",
" out[j].data .+= sm(i).data .* Gpi[j]\n",
" end\n",
" end\n",
" out\n",
"end\n",
"\n",
"function Intensity(r,N)\n",
" E = Ep(r,N)\n",
" sum(dagger(e)*e for e=E)\n",
"end\n",
"\n",
"I_exp(r, ψ) = real(expect(Intensity(r,N2), ψ))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Finally a grid in real space for the radiation field is defined and the instances of time are chosen such that the excitation is either almost completely in ring one or two."
]
},
{
"cell_type": "code",
"execution_count": 101,
"metadata": {
"scrolled": true
},
"outputs": [],
"source": [
"# Define grid\n",
"xmin = -0.75+R ; xmax = 0.75-R;\n",
"x = collect(range(-0.5, stop=0.5, length=61))\n",
"y = collect(range(-0.5, stop=2R+0.5, length=61))\n",
"z = 1.5R\n",
"grid = [[x1, y1, z1] for x1=x, y1=y, z1=z];\n",
"\n",
"figure(figsize=(9, 3))\n",
"kspan = [1,81,181]\n",
"for k = 1:3\n",
" t = kspan[k];t2 = tout[t]\n",
" I_profile = [I_exp(gr, psit[t]) for gr=grid];\n",
" subplot(1,3,k)\n",
" title(\"t/Γ₀ = $t2\")\n",
" contourf(y, x, I_profile, 100,cmap=\"viridis\")\n",
" for i = 1:N2\n",
" plot(pos[i][2],pos[i][1],\"wo\",markersize=6,alpha=0.5)\n",
" end\n",
" xlabel(\"x\")\n",
" ylabel(\"y\")\n",
" tight_layout()\n",
"end\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"To look more closely at the transport behaviour of the two rings we define a projector which counts only excitations in the second ring and look at its time evolution for the 5 most subradiant states starting in the first ring."
]
},
{
"cell_type": "code",
"execution_count": 102,
"metadata": {
"scrolled": false
},
"outputs": [],
"source": [
"function projector(psi)\n",
" acc = 0.0\n",
" for i = N1+1:N2\n",
" acc += expect(sp(i)*sm(i),psi)\n",
" end\n",
" return acc\n",
"end\n",
"\n",
"figure(figsize=(7, 4))\n",
"for i = 1:4\n",
" psi0 = λ[2][end+1-i]\n",
" tout,psit = timeevolution.schroedinger(tspan,psi0,Heff)\n",
" subplot(2,2,i)\n",
" plot(tout,projector.(psit))\n",
" xlabel(\"Γ₀t\")\n",
" ylabel(\"2.Ring population\")\n",
" tight_layout()\n",
" ylim(0,1)\n",
"end\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Julia 1.6.1",
"language": "julia",
"name": "julia-1.6"
},
"language_info": {
"file_extension": ".jl",
"mimetype": "application/julia",
"name": "julia",
"version": "1.6.1"
}
},
"nbformat": 4,
"nbformat_minor": 2
}