Setting up a calculation
Contents
Index
BialkaliSpectrum.K40Rb87.KRb_Nuclear_Neyenhuis
BialkaliSpectrum.K40Rb87.KRb_Nuclear_Ospelkaus
BialkaliSpectrum.K40Rb87.KRb_Parameters_Neyenhuis
BialkaliSpectrum.K40Rb87.KRb_Parameters_Ospelkaus
BialkaliSpectrum.K40Rb87.KRb_Polarizability
BialkaliSpectrum.K40Rb87.KRb_Zeeman
BialkaliSpectrum.Toy.TOY_PARAMETERS
BialkaliSpectrum.ExternalFields
BialkaliSpectrum.HamiltonianParts
BialkaliSpectrum.MolecularParameters
BialkaliSpectrum.NuclearParameters
BialkaliSpectrum.Polarizability
BialkaliSpectrum.SphericalVector
BialkaliSpectrum.State
BialkaliSpectrum.State
BialkaliSpectrum.ZeemanParameters
Base.convert
Base.show
BialkaliSpectrum.K40Rb87.KRbState
BialkaliSpectrum.K40Rb87.make_krb_hamiltonian_parts
BialkaliSpectrum.Toy.make_toy_hamiltonian_parts
BialkaliSpectrum.VectorX
BialkaliSpectrum.VectorY
BialkaliSpectrum.VectorZ
BialkaliSpectrum.basis_index
BialkaliSpectrum.basis_state
BialkaliSpectrum.generate_fields_scan
BialkaliSpectrum.hamiltonian
BialkaliSpectrum.make_hamiltonian_parts
Molecular parameters
Structs for organizing the coupling constants, dipole moment, and nuclear spins of the molecule.
Types
BialkaliSpectrum.MolecularParameters
— TypeMolecularParameters
Contains the coupling constants for the molecular Hamiltonian and the nuclear angular momenta (needed to construct the basis states).
BialkaliSpectrum.NuclearParameters
— TypeNuclearParameters
Contains the nuclear electric quadrupole moments, nuclear spin-rotation couplings, and nuclear spin-spin coupling for calculating the hyperfine Hamiltonian.
BialkaliSpectrum.Polarizability
— TypePolarizability
Contains the parallel and perpendicular ac polarizabilities at a particular optical wavelength λ.
So far, it is assumed that λ is far-detuned from any electronic transitions, such that the polarizability does not depend on rotational state $N$.
BialkaliSpectrum.ZeemanParameters
— TypeZeemanParameters
Contains the g-factors and nuclear shielding factors for computing Zeeman shifts.
Constants
Hamiltonian
Functions for building the molecular Hamiltonian.
The HamiltonianParts
struct is used as an input to several analysis methods, e.g. computing dipole matrix elements, since it holds the individual terms in the Hamiltonian separately.
Types
BialkaliSpectrum.HamiltonianParts
— TypeHamiltonianParts
Contains all parts of the Hamiltonian except external fields.
Zeeman, dc Stark, and ac Stark terms are stored as vectors of matrices that can be contracted with the appropriate external field tensors. This allows the full Hamiltonian to be re-constructed at various field values and orientations without recalculating all of the matrix elements.
Should be created by make_hamiltonian_parts
or make_krb_hamiltonian_parts
. Used as an input to get_spectrum
and hamiltonian
.
See also make_hamiltonian_parts
, make_krb_hamiltonian_parts
.
Methods
BialkaliSpectrum.hamiltonian
— Methodhamiltonian(parts, external_fields)
Construct the full Hamiltonian including magnetic, electric, and optical fields.
The field-independent building blocks in parts
can be reused over calls to hamiltonian
to avoid recalculating the matrix elements each time.
See also make_hamiltonian_parts
, make_krb_hamiltonian_parts
, ExternalFields
.
BialkaliSpectrum.make_hamiltonian_parts
— Methodmake_hamiltonian_parts(molecular_parameters, N_max)
Construct all parts of the Hamiltonian that do not depend on external fields.
The size of the basis is determined by molecular_parameters
, which contains the nuclear spin quantum numbers molecular_parameters.I
, and the rotational states 0:N_max
to include.
See also make_krb_hamiltonian_parts
.
Fields
Utilities for defining the external magnetic, electric, and optical fields experienced by the molecules.
The ExternalFields
holds the external fields, which are each defined by a SphericalVector
with a magnitude
, polar angle $\theta$, and azimuthal angle $\phi$. Note that the molecular basis is defined with respect to the $z$-axis ($\theta = 0$).
A Vector{ExternalFields}
is expected by get_spectra
, and can be generated manually or with generate_fields_scan
.
Types
BialkaliSpectrum.ExternalFields
— TypeExternalFields(B::SphericalVector, E::SphericalVector, Optical::Vector{SphericalVector})
ExternalFields(B::Float64, E::Float64)
External magnetic, electric, optical fields to use in constructing the Hamiltonian.
If B
and E
are provided as Float64
s, then the fields are assumed to be along z
. The Optical
argument can also be left as an empty vector []
.
See also get_spectrum
, hamiltonian
, SphericalVector
.
TODO: add the other signatures
Examples
julia> ExternalFields(VectorZ(545.9), VectorX(1020.0), [])
ExternalFields(SphericalVector(545.9, 0.0, 0.0), SphericalVector(1020.0, 1.5707963267948966, 0.0), SphericalVector[])
julia> ExternalFields(545.9, 1020.0)
ExternalFields(SphericalVector(545.9, 0.0, 0.0), SphericalVector(1020.0, 0.0, 0.0), SphericalVector[])
BialkaliSpectrum.SphericalVector
— TypeSphericalVector(magnitude, θ, φ)
Construct a vector with magnitude
, polar angle θ
, and azimuthal angle φ
.
Represents an external field vector in spherical coordinates used to construct ExternalFields
. Currently, SphericalVector
s may be negated but no other mathematical operations are implemented.
Vectors along $x$, $y$, or $z$ can be quickly constructed using VectorX
, VectorY
, and VectorZ
, respectively.
See also SphericalUnitVector
, ExternalFields
.
Methods
BialkaliSpectrum.VectorX
— MethodVectorX(magnitude)
Construct a SphericalVector
with magnitude
along x
.
Examples
julia> VectorX(5.0)
SphericalVector(5.0, 1.5707963267948966, 0.0)
BialkaliSpectrum.VectorY
— MethodVectorY(magnitude)
Construct a SphericalVector
with magnitude
along y
.
Examples
julia> VectorY(2.0)
SphericalVector(2.0, 1.5707963267948966, 1.5707963267948966)
BialkaliSpectrum.VectorZ
— MethodVectorZ(magnitude)
Construct a SphericalVector
with magnitude
along z
.
Examples
julia> VectorZ(10.0)
SphericalVector(10.0, 0.0, 0.0)
BialkaliSpectrum.generate_fields_scan
— Methodgenerate_fields_scan(Bs, Es, Opticals)
Produce a vector of ExternalFields
for creating a scan of spectra as a function of fields.
TODO: explain the behavior as a product of the lists (not zipped)
Basis states
Utilities for working with molecular states in the uncoupled basis $|N, m_N, m_{i1}, m_{i2} \rangle$.
The basis states are indexed with the quantum numbers changing fastest on the right. For example, with $I_1 = 1/2$ and $I_2 = 1/2$, the first few states are:
- $|0, 0, -1/2, -1/2⟩$
- $|0, 0, -1/2, +1/2⟩$
- $|0, 0, +1/2, -1/2⟩$
- $|0, 0, +1/2, +1/2⟩$
- $|1, -1, -1/2, -1/2⟩$
etc.
The methods basis_state
and basis_index
are used to switch between an index
and the corresponding basis State
.
Types
BialkaliSpectrum.State
— TypeState(N, mₙ, I, mᵢ)
State(N, mₙ, I₁, mᵢ₁, I₂, mᵢ₂)
Represents a molecular state in the uncoupled basis $|N, mₙ, I₁, mᵢ₁, I₂, mᵢ₂⟩$.
Write something more descriptive here
Standard library interfaces
Base.convert
— FunctionBase.convert(t::Type{NamedTuple}, s::State)
Returns a named tuple with fields N
, m_n
, I_1
, m_i1
, I_2
, m_i2
from s
.
This is a utility function to simplify outputting to a DataFrame
.
Base.show
— FunctionBase.show(io::IO, s::State)
Pretty prints a State
in ket notation, $|N, m_N, m_{i1}, m_{i2}⟩$.
Methods
BialkaliSpectrum.basis_index
— Methodbasis_index(s::State)
Returns index of state s
in the basis.
The uncoupled basis $|N, mₙ, I₁, mᵢ₁, I₂, mᵢ₂⟩$ is ordered with the quantum numbers on the left changing the slowest.
See also State
, basis_index
.
Examples
julia> import BialkaliSpectrum.K40Rb87: KRbState
julia> basis_index(KRbState(1, 1, -4, 1/2))
111
julia> basis_index(basis_state(42, 4, 3/2))
42
BialkaliSpectrum.basis_state
— Methodbasis_state(i, I₁, I₂)
Returns the State
corresponding to the i
th member of the basis.
The uncoupled basis $|N, mₙ, I₁, mᵢ₁, I₂, mᵢ₂⟩$ is ordered with the quantum numbers on the left changing the slowest.
See also State
, basis_index
.
Examples
julia> s = basis_state(1, 4, 3/2)
BialkaliSpectrum.State basis state:
|0, 0, -4, -3/2⟩
julia> s = basis_state(37, 4, 3/2)
BialkaliSpectrum.State basis state:
|1, -1, -4, -3/2⟩
Molecule-specific definitions
Constants and method defs for working with specific molecule species.
${}^{40}\text{K}^{87}\text{Rb}$ (K40Rb87
)
BialkaliSpectrum.K40Rb87.KRb_Nuclear_Neyenhuis
— ConstantKRb_Nuclear_Neyenhuis
Experimental values from Neyenhuis et al., PRL 109, 230403 (2012)
BialkaliSpectrum.K40Rb87.KRb_Nuclear_Ospelkaus
— ConstantKRb_Nuclear_Ospelkaus
Experimental values from Ospelkaus et al., PRL 104, 030402 (2010)
BialkaliSpectrum.K40Rb87.KRb_Parameters_Neyenhuis
— ConstantKRb_Parameters_Neyenhuis
Experimental values from Neyenhuis et al., PRL 109, 230403 (2012)
BialkaliSpectrum.K40Rb87.KRb_Parameters_Ospelkaus
— ConstantKRb_Parameters_Ospelkaus
Experimental values from Ospelkaus et al., PRL 104, 030402 (2010)
BialkaliSpectrum.K40Rb87.KRb_Polarizability
— ConstantKRb_Polarizability
Experimental values from Neyenhuis et al., PRL 109, 230403 (2012)
BialkaliSpectrum.K40Rb87.KRb_Zeeman
— ConstantKRb_Zeeman
Theoretical values from Aldegunde et al. PRA (2008)
BialkaliSpectrum.K40Rb87.KRbState
— MethodKRbState(N, mₙ, mK, mRb)
Creates a basis state $|N, m_n, m_{\text{K}}, m_{\text{Rb}}⟩$ for ${}^{40}\text{K}^{87}\text{Rb}$.
This is a wrapper around State
to avoid having to specify the nuclear spins $I_k$ each time.
See also State
.
BialkaliSpectrum.K40Rb87.make_krb_hamiltonian_parts
— Methodmake_krb_hamiltonian_parts(N_max::Int)
Construct all parts of the ${}^{40}\text{K}^{87}\text{Rb}$ Hamiltonian that do not depend on external fields.
The rotational states 0:N_max
are included. This is a shortcut method that replaces make_hamiltonian_parts
for KRb.
See also make_hamiltonian_parts
.
Toy molecule without hyperfine structure (Toy
)
BialkaliSpectrum.Toy.TOY_PARAMETERS
— ConstantTOY_PARAMETERS
Toy model values with dipole = 1 D and no hyperfine structure. Intended for testing.
BialkaliSpectrum.State
— MethodState(N, mₙ)
Creates a basis state $|N, m_n⟩$ for the toy molecule.
This is a wrapper around State
to avoid having to specify the nuclear spins $I_k$ each time.
See also State
.
BialkaliSpectrum.Toy.make_toy_hamiltonian_parts
— Methodmake_toy_hamiltonian_parts(N_max::Int)
Construct all parts of the toy molecule Hamiltonian that do not depend on external fields.
The rotational states 0:N_max
are included. This is a shortcut method that replaces make_hamiltonian_parts
for the toy molecule.
See also make_hamiltonian_parts
.