Lamb-Dicke orders cut off if timescale in hamiltonian is too large

[neil-glikin]

Setup:

ion = Ca40([("S1/2", -1/2),
             ("D5/2", -1/2)])
set_sublevel_alias!(ion, Dict("S" => ("S1/2", -1/2),
                              "D" => ("D5/2", -1/2)))

L = Laser(k=(x̂+ẑ)/√2, ϵ=(x̂-ẑ)/√2)
chain = LinearChain(ions=[ion], com_frequencies=(x=3e6,y=3e6,z=1e6), vibrational_modes=(;z=[1]))
mode = chain.vibrational_modes.z[1]
mode.N = 2  #for easier interpretation of hamiltonian printout

T = Trap(configuration=chain, B=4e-4, Bhat=ẑ, δB=0, lasers=[L])
Efield_from_pi_time!(20e-6, T, 1, 1, ("S", "D"))  #not too fast to not cause too much carrier stark shift
L.λ = transitionwavelength(ion, ("S", "D"), T)
L.Δ = mode.ν

tscale = 1e-6
h = hamiltonian(T, timescale=tscale, lamb_dicke_order=2)

tspan = (0:1:1000)*1e-6/tscale
tout, sol = timeevolution.schroedinger_dynamic(tspan, ion["S"] ⊗ mode[0], h);

This should drive the first order motional sideband; I have explicitly made lamb_dicke_order equal to 2. I have defined a variable tscale which in theory should change nothing if it takes on a different value, since tspan is defined to scale with tscale; even though tscale will scale the matrix elements of the hamiltonian, the dynamics should be unaffected since the time will scale appropriately. With tscale = 1e-6, we get the expected results (slightly off resonance due to small Stark shift from the carrier):

Plot:

ex = expect(ionprojector(T, "D"), sol)
plt.plot(1*tout, ex)
plt.xlim(1*tout[1], 1*tout[end])
plt.ylim(0, 1)
plt.ylabel("Excitation")
plt.xlabel("Time (μs)");

And for a peek into the hamiltonian (don't know of a more elegant way to see this; the sparse matrix print is cuts off too much):

for el in h(0, 0).data
    print(el)
    print("\n")
end

0.0 + 0.0im
0.07835564875489125 + 0.0im
0.0 + 0.0im
3.287624918714476e-19 + 0.005369098945105558im
0.0 + 0.0im
-0.0002604522962494028 + 3.1896207093240934e-20im
0.07835564875489125 - 0.0im
0.0 + 0.0im
-3.287624918714476e-19 - 0.005369098945105558im
0.0 + 0.0im
-0.0002604522962494028 + 3.1896207093240934e-20im
0.0 + 0.0im
0.0 + 0.0im
-3.287624918714476e-19 + 0.005369098945105558im
0.0 + 0.0im
0.07798731358518414 + 0.0im
0.0 + 0.0im
4.649403748041756e-19 + 0.007593052545891358im
3.287624918714476e-19 - 0.005369098945105558im
0.0 + 0.0im
0.07798731358518414 - 0.0im
0.0 + 0.0im
-4.649403748041756e-19 - 0.007593052545891358im
0.0 + 0.0im
0.0 + 0.0im
-0.0002604522962494028 - 3.1896207093240934e-20im
0.0 + 0.0im
-4.649403748041756e-19 + 0.007593052545891358im
0.0 + 0.0im
0.0781714811700377 + 0.0im
-0.0002604522962494028 - 3.1896207093240934e-20im
0.0 + 0.0im
4.649403748041756e-19 - 0.007593052545891358im
0.0 + 0.0im
0.0781714811700377 - 0.0im
0.0 + 0.0im

The exact values here aren't important, I just want to keep track of what elements are nonzero.

Now, if I instead set tscale to, say, 1e-5, we get (almost) the same results, except with the time axis scaled down by a factor of 10, and the values of the hamiltonian matrix elements similarly scaled up by a factor of 10:

Plot:

Hamiltonian matrix elements:

0.0 + 0.0im
0.7835564875489125 + 0.0im
0.0 + 0.0im
3.2876249187144764e-18 + 0.05369098945105559im
0.0 + 0.0im
0.0 + 0.0im
0.7835564875489125 - 0.0im
0.0 + 0.0im
-3.2876249187144764e-18 - 0.05369098945105559im
0.0 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
-3.2876249187144764e-18 + 0.05369098945105559im
0.0 + 0.0im
0.7798731358518415 + 0.0im
0.0 + 0.0im
4.649403748041757e-18 + 0.07593052545891359im
3.2876249187144764e-18 - 0.05369098945105559im
0.0 + 0.0im
0.7798731358518415 - 0.0im
0.0 + 0.0im
-4.649403748041757e-18 - 0.07593052545891359im
0.0 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
-4.649403748041757e-18 + 0.07593052545891359im
0.0 + 0.0im
0.781714811700377 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
4.649403748041757e-18 - 0.07593052545891359im
0.0 + 0.0im
0.781714811700377 - 0.0im
0.0 + 0.0im

(Actually, some terms are now missing! These are the second order Lamb-Dicke terms, which don't affect the dynamics of driving the first order sideband which is why the plot looks the same. Back to this later.)

But now, if I make tscale = 1e-4, the first order sideband disappears.

Plot:

Hamiltonian matrix elements:

0.0 + 0.0im
7.835564875489124 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
7.835564875489124 - 0.0im
0.0 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
7.798731358518413 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
7.798731358518413 - 0.0im
0.0 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
7.817148117003769 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
0.0 + 0.0im
7.817148117003769 - 0.0im
0.0 + 0.0im

The first order sideband is now gone, and we can see in the hamiltonian that it's because some terms have vanished now. If I were to drive the carrier instead, that would be always unchanged (not putting all the plots here so there's not too much in this post).

It seems that the (effective) Lamb-Dicke order is changing depending on the value of timescale in the hamiltonian; indeed, the above matrix elements are what you would get with hamiltonian(T, timescale=1e-6, lamb_dicke_order=0) (except scaled by the expected factor of 100).

I tried this with the second order sideband and got similar results. More precisely, it seems that the 2nd order Lamb-Dicke terms are cut off if timescale is any greater than 8.66e-6, and the 1st order terms are cut off if timescale is any greater than 1.73e-5. (1e-5 is in between these two values which is why I got the results from the second plot/hamiltonian output print.) I have no idea why; this is independent of the Rabi frequency but I didn't play with other parameters. The ratio of these two values is about 5. This is probably caused by something going on in the construction of the hamiltonian itself in hamiltonians.jl.