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 style #106

Merged
merged 39 commits into from
Jan 30, 2023
Merged
Changes from 1 commit
Commits
Show all changes
39 commits
Select commit Hold shift + click to select a range
8ae98f2
Rename objects
marwahaha Jul 19, 2022
dfc22be
Format
marwahaha Jul 19, 2022
2481b5d
Remove PhysicalConstants import in species files
neil-glikin Nov 23, 2022
9e2377c
Add some convenience functions and change some function names for imp…
neil-glikin Nov 23, 2022
86a6ea3
Add functions wavelength! and wavelength_from_transition!
neil-glikin Nov 30, 2022
45a7194
Change name of starkshift/stark_shift to manualshift/manual_shift
neil-glikin Nov 30, 2022
d4e1b75
Update tests to accommodate name changes so far
neil-glikin Nov 30, 2022
55a5992
Merge branch 'new-style' into rename
neil-glikin Dec 1, 2022
2b25d7a
Merge pull request #103 from HaeffnerLab/rename
neil-glikin Dec 1, 2022
10999a9
Rename Chamber.configuration to Chamber.iontrap
neil-glikin Dec 1, 2022
3834fb9
Change default timescale from 1e-6 to 1 in hamiltonian
neil-glikin Dec 2, 2022
1e6c341
Add optional argument in LinearChain constructor for setting mode siz…
neil-glikin Jan 5, 2023
e36a307
Rename a few things related to ions
neil-glikin Jan 7, 2023
6545e65
Add getters and setters for Laser struct
neil-glikin Jan 10, 2023
2a054b3
Add getters and setters for VibrationalMode struct
neil-glikin Jan 10, 2023
f55b111
Add getters for IonTrap struct
neil-glikin Jan 10, 2023
9c2c09a
Add getters and setters for Chamber struct
neil-glikin Jan 10, 2023
80541d6
Small tweaks to Ion struct getters
neil-glikin Jan 11, 2023
0adb40d
Remove redundant basis(Chamber) function
neil-glikin Jan 11, 2023
e904150
Change some functions' names and/or arguments
neil-glikin Jan 17, 2023
7360b17
Fix chamber basis problem that affected last 2 commits
neil-glikin Jan 17, 2023
531ac4d
Change more function names
neil-glikin Jan 18, 2023
26a55f2
Naming changes to IonTraps
neil-glikin Jan 19, 2023
8ac253b
Change name of iontraps.jl in IonSim.jl
neil-glikin Jan 19, 2023
cd5e861
Update internal source code to follow new style
neil-glikin Jan 19, 2023
92fa1b1
Allow Ion constructor to specify sublevel aliases directly
neil-glikin Jan 20, 2023
0411d43
Fix Ion print
neil-glikin Jan 20, 2023
862fc2c
Update Ion print to include sublevel aliases
neil-glikin Jan 20, 2023
f2351c8
Rename files
neil-glikin Jan 20, 2023
32d75e6
Change filenames in IonSim.jl
neil-glikin Jan 20, 2023
feab0e0
Change Laser strength paramter from Efield to intensity
neil-glikin Jan 20, 2023
9d45fee
Fix efield(Laser) function
neil-glikin Jan 21, 2023
e725234
Change behavior of ionstate(IonTrap) and ionstate(Chamber)
neil-glikin Jan 24, 2023
382ce9c
Change a small amount of documentation for tutorial notebook
neil-glikin Jan 24, 2023
5dbf61b
Fix dynamics tests for change in ionstate function
neil-glikin Jan 24, 2023
bd56fe1
Fix hamiltonians test
neil-glikin Jan 24, 2023
9c7a6dd
Update hamiltonians tests
neil-glikin Jan 24, 2023
1b9ac5b
Add basis(IonTrap) function
neil-glikin Jan 25, 2023
f6f9e1c
Change a few more things in hamiltonians.jl for new style
neil-glikin Jan 25, 2023
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
Prev Previous commit
Next Next commit
Update internal source code to follow new style
neil-glikin committed Jan 19, 2023
commit cd5e8613767bfa1214d33f7928fa6a43e8ad71dd
14 changes: 7 additions & 7 deletions src/chambers.jl
Original file line number Diff line number Diff line change
@@ -285,7 +285,7 @@ function efield_from_pitime(
transition::Tuple,
Bhat::NamedTuple{(:x, :y, :z)}
)
p = laser.pointing
p = pointing(laser)
s_indx = findall(x -> x[1] == ionnumber(ion), p)
@assert length(s_indx) > 0 "This laser doesn't shine on this ion"
s = p[s_indx[1]][2]
@@ -486,7 +486,7 @@ Retuns the value of the magnetic field in `T` at the location of `ion`, includin
"""
function bfield(chamber::Chamber, ion::Ion)
@assert ionintrap(chamber, ion) "trap does not contain ion"
return chamber.B + chamber.∇B * ionposition(ion)
return bfield(chamber) + bgradient(chamber) * ionposition(ion)
end
function bfield(chamber::Chamber, ion_index::Int)
return bfield(chamber, ions(chamber)[ion_index])
@@ -604,8 +604,8 @@ ions, which are assumed to be of the same species. `ion_indxs` refer to the
ordering of the ions in the chain.
"""
function bgradient!(T::Chamber, ion_indxs::Tuple{Int, Int}, transition::Tuple, d::Real)
ionA = T.iontrap.ions[ion_indxs[1]]
ionB = T.iontrap.ions[ion_indxs[2]]
ionA = ions(T)[ion_indxs[1]]
ionB = ions(T)[ion_indxs[2]]
separation = abs(ionposition(ionA) - ionposition(ionB))

(SL1, SL2) = transition
@@ -618,7 +618,7 @@ function bgradient!(T::Chamber, ion_indxs::Tuple{Int, Int}, transition::Tuple, d
# Calculate Zeeman shifts with a unit B-field using a method of zeemanshift that ensures a nonlinear term is not used
E1 = zeemanshift(1.0, g1, m1)
E2 = zeemanshift(1.0, g2, m2)
return T.∇B = d / (abs(E2 - E1) * separation)
bgradient!(T, d / (abs(E2 - E1) * separation))
end

# In QunatumOptics.jl, this method will return true whenever the shapes of b1 and b2 match,
@@ -653,7 +653,7 @@ function lambdicke(V::VibrationalMode, I::Ion, L::Laser; scaled = false)
k = 2π / L.λ
scaled ? ν = 1 : ν = V.ν
x0 = √(ħ / (2 * mass(I) * 2π * ν))
cosθ = ndot(L.k, V.axis)
k * x0 * cosθ * V.modestructure[ionnumber(I)]
cosθ = ndot(wavevector(L), axis(V))
k * x0 * cosθ * modestructure(V)[ionnumber(I)]
end
end
84 changes: 42 additions & 42 deletions src/hamiltonians.jl
Original file line number Diff line number Diff line change
@@ -187,12 +187,12 @@ function _setup_base_hamiltonian(
rwa_cutoff *= timescale
allmodes = reverse(modes(T))
L = length(allmodes)
νlist = Tuple([mode for mode in allmodes])
mode_dims = [mode.N + 1 for mode in allmodes]
νlist = Tuple([frequency(mode) for mode in allmodes])
mode_dims = [modecutoff(mode) + 1 for mode in allmodes]

ions = reverse(T.iontrap.ions)
Q = prod([ion.shape[1] for ion in ions])
ion_arrays = [spdiagm(0 => [true for _ in 1:ion.shape[1]]) for ion in ions]
all_ions = reverse(ions(T))
Q = prod([shape(ion)[1] for ion in all_ions])
ion_arrays = [spdiagm(0 => [true for _ in 1:shape(ion)[1]]) for ion in all_ions]

ηm, Δm, Ωm = _ηmatrix(T), _Δmatrix(T, timescale), _Ωmatrix(T, timescale)
lamb_dicke_order = _check_lamb_dicke_order(lamb_dicke_order, L)
@@ -209,17 +209,17 @@ function _setup_base_hamiltonian(
local ts, ion_rows, ion_cols, ion_idxs, ion_reps, rn

# iterate over ions and lasers
for n in eachindex(ions), m in eachindex(T.lasers)
for n in eachindex(all_ions), m in eachindex(T.lasers)
if m ≡ 1
rn = length(ions) - n + 1
ts = subleveltransitions(ions[n])
rn = length(all_ions) - n + 1
ts = subleveltransitions(all_ions[n])
C = sum([
i * real.(sigma(ions[n], reverse(ts[i])...).data) for i in 1:length(ts)
i * real.(sigma(all_ions[n], reverse(ts[i])...).data) for i in 1:length(ts)
])
if length(ions) == 1
if length(all_ions) == 1
K = C
else
K = kron(ion_arrays[1:(n - 1)]..., C, ion_arrays[(n + 1):length(ions)]...)
K = kron(ion_arrays[1:(n - 1)]..., C, ion_arrays[(n + 1):length(all_ions)]...)
end
ion_rows, ion_cols, ion_vals = findnz(K)
ion_idxs = sortperm(real.(ion_vals))
@@ -444,24 +444,24 @@ end
# the time-dependence of δB. When the system is integrated, the Hamiltonian terms will be
# updated at each time step by by bfunc(t) times the individual susceptibilities.
function _setup_global_B_hamiltonian(T, timescale)
ions = T.iontrap.ions
all_ions = ions(T)
global_B_indices = Vector{Vector{Int64}}(undef, 0)
global_B_scales = Vector{Float64}(undef, 0)
δB = T.δB
δB = bfield_fluctuation(T)
τ = timescale
bfunc = FunctionWrapper{Float64, Tuple{Float64}}(t -> 2π * δB(t * τ))
if T._cnst_δB && δB(0) == 0
return global_B_indices, global_B_scales, bfunc
end
for n in eachindex(ions)
for sublevel in sublevels(ions[n])
ion_op = sigma(ions[n], sublevel)
for n in eachindex(all_ions)
for sublevel in sublevels(all_ions[n])
ion_op = sigma(all_ions[n], sublevel)
A = embed(basis(T), [n], [ion_op]).data
indices = [x[1] for x in getfield.(findall(x -> x .== complex(1, 0), A), :I)]
push!(global_B_indices, indices)
# zeemanshift(ions[n], sublevel, 1]) is the Zeeman shift of
# sublevel in units of δB.
push!(global_B_scales, τ * zeemanshift(ions[n], sublevel, 1))
push!(global_B_scales, τ * zeemanshift(all_ions[n], sublevel, 1))
end
end
return global_B_indices, global_B_scales, bfunc
@@ -493,15 +493,15 @@ end
# A 3D array of Lamb-Dicke parameters for each combination of ion, laser and mode. Modes are
# populated in reverse order.
function _ηmatrix(T)
ions = T.iontrap.ions
all_ions = ions(T)
vms = modes(T)
lasers = T.lasers
(N, M, L) = map(x -> length(x), [ions, lasers, vms])
all_lasers = lasers(T)
(N, M, L) = map(x -> length(x), [all_ions, all_lasers, vms])
ηnml = Array{Any}(undef, N, M, L)
for n in 1:N, m in 1:M, l in 1:L
δν = vms[l].δν
ν = vms[l]
eta = lambdicke(vms[l], ions[n], lasers[m], scaled = true)
δν = frequency_fluctuation(vms[l])
ν = frequency(vms[l])
eta = lambdicke(vms[l], all_ions[n], all_lasers[m], scaled = true)
if eta == 0
ηnml[n, m, L - l + 1] = 0
else
@@ -517,18 +517,18 @@ end
# for each ion transition. We need to separate this calculation from _Ωmatrix to implement
# RWA easily.
function _Δmatrix(T, timescale)
ions = T.iontrap.ions
lasers = T.lasers
(N, M) = length(ions), length(lasers)
B = T.B
∇B = T.∇B
all_ions = ions(T)
all_lasers = lasers(T)
(N, M) = length(all_ions), length(all_lasers)
B = bfield(T)
∇B = bgradient(T)
Δnmkj = Array{Vector}(undef, N, M)
for n in 1:N, m in 1:M
Btot = B + ∇B * ionposition(ions[n])
Btot = bfield(T, all_ions[n])
v = Vector{Float64}(undef, 0)
for transition in subleveltransitions(ions[n])
ωa = transitionfrequency(ions[n], transition, B = Btot)
push!(v, 2π * timescale * ((c / lasers[m].λ) + lasers[m] - ωa))
for transition in subleveltransitions(all_ions[n])
ωa = transitionfrequency(all_ions[n], transition, B = Btot)
push!(v, 2π * timescale * ((c / wavelength(all_lasers[m])) + detuning(all_lasers[m]) - ωa))
end
Δnmkj[n, m] = v
end
@@ -539,35 +539,35 @@ end
# respectively. For each row/column we have a vector of coupling strengths between the laser
# and all allowed electronic ion transitions.
function _Ωmatrix(T, timescale)
ions = T.iontrap.ions
lasers = T.lasers
(N, M) = length(ions), length(lasers)
all_ions = ions(T)
all_lasers = lasers(T)
(N, M) = length(all_ions), length(all_lasers)
Ωnmkj = Array{Vector}(undef, N, M)
for n in 1:N, m in 1:M
E = lasers[m].E
phase = lasers[m]
transitions = subleveltransitions(ions[n])
s_indx = findall(x -> x[1] == n, lasers[m].pointing)
E = efield(all_lasers[m])
ϕ = phase(all_lasers[m])
transitions = subleveltransitions(all_ions[n])
s_indx = findall(x -> x[1] == n, all_lasers[m].pointing)
if length(s_indx) == 0
Ωnmkj[n, m] = [0 for _ in 1:length(transitions)]
continue
else
s = lasers[m].pointing[s_indx[1]][2]
s = pointing(all_lasers[m])[s_indx[1]][2]
end
v = []
for t in transitions
Ω0 =
2π *
timescale *
s *
matrixelement(ions[n], t, 1.0, lasers[m].ϵ, lasers[m].k, T.Bhat) / 2.0
matrixelement(all_ions[n], t, 1.0, polarization(all_lasers[m]), wavevector(all_lasers[m]), bfield_unitvector(T)) / 2.0
if Ω0 == 0
push!(v, 0)
else
push!(
v,
FunctionWrapper{ComplexF64, Tuple{Float64}}(
t -> Ω0 * E(t) * exp(-im * phase(t))
t -> Ω0 * E(t) * exp(-im * ϕ(t))
)
)
end
2 changes: 1 addition & 1 deletion src/ions.jl
Original file line number Diff line number Diff line change
@@ -172,7 +172,7 @@ Sets the manual shift of all sublevels of `I` to zero.
"""
function zeromanualshift!(I::Ion)
for sublevel in keys(manualshift(I))
I.manualshift[sublevel] = 0.0
manualshift!(I, sublevel, 0.0)
end
end

6 changes: 3 additions & 3 deletions src/iontraps.jl
Original file line number Diff line number Diff line change
@@ -224,17 +224,17 @@ end
xmodes(lc::LinearChain)
Returns an array of all of the selected `VibrationalModes` in the x-direction in the `LinearChain`.
"""
xmodes(lc::LinearChain) = lc.selectedmodes.x
xmodes(lc::LinearChain) = selectedmodes(lc).x
"""
ymodes(lc::LinearChain)
Returns an array of all of the selected `VibrationalModes` in the y-direction in the `LinearChain`.
"""
ymodes(lc::LinearChain) = lc.selectedmodes.y
ymodes(lc::LinearChain) = selectedmodes(lc).y
"""
zmodes(lc::LinearChain)
Returns an array of all of the selected `VibrationalModes` in the z-direction in the `LinearChain`.
"""
zmodes(lc::LinearChain) = lc.selectedmodes.z
zmodes(lc::LinearChain) = selectedmodes(lc).z

"""
modecutoff!(lc::LinearChain, N::Int)
40 changes: 20 additions & 20 deletions src/operators.jl
Original file line number Diff line number Diff line change
@@ -22,7 +22,7 @@ export create,
create(v::VibrationalMode)
returns the creation operator for `v` such that: `create(v) * v[i] = √(i+1) * v[i+1]`.
"""
create(v::VibrationalMode) = SparseOperator(v, diagm(-1 => sqrt.(1:(v.N))))
create(v::VibrationalMode) = SparseOperator(v, diagm(-1 => sqrt.(1:(modecutoff(v)))))

"""
destroy(v::VibrationalMode)
@@ -34,15 +34,15 @@ destroy(v::VibrationalMode) = create(v)'
number(v::VibrationalMode)
Returns the number operator for `v` such that: `number(v) * v[i] = i * v[i]`.
"""
number(v::VibrationalMode) = SparseOperator(v, diagm(0 => 0:(v.N)))
number(v::VibrationalMode) = SparseOperator(v, diagm(0 => 0:(modecutoff(v))))

"""
displace(v::VibrationalMode, α::Number; method="truncated")
Returns the displacement operator ``D(α)`` corresponding to `v`.

If `method="truncated"` (default), the matrix elements are computed according to
``D(α) = exp[αa^† - α^*a]`` where ``a`` and ``a^†`` live in a truncated Hilbert space of
dimension `v.N+1`.
dimension `modecutoff(v)+1`.
Otherwise if `method="analytic"`, the matrix elements are computed assuming an
infinite-dimension Hilbert space. In general, this option will not return a unitary operator.
"""
@@ -51,13 +51,13 @@ function displace(v::VibrationalMode, α::Number; method = "truncated")
# Above line commented out to allow for Hamiltonian construction even if vibrational mode N = 0.
# May want to think of a different way to perform this check in the future.
@assert method in ["truncated", "analytic"] "method ∉ [truncated, analytic]"
D = zeros(ComplexF64, v.N + 1, v.N + 1)
D = zeros(ComplexF64, modecutoff(v) + 1, modecutoff(v) + 1)
if α == 0
return one(v)
elseif method ≡ "analytic"
@inbounds begin
@simd for n in 1:(v.N + 1)
@simd for m in 1:(v.N + 1)
@simd for n in 1:(modecutoff(v) + 1)
@simd for m in 1:(modecutoff(v) + 1)
D[n, m] = _Dnm(α, n, m)
end
end
@@ -80,15 +80,15 @@ assuming an infinite-dimensional Hilbert space, is used:
``[ρ_{th}]_{ij} = δ_{ij} \\frac{nⁱ}{(n+1)^{i+1}}.``
"""
function thermalstate(v::VibrationalMode, n̄::Real; method = "truncated")
@assert v.N ≥ n̄ "`n̄` must be less than `v.N`"
@assert modecutoff(v) ≥ n̄ "`n̄` must be less than `modecutoff(v)`"
@assert method in ["truncated", "analytic"] "method ∉ [truncated, analytic]"
if n̄ == 0
return v[0] ⊗ v[0]'
elseif method ≡ "truncated"
d = [(n̄ / (n̄ + 1))^i for i in 0:(v.N)]
d = [(n̄ / (n̄ + 1))^i for i in 0:(modecutoff(v))]
return DenseOperator(v, diagm(0 => d) ./ sum(d))
elseif method ≡ "analytic"
return DenseOperator(v, diagm(0 => [(n̄ / (n̄ + 1))^i / (n̄ + 1) for i in 0:(v.N)]))
return DenseOperator(v, diagm(0 => [(n̄ / (n̄ + 1))^i / (n̄ + 1) for i in 0:(modecutoff(v))]))
end
end

@@ -99,10 +99,10 @@ Returns a coherent state on `v` with complex amplitude ``α``.
function coherentstate(v::VibrationalMode, α::Number)
# this implementation is the same as in QuantumOptics.jl, but there the function is
# restricted to v::FockBasis, so we must reimplement here
@assert v.N ≥ abs(α) "`α` must be less than `v.N`"
k = zeros(ComplexF64, v.N + 1)
@assert modecutoff(v) ≥ abs(α) "`α` must be less than `modecutoff(v)`"
k = zeros(ComplexF64, modecutoff(v) + 1)
k[1] = exp(-abs2(α) / 2)
@inbounds for n in 1:(v.N)
@inbounds for n in 1:(modecutoff(v))
k[n + 1] = k[n] * α / √n
end
return Ket(v, k)
@@ -118,7 +118,7 @@ complex amplitude of the displacement.
displacement operator is computed (see: [`displace`](@ref)) .
"""
function coherentthermalstate(v::VibrationalMode, n̄::Real, α::Number; method = "truncated")
@assert (v.N ≥ n̄ && v.N ≥ abs(α)) "`n̄`, `α` must be less than `v.N`"
@assert (modecutoff(v) ≥ n̄ && modecutoff(v) ≥ abs(α)) "`n̄`, `α` must be less than `modecutoff(v)`"
@assert method in ["truncated", "analytic"] "method ∉ [truncated, analytic]"
if method ≡ "truncated"
d = displace(v, α)
@@ -159,10 +159,10 @@ end
ionstate(I::Ion, sublevelalias::String) = ionstate(I, sublevel(I, sublevelalias))
ionstate(I::Ion, sublevel::Int) = basisstate(I, sublevel)
function ionstate(IC::IonTrap, states::Union{Tuple{String, Real}, String, Int}...)
ions = IC.ions
L = length(ions)
allions = ions(IC)
L = length(allions)
@assert L ≡ length(states) "wrong number of states"
return tensor([ionstate(ions[i], states[i]) for i in 1:L]...)
return tensor([ionstate(allions[i], states[i]) for i in 1:L]...)
end
ionstate(T::Chamber, states::Union{Tuple{String, Real}, String, Int}...) =
ionstate(T.iontrap, states...)
@@ -194,11 +194,11 @@ function ionprojector(
sublevels::Union{Tuple{String, Real}, String, Int}...;
only_ions = false
)
ions = IC.ions
L = length(ions)
allions = ions(IC)
L = length(allions)
@assert L ≡ length(sublevels) "wrong number of sublevels"
allmodes = modes(IC)
observable = tensor([projector(ions[i][sublevels[i]]) for i in 1:L]...)
observable = tensor([projector(allions[i][sublevels[i]]) for i in 1:L]...)
if !only_ions
for mode in allmodes
observable = observable ⊗ one(mode)
@@ -211,7 +211,7 @@ function ionprojector(
sublevels::Union{Tuple{String, Real}, String, Int}...;
only_ions = false
)
return ionprojector(T.iontrap, sublevels..., only_ions = only_ions)
return ionprojector(iontrap(T), sublevels..., only_ions = only_ions)
end

#############################################################################################
2 changes: 1 addition & 1 deletion src/vibrational_modes.jl
Original file line number Diff line number Diff line change
@@ -118,6 +118,6 @@ Base.show(io::IO, V::VibrationalMode) = print(
)

function Base.getindex(V::VibrationalMode, n::Int)
@assert 0 <= n <= V.N "n ∉ [0, $(V.N+1)]"
@assert 0 <= n <= modecutoff(V) "n ∉ [0, $(V.N+1)]"
return basisstate(V, n + 1)
end