Divalent atom functions

Note

This is a completely new module added in ARC 3.0.0 version. See more at E. J. Robertson, N. Šibalić, R. M. Potvliege and M. P. A. Jones, arXiv:2007.12016 .

Overview

DivalentAtom Methods

DivalentAtom.getDipoleMatrixElement(n1, l1, …) Dipole matrix element \(\langle n_1 l_1 j_1 m_{j_1} |e\mathbf{r}|\ n_2 l_2 j_2 m_{j_2}\rangle\) in units of \(a_0 e\)
DivalentAtom.getTransitionWavelength(n1, l1, …) Calculated transition wavelength (in vacuum) in m.
DivalentAtom.getTransitionFrequency(n1, l1, …) Calculated transition frequency in Hz
DivalentAtom.getRabiFrequency(n1, l1, j1, …) Returns a Rabi frequency for resonantly driven atom in a center of TEM00 mode of a driving field
DivalentAtom.getRabiFrequency2(n1, l1, j1, …) Returns a Rabi frequency for resonant excitation with a given electric field amplitude
DivalentAtom.getStateLifetime(n, l, j[, …]) Returns the lifetime of the state (in s)
DivalentAtom.getTransitionRate(n1, l1, j1, …) Transition rate due to coupling to vacuum modes (black body included)
DivalentAtom.getReducedMatrixElementJ_asymmetric(n1, …) Reduced matrix element in \(J\) basis, defined in asymmetric notation.
DivalentAtom.getReducedMatrixElementJ(n1, …) Reduced matrix element in \(J\) basis (symmetric notation)
DivalentAtom.getReducedMatrixElementL(n1, …) Reduced matrix element in \(L\) basis (symmetric notation)
DivalentAtom.getRadialMatrixElement(n1, l1, …) Radial part of the dipole matrix element
DivalentAtom.getQuadrupoleMatrixElement(n1, …) Radial part of the quadrupole matrix element
DivalentAtom.getPressure(temperature) Vapour pressure (in Pa) at given temperature
DivalentAtom.getNumberDensity(temperature) Atom number density at given temperature
DivalentAtom.getAverageInteratomicSpacing(…) Returns average interatomic spacing in atomic vapour
DivalentAtom.getEnergy(n, l, j[, s]) Energy of the level relative to the ionisation level (in eV)
DivalentAtom.getZeemanEnergyShift(l, j, mj, …) Retuns linear (paramagnetic) Zeeman shift.
DivalentAtom.getQuantumDefect(n, l, j[, s]) Quantum defect of the level.
DivalentAtom.getC6term(n, l, j, n1, l1, j1, …) C6 interaction term for the given two pair-states
DivalentAtom.getC3term(n, l, j, n1, l1, j1, …) C3 interaction term for the given two pair-states
DivalentAtom.getEnergyDefect(n, l, j, n1, …) Energy defect for the given two pair-states (one of the state has two atoms in the same state)
DivalentAtom.getEnergyDefect2(n, l, j, nn, …) Energy defect for the given two pair-states
DivalentAtom.updateDipoleMatrixElementsFile() Updates the file with pre-calculated dipole matrix elements.
DivalentAtom.getRadialCoupling(n, l, j, n1, …) Returns radial part of the coupling between two states (dipole and quadrupole interactions only)
DivalentAtom.getAverageSpeed(temperature) Average (mean) speed at a given temperature
DivalentAtom.getLiteratureDME(n1, l1, j1, …) Returns literature information on requested transition.

Detailed documentation

class arc.divalent_atom_functions.DivalentAtom(preferQuantumDefects=True, cpp_numerov=True)[source]

Bases: arc.alkali_atom_functions.AlkaliAtom

Implements general calculations for Alkaline Earths, and other divalent atoms.

This class inherits arc.alkali_atom_functions.AlkaliAtom . Most of the methods can be directly used from there, and the source for them is provided in the base class. Few methods that are implemented differently for Alkaline Earths are defined here.

Parameters:
  • preferQuantumDefects (bool) – Use quantum defects for energy level calculations. If False, uses NIST ASD values where available. If True, uses quantum defects for energy calculations for principal quantum numbers within the range specified in defectFittingRange which is specified for each element and series separately. For principal quantum numbers below this value, NIST ASD values are used if existing, since quantum defects. Default is True.
  • cpp_numerov (bool) – This switch for Alkaline Earths at the moment doesn’t have any effect since wavefunction calculation function is not implemented (d.m.e. and quadrupole matrix elements are calculated directly semiclassically)
corePotential(l, r)[source]

Not implemented for Alkaline earths

defectFittingRange = {}

Used for AlkalineEarths to define minimum and maximum principal quantum number for which quantum defects are valid. Ranges are stored under keys defined as state terms ({‘stateLabel’:[minN, maxN]}, e.g. ‘1S0’). Dictionary returns array stating minimal and maximal principal quantum number for which quantum defects were fitted. For example:

limits = self.defectFittingRange['1S0']
print("Minimal n = %d" % limits[0])
print("Maximal n = %d" % limits[1]) 1
effectiveCharge(l, r)[source]

Not implemented for Alkaline earths

energyLevelsExtrapolated = False

flag that is turned to True if the energy levels of this atom were calculated by extrapolating with quantum defects values outside the quantum defect fitting range.

getAverageInteratomicSpacing(temperature)

Returns average interatomic spacing in atomic vapour

See calculation of basic properties example snippet.

Parameters:temperature (float) – temperature of the atomic vapour
Returns:average interatomic spacing in m
Return type:float
getAverageSpeed(temperature)

Average (mean) speed at a given temperature

Parameters:temperature (float) – temperature (K)
Returns:mean speed (m/s)
Return type:float
getC3term(n, l, j, n1, l1, j1, n2, l2, j2, s=0.5)

C3 interaction term for the given two pair-states

Calculates \(C_3\) intaraction term for
\(|n,l,j,n,l,j\rangle \leftrightarrow |n_1,l_1,j_1,n_2,l_2,j_2\rangle\)
Parameters:
  • n (int) – principal quantum number
  • l (int) – orbital angular momenutum
  • j (float) – total angular momentum
  • n1 (int) – principal quantum number
  • l1 (int) – orbital angular momentum
  • j1 (float) – total angular momentum
  • n2 (int) – principal quantum number
  • l2 (int) – orbital angular momentum
  • j2 (float) – total angular momentum
  • s (float) – optional, total spin angular momentum of state. By default 0.5 for Alkali atoms.
Returns:

\(C_3 = \frac{\langle n,l,j |er |n_1,l_1,j_1\rangle \langle n,l,j |er|n_2,l_2,j_2\rangle}{4\pi\varepsilon_0}\) (\(h\) Hz m \({}^3\)).

Return type:

float

getC6term(n, l, j, n1, l1, j1, n2, l2, j2, s=0.5)

C6 interaction term for the given two pair-states

Calculates \(C_6\) intaraction term for \(|n,l,j,n,l,j \rangle \leftrightarrow |n_1,l_1,j_1,n_2,l_2,j_2\rangle\). For details of calculation see Ref. [2].

Parameters:
  • n (int) – principal quantum number
  • l (int) – orbital angular momenutum
  • j (float) – total angular momentum
  • n1 (int) – principal quantum number
  • l1 (int) – orbital angular momentum
  • j1 (float) – total angular momentum
  • n2 (int) – principal quantum number
  • l2 (int) – orbital angular momentum
  • j2 (float) – total angular momentum
  • s (float) – optional, total spin angular momentum of state. By default 0.5 for Alkali atoms.
Returns:

\(C_6 = \frac{1}{4\pi\varepsilon_0} \frac{|\langle n,l,j |er|n_1,l_1,j_1\rangle|^2| \langle n,l,j |er|n_2,l_2,j_2\rangle|^2} {E(n_1,l_1,j_2,n_2,j_2,j_2)-E(n,l,j,n,l,j)}\) (\(h\) Hz m \({}^6\)).

Return type:

float

Example

We can reproduce values from Ref. [2] for C3 coupling to particular channels. Taking for example channels described by the Eq. (50a-c) we can get the values:

from arc import *

channels = [[70,0,0.5, 70, 1,1.5, 69,1, 1.5],\
            [70,0,0.5, 70, 1,1.5, 69,1, 0.5],\
            [70,0,0.5, 69, 1,1.5, 70,1, 0.5],\
            [70,0,0.5, 70, 1,0.5, 69,1, 0.5]]

print(" = = = Caesium = = = ")
atom = Caesium()
for channel in channels:
    print("%.0f  GHz (mu m)^6" % ( atom.getC6term(*channel)
                                  / C_h * 1.e27 ))

print("\n = = = Rubidium  = = =")
atom = Rubidium()
for channel in channels:
    print("%.0f  GHz (mu m)^6" % ( atom.getC6term(*channel)
                                  / C_h * 1.e27 ))

Returns:

 = = = Caesium = = =
722  GHz (mu m)^6
316  GHz (mu m)^6
383  GHz (mu m)^6
228  GHz (mu m)^6

 = = = Rubidium  = = =
799  GHz (mu m)^6
543  GHz (mu m)^6
589  GHz (mu m)^6
437  GHz (mu m)^6

which is in good agreement with the values cited in the Ref. [2]. Small discrepancies for Caesium originate from slightly different quantum defects used in calculations.

References

[2](1, 2, 3) T. G. Walker, M. Saffman, PRA 77, 032723 (2008) https://doi.org/10.1103/PhysRevA.77.032723
getDipoleMatrixElement(n1, l1, j1, mj1, n2, l2, j2, mj2, q, s=0.5)

Dipole matrix element \(\langle n_1 l_1 j_1 m_{j_1} |e\mathbf{r}|\ n_2 l_2 j_2 m_{j_2}\rangle\) in units of \(a_0 e\)

Parameters:
  • l1, j1, mj1 (n1.) – principal, orbital, total angular momentum, and projection of total angular momenutum for state 1
  • l2, j2, mj2 (n2.) – principal, orbital, total angular momentum, and projection of total angular momenutum for state 2
  • q (int) – specifies transition that the driving field couples to, +1, 0 or -1 corresponding to driving \(\sigma^+\), \(\pi\) and \(\sigma^-\) transitions respectively.
  • s (float) – optional, total spin angular momentum of state. By default 0.5 for Alkali atoms.
Returns:

dipole matrix element( \(a_0 e\))

Return type:

float

Example

For example, calculation of \(5 S_{1/2}m_j=-\frac{1}{2}\ \rightarrow 5 P_{3/2}m_j=-\frac{3}{2}\) transition dipole matrix element for laser driving \(\sigma^-\) transition:

from arc import *
atom = Rubidium()
# transition 5 S_{1/2} m_j=-0.5 -> 5 P_{3/2} m_j=-1.5
# for laser driving sigma- transition
print(atom.getDipoleMatrixElement(5,0,0.5,-0.5,5,1,1.5,-1.5,-1))
getDipoleMatrixElementHFS(n1, l1, j1, f1, mf1, n2, l2, j2, f2, mf2, q, s=0.5)

Dipole matrix element for hyperfine structure resolved transitions \(\langle n_1 l_1 j_1 f_1 m_{f_1} |e\mathbf{r}|\ n_2 l_2 j_2 f_2 m_{f_2}\rangle\) in units of \(a_0 e\)

For hyperfine resolved transitions, the dipole matrix element is \(\langle n_1,\ell_1,j_1,f_1,m_{f1} | \ \mathbf{\hat{r}}\cdot \mathbf{\varepsilon}_q \ | n_2,\ell_2,j_2,f_2,m_{f2} \rangle = (-1)^{f_1-m_{f1}} \ \left( \ \begin{matrix} \ f_1 & 1 & f_2 \\ \ -m_{f1} & q & m_{f2} \ \end{matrix}\right) \ \langle n_1 \ell_1 j_1 f_1|| r || n_2 \ell_2 j_2 f_2 \rangle,\) where \(\langle n_1 \ell_1 j_1 f_1 ||r|| n_2 \ell_2 j_2 f_2 \rangle \ = (-1)^{j_1+I+F_2+1}\sqrt{(2f_1+1)(2f_2+1)} ~ \ \left\{ \begin{matrix}\ F_1 & 1 & F_2 \\ \ j_2 & I & j_1 \ \end{matrix}\right\}~ \ \langle n_1 \ell_1 j_1||r || n_2 \ell_2 j_2 \rangle.\)

Parameters:
  • l1, j1, f1, mf1 (n1.) – principal, orbital, total orbital, fine basis (total atomic) angular momentum, and projection of total angular momenutum for state 1
  • l2, j2, f2, mf2 (n2.) – principal, orbital, total orbital, fine basis (total atomic) angular momentum, and projection of total angular momenutum for state 2
  • q (int) – specifies transition that the driving field couples to, +1, 0 or -1 corresponding to driving \(\sigma^+\), \(\pi\) and \(\sigma^-\) transitions respectively.
  • s (float) – optional, total spin angular momentum of state. By default 0.5 for Alkali atoms.
Returns:

dipole matrix element( \(a_0 e\))

Return type:

float

getEnergy(n, l, j, s=None)[source]

Energy of the level relative to the ionisation level (in eV)

Returned energies are with respect to the center of gravity of the hyperfine-split states. If preferQuantumDefects =False (set during initialization) program will try use NIST energy value, if such exists, falling back to energy calculation with quantum defects if the measured value doesn’t exist. For preferQuantumDefects =True, program will always calculate energies from quantum defects (useful for comparing quantum defect calculations with measured energy level values).

Parameters:
  • n (int) – principal quantum number
  • l (int) – orbital angular momentum
  • j (float) – total angular momentum
  • s (float) – optional, total spin angular momentum. Default value of 0.5 is correct for Alkali atoms, and has to be specified explicitly for divalent atoms.
Returns:

state energy (eV)

Return type:

float

getEnergyDefect(n, l, j, n1, l1, j1, n2, l2, j2, s=0.5)

Energy defect for the given two pair-states (one of the state has two atoms in the same state)

Energy difference between the states \(E(n_1,l_1,j_1,n_2,l_2,j_2) - E(n,l,j,n,l,j)\)

Parameters:
  • n (int) – principal quantum number
  • l (int) – orbital angular momenutum
  • j (float) – total angular momentum
  • n1 (int) – principal quantum number
  • l1 (int) – orbital angular momentum
  • j1 (float) – total angular momentum
  • n2 (int) – principal quantum number
  • l2 (int) – orbital angular momentum
  • j2 (float) – total angular momentum
  • s (float) – optional. Spin angular momentum (default 0.5 for Alkali)
Returns:

energy defect (SI units: J)

Return type:

float

getEnergyDefect2(n, l, j, nn, ll, jj, n1, l1, j1, n2, l2, j2, s=0.5)

Energy defect for the given two pair-states

Energy difference between the states \(E(n_1,l_1,j_1,n_2,l_2,j_2) - E(n,l,j,nn,ll,jj)\)

See pair-state energy defects example snippet.

Parameters:
  • n (int) – principal quantum number
  • l (int) – orbital angular momenutum
  • j (float) – total angular momentum
  • nn (int) – principal quantum number
  • ll (int) – orbital angular momenutum
  • jj (float) – total angular momentum
  • n1 (int) – principal quantum number
  • l1 (int) – orbital angular momentum
  • j1 (float) – total angular momentum
  • n2 (int) – principal quantum number
  • l2 (int) – orbital angular momentum
  • j2 (float) – total angular momentum
  • s (float) – optional. Spin angular momentum (default 0.5 for Alkali)
Returns:

energy defect (SI units: J)

Return type:

float

getLiteratureDME(n1, l1, j1, n2, l2, j2, s=0)[source]

Returns literature information on requested transition.

Parameters:
  • n1,l1,j1 – one of the states we are coupling
  • n2,l2,j2 – the other state to which we are coupling
  • s – (optional) spin of the state. Default s=0.
Returns:

hasLiteratureValue?, dme, referenceInformation

If Boolean value is True, a literature value for dipole matrix element was found and reduced DME in J basis is returned as the number. The third returned argument (array) contains additional information about the literature value in the following order [ typeOfSource, errorEstimate , comment , reference, reference DOI] upon success to find a literature value for dipole matrix element:

  • typeOfSource=1 if the value is theoretical calculation; otherwise, if it is experimentally obtained value typeOfSource=0
  • comment details where within the publication the value can be found
  • errorEstimate is absolute error estimate
  • reference is human-readable formatted reference
  • reference DOI provides link to the publication.

Boolean value is False, followed by zero and an empty array if no literature value for dipole matrix element is found.

Return type:

bool, float, [int,float,string,string,string]

Note

The literature values are stored in /data folder in <element name>_literature_dme.csv files as a ; separated values. Each row in the file consists of one literature entry, that has information in the following order:

  • n1
  • l1
  • j1
  • n2
  • l2
  • j2
  • s
  • dipole matrix element reduced l basis (a.u.)
  • comment (e.g. where in the paper value appears?)
  • value origin: 1 for theoretical; 0 for experimental values
  • accuracy
  • source (human readable formatted citation)
  • doi number (e.g. 10.1103/RevModPhys.82.2313 )

If there are several values for a given transition, program outputs the value that has smallest error (under column accuracy). The list of values can be expanded - every time program runs this file is read and the list is parsed again for use in calculations.

getNumberDensity(temperature)

Atom number density at given temperature

See calculation of basic properties example snippet.

Parameters:temperature (float) – temperature in K
Returns:atom concentration in \(1/m^3\)
Return type:float
getPressure(temperature)

Vapour pressure (in Pa) at given temperature

Parameters:temperature (float) – temperature in K
Returns:vapour pressure in Pa
Return type:float
getQuadrupoleMatrixElement(n1, l1, j1, n2, l2, j2, s=0.5)[source]

Radial part of the quadrupole matrix element

Calculates \(\int \mathbf{d}r~R_{n_1,l_1,j_1}(r)\cdot R_{n_1,l_1,j_1}(r) \cdot r^4\). See Quadrupole calculation example snippet .

Parameters:
  • n1 (int) – principal quantum number of state 1
  • l1 (int) – orbital angular momentum of state 1
  • j1 (float) – total angular momentum of state 1
  • n2 (int) – principal quantum number of state 2
  • l2 (int) – orbital angular momentum of state 2
  • j2 (float) – total angular momentum of state 2
  • s (float) – optional. Spin of the state. Default 0.5 is for Alkali
Returns:

quadrupole matrix element (\(a_0^2 e\)).

Return type:

float

getQuantumDefect(n, l, j, s=0.5)

Quantum defect of the level.

For an example, see Rydberg energy levels example snippet.

Parameters:
  • n (int) – principal quantum number
  • l (int) – orbital angular momentum
  • j (float) – total angular momentum
  • s (float) – (optional). Total spin angular momentum. Default value of 0.5 correct for Alkali atoms. For divalent atoms it has to be explicitly defined.
Returns:

quantum defect

Return type:

float

getRabiFrequency(n1, l1, j1, mj1, n2, l2, j2, q, laserPower, laserWaist, s=0.5)

Returns a Rabi frequency for resonantly driven atom in a center of TEM00 mode of a driving field

Parameters:
  • n1,l1,j1,mj1 – state from which we are driving transition
  • n2,l2,j2 – state to which we are driving transition
  • q – laser polarization (-1,0,1 correspond to \(\sigma^-\), \(\pi\) and \(\sigma^+\) respectively)
  • laserPower – laser power in units of W
  • laserWaist – laser \(1/e^2\) waist (radius) in units of m
  • s (float) – optional, total spin angular momentum of state. By default 0.5 for Alkali atoms.
Returns:

Frequency in rad \(^{-1}\). If you want frequency in Hz, divide by returned value by \(2\pi\)

Return type:

float

getRabiFrequency2(n1, l1, j1, mj1, n2, l2, j2, q, electricFieldAmplitude, s=0.5)

Returns a Rabi frequency for resonant excitation with a given electric field amplitude

Parameters:
  • n1,l1,j1,mj1 – state from which we are driving transition
  • n2,l2,j2 – state to which we are driving transition
  • q – laser polarization (-1,0,1 correspond to \(\sigma^-\), \(\pi\) and \(\sigma^+\) respectively)
  • electricFieldAmplitude – amplitude of electric field driving (V/m)
  • s (float) – optional, total spin angular momentum of state. By default 0.5 for Alkali atoms.
Returns:

Frequency in rad \(^{-1}\). If you want frequency in Hz, divide by returned value by \(2\pi\)

Return type:

float

getRadialCoupling(n, l, j, n1, l1, j1, s=0.5)

Returns radial part of the coupling between two states (dipole and quadrupole interactions only)

Parameters:
  • n1 (int) – principal quantum number
  • l1 (int) – orbital angular momentum
  • j1 (float) – total angular momentum
  • n2 (int) – principal quantum number
  • l2 (int) – orbital angular momentum
  • j2 (float) – total angular momentum
  • s (float) – optional, total spin angular momentum of state. By default 0.5 for Alkali atoms.
Returns:

radial coupling strength (in a.u.), or zero for forbidden transitions in dipole and quadrupole approximation.

Return type:

float

getRadialMatrixElement(n1, l1, j1, n2, l2, j2, s=None, useLiterature=True)[source]

Radial part of the dipole matrix element

Calculates \(\int \mathbf{d}r~R_{n_1,l_1,j_1}(r)\cdot R_{n_1,l_1,j_1}(r) \cdot r^3\).

Parameters:
  • n1 (int) – principal quantum number of state 1
  • l1 (int) – orbital angular momentum of state 1
  • j1 (float) – total angular momentum of state 1
  • n2 (int) – principal quantum number of state 2
  • l2 (int) – orbital angular momentum of state 2
  • j2 (float) – total angular momentum of state 2
  • s (float) – is required argument, total spin angular momentum of state. Specify s=0 for singlet state or s=1 for triplet state.
Returns:

dipole matrix element (\(a_0 e\)).

Return type:

float

getReducedMatrixElementJ(n1, l1, j1, n2, l2, j2, s=0.5)

Reduced matrix element in \(J\) basis (symmetric notation)

Parameters:
  • n1 (int) – principal quantum number of state 1
  • l1 (int) – orbital angular momentum of state 1
  • j1 (float) – total angular momentum of state 1
  • n2 (int) – principal quantum number of state 2
  • l2 (int) – orbital angular momentum of state 2
  • j2 (float) – total angular momentum of state 2
  • s (float) – optional, total spin angular momentum of state. By default 0.5 for Alkali atoms.
Returns:

reduced dipole matrix element in \(J\) basis \(\langle j || er || j' \rangle\) (\(a_0 e\)).

Return type:

float

getReducedMatrixElementJ_asymmetric(n1, l1, j1, n2, l2, j2, s=0.5)

Reduced matrix element in \(J\) basis, defined in asymmetric notation.

Note that notation for symmetric and asymmetricly defined reduced matrix element is not consistent in the literature. For example, notation is used e.g. in Steck [1] is precisely the oposite.

Note

Note that this notation is asymmetric: \(( j||e r ||j' ) \neq ( j'||e r ||j )\). Relation between the two notation is \(\langle j||er||j' \rangle=\sqrt{2j+1} ( j ||er ||j')\). This function always returns value for transition from lower to higher energy state, independent of the order of states entered in the function call.

Parameters:
  • n1 (int) – principal quantum number of state 1
  • l1 (int) – orbital angular momentum of state 1
  • j1 (float) – total angular momentum of state 1
  • n2 (int) – principal quantum number of state 2
  • l2 (int) – orbital angular momentum of state 2
  • j2 (float) – total angular momentum of state 2
  • s (float) – optional, total spin angular momentum of state. By default 0.5 for Alkali atoms.
Returns:

reduced dipole matrix element in Steck notation \(( j || er || j' )\) (\(a_0 e\)).

Return type:

float

[1]Daniel A. Steck, “Cesium D Line Data,” (revision 2.0.1, 2 May 2008). http://steck.us/alkalidata
getReducedMatrixElementL(n1, l1, j1, n2, l2, j2, s=0.5)

Reduced matrix element in \(L\) basis (symmetric notation)

Parameters:
  • n1 (int) – principal quantum number of state 1
  • l1 (int) – orbital angular momentum of state 1
  • j1 (float) – total angular momentum of state 1
  • n2 (int) – principal quantum number of state 2
  • l2 (int) – orbital angular momentum of state 2
  • j2 (float) – total angular momentum of state 2
Returns:

reduced dipole matrix element in \(L\) basis \(\langle l || er || l' \rangle\) (\(a_0 e\)).

Return type:

float

getStateLifetime(n, l, j, temperature=0, includeLevelsUpTo=0, s=0)[source]

Returns the lifetime of the state (in s)

For non-zero temperatures, user must specify up to which principal quantum number levels, that is above the initial state, should be included in order to account for black-body induced transitions to higher lying states. See Rydberg lifetimes example snippet.

Parameters:
  • l, j (n,) – specifies state whose lifetime we are calculating
  • temperature – optional. Temperature at which the atom environment is, measured in K. If this parameter is non-zero, user has to specify transitions up to which state (due to black-body decay) should be included in calculation.
  • includeLevelsUpTo (int) – optional and not needed for atom lifetimes calculated at zero temperature. At non zero temperatures, this specify maximum principal quantum number of the state to which black-body induced transitions will be included. Minimal value of the parameter in that case is \(n+1\)
  • s (float) – optional, total spin angular momentum of state. By default 0.5 for Alkali atoms.
Returns:

State lifetime in units of s (seconds)

Return type:

float

See also

getTransitionRate for calculating rates of individual transition rates between the two states

getTransitionFrequency(n1, l1, j1, n2, l2, j2, s=0.5, s2=None)

Calculated transition frequency in Hz

Returned values is given relative to the centre of gravity of the hyperfine-split states.

Parameters:
  • n1 (int) – principal quantum number of the state from which we are going
  • l1 (int) – orbital angular momentum of the state from which we are going
  • j1 (float) – total angular momentum of the state from which we are going
  • n2 (int) – principal quantum number of the state to which we are going
  • l2 (int) – orbital angular momentum of the state to which we are going
  • j2 (float) – total angular momentum of the state to which we are going
  • s (float) – optional, spin of the intial state (for Alkali this is fixed to 0.5)
  • s2 (float) – optional, spin of the final state If not set, defaults to the same value as s
Returns:

transition frequency (in Hz). If the returned value is negative, level from which we are going is above the level to which we are going.

Return type:

float

getTransitionRate(n1, l1, j1, n2, l2, j2, temperature=0.0, s=0.5)

Transition rate due to coupling to vacuum modes (black body included)

Calculates transition rate from the first given state to the second given state \(|n_1,l_1,j_1\rangle \rightarrow |n_2,j_2,j_2\rangle\) at given temperature due to interaction with the vacuum field. For zero temperature this returns Einstein A coefficient. For details of calculation see Ref. [3] and Ref. [4]. See Black-body induced population transfer example snippet.

Parameters:
  • n1 (int) – principal quantum number
  • l1 (int) – orbital angular momentum
  • j1 (float) – total angular momentum
  • n2 (int) – principal quantum number
  • l2 (int) – orbital angular momentum
  • j2 (float) – total angular momentum
  • [temperature] (float) – temperature in K
  • s (float) – optional, total spin angular momentum of state. By default 0.5 for Alkali atoms.
Returns:

transition rate in s \({}^{-1}\) (SI)

Return type:

float

References

[3]C. E. Theodosiou, PRA 30, 2881 (1984) https://doi.org/10.1103/PhysRevA.30.2881
[4]I. I. Beterov, I. I. Ryabtsev, D. B. Tretyakov, and V. M. Entin, PRA 79, 052504 (2009) https://doi.org/10.1103/PhysRevA.79.052504
getTransitionWavelength(n1, l1, j1, n2, l2, j2, s=0.5, s2=None)

Calculated transition wavelength (in vacuum) in m.

Returned values is given relative to the centre of gravity of the hyperfine-split states.

Parameters:
  • n1 (int) – principal quantum number of the state from which we are going
  • l1 (int) – orbital angular momentum of the state from which we are going
  • j1 (float) – total angular momentum of the state from which we are going
  • n2 (int) – principal quantum number of the state to which we are going
  • l2 (int) – orbital angular momentum of the state to which we are going
  • j2 (float) – total angular momentum of the state to which we are going
  • s (float) – optional, spin of the intial state (for Alkali this is fixed to 0.5)
  • s2 (float) – optional, spin of the final state. If not set, defaults to same value as s
Returns:

transition wavelength (in m). If the returned value is negative, level from which we are going is above the level to which we are going.

Return type:

float

getZeemanEnergyShift(l, j, mj, magneticFieldBz, s=0.5)

Retuns linear (paramagnetic) Zeeman shift.

\(\mathcal{H}_P=\frac{\mu_B B_z}{\hbar}(\hat{L}_{\rm z}+\ g_{\rm S}S_{\rm z})\)

Parameters:
  • l (int) – orbital angular momentum
  • j (float) – total angular momentum
  • mj (float) – projection of total angular momentum alon z-axis
  • magneticFieldBz (float) – applied magnetic field (alon z-axis only) in units of T (Tesla)
  • s (float) – optional, total spin angular momentum of state. By default 0.5 for Alkali atoms.
Returns:

energy offset of the state (in J)

Return type:

float

levelDataFromNIST = ''

file with .csv data, each row is [n, l, s, j, energy, source, absolute uncertanty]

minQuantumDefectN = None

Not used with DivalentAtom, see defectFittingRange instead.

modelPotential_coef = {}

Model potential parameters fitted from experimental observations for different l (electron angular momentum)

potential(l, s, j, r)[source]

Not implemented for Alkaline earths

quantumDefect = [[[0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0]], [[0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0]], [[0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0]], [[0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0]]]

Contains list of modified Rydberg-Ritz coefficients for calculating quantum defects for [[ \(^1S_{0},^1P_{1},^1D_{2},^1F_{3}\)], [ \(^3S_{0},^3P_{0},^3D_{1},^3F_{2}\)], [ \(^3S_{0},^3P_{1},^3D_{2},^3F_{3}\)], [ \(^3S_{1},^3P_{2},^3D_{3},^3F_{4}\)]].

radialWavefunction(l, s, j, stateEnergy, innerLimit, outerLimit, step)[source]

Not implemented for Alkaline earths

updateDipoleMatrixElementsFile()

Updates the file with pre-calculated dipole matrix elements.

This function will add the the file all the elements that have been calculated in the previous run, allowing quick access to them in the future calculations.