Governing Equations

The classes here correspond to the three governing equations.

Overview

governingeq.governingeq(laserBeams, magField)

Governing equation base class

heuristiceq(laserBeams, magField[, a, mass, ...])

Heuristic force equation

rateeq(laserBeams, magField, hamitlonian[, ...])

The rate equations

obe(laserBeams, magField, hamitlonian[, a, ...])

The optical Bloch equations

Details

class pylcp.governingeq.governingeq(laserBeams, magField, hamiltonian=None, a=array([0.0, 0.0, 0.0]), r0=array([0.0, 0.0, 0.0]), v0=array([0.0, 0.0, 0.0]))[source]

Governing equation base class

This class is the basis for making all the governing equations in pylcp, including the rate equations, heuristic equation, and the optical Bloch equations. Its methods are available to other governing equations.

Parameters
  • laserBeams (dictionary of pylcp.laserBeams, pylcp.laserBeams, or list of pylcp.laserBeam) –

    The laserBeams that will be used in constructing the optical Bloch equations. which transitions in the block diagonal hamiltonian. It can be any of the following:

    • A dictionary of pylcp.laserBeams: if this is the case, the keys of the dictionary should match available \(d^{nm}\) matrices in the pylcp.hamiltonian object. The key structure should be n->m.

    • pylcp.laserBeams: a single set of laser beams is assumed to address the transition g->e.

    • a list of pylcp.laserBeam: automatically promoted to a pylcp.laserBeams object assumed to address the transtion g->e.

  • magField (pylcp.magField or callable) – The function or object that defines the magnetic field.

  • hamiltonian (pylcp.hamiltonian or None) – The internal hamiltonian of the particle.

  • a (array_like, shape (3,), optional) – A default acceleraiton to apply to the particle’s motion, usually gravity. Default: [0., 0., 0.]

  • r0 (array_like, shape (3,)) – Initial position. Default: [0.,0.,0.]

  • v0 (array_like, shape (3,)) – Initial velocity. Default: [0.,0.,0.]

damping_coeff(axes, r=None, eps=0.01, **kwargs)[source]

Find the damping coefficent

Uses the find_equilibrium force() method to calculate the damping coefficient for the particular configuration.

Parameters
  • axes (array_like) – A list of axis indices to compute the damping coefficient(s) along. Here, \(\hat{x}\) is index 0, \(\hat{y}\) is index 1, and \(\hat{z}\) is index 2. For example, axes=[2] calculates the damping parameter along \(\hat{z}\).

  • r (array_like, optional) – The position at which to calculate the damping coefficient. By default r=None, which defaults to calculating at the equilibrium position as found by the find_equilibrium_position() method. If this method has not been run, it defaults to the origin.

  • eps (float) – The small numerical \(\epsilon\) parameter used for calculating the \(df/dv\) derivative. Default: 0.01

  • kwargs – Any additional keyword arguments to pass to find_equilibrium_force()

Returns

beta – The damping coefficients along the selected axes.

Return type

list or float

find_equilibrium_force()[source]

Find the equilibrium force at the initial conditions

Returns

force – Equilibrium force experienced by the atom

Return type

array_like

find_equilibrium_position(axes, **kwargs)[source]

Find the equilibrium position

Uses the find_equilibrium force() method to calculate the where the \(\mathbf{f}(\mathbf{r}, \mathbf{v}=0)=0\).

Parameters
  • axes (array_like) – A list of axis indices to compute the trapping frequencies along. Here, \(\hat{x}\) is index 0, \(\hat{y}\) is index 1, and \(\hat{z}\) is index 2. For example, axes=[2] calculates the trapping frquency along \(\hat{z}\).

  • kwargs – Any additional keyword arguments to pass to find_equilibrium_force()

Returns

r_eq – The equilibrium positions along the selected axes.

Return type

list or float

force()[source]

Find the instantaneous force

Returns

force – Force experienced by the atom

Return type

array_like

generate_force_profile()[source]

Map out the equilibrium force vs. position and velocity

Parameters
  • R (array_like, shape(3, ...)) – Position vector. First dimension of the array must be length 3, and corresponds to \(x\), \(y\), and \(z\) components, repsectively.

  • V (array_like, shape(3, ...)) – Velocity vector. First dimension of the array must be length 3, and corresponds to \(v_x\), \(v_y\), and \(v_z\) components, repsectively.

  • name (str, optional) – Name for the profile. Stored in profile dictionary in this object. If None, uses the next integer, cast as a string, (i.e., ‘0’) as the name.

  • progress_bar (boolean, optional) – Displays a progress bar as the proceeds. Default: False

Returns

profile – Resulting force profile.

Return type

pylcp.common.base_force_profile

set_initial_position(r0)[source]

Sets the initial position

Parameters

r0 (array_like, shape (3,)) – Initial position. Default: [0.,0.,0.]

set_initial_position_and_velocity(r0, v0)[source]

Sets the initial position and velocity

Parameters
  • r0 (array_like, shape (3,)) – Initial position. Default: [0.,0.,0.]

  • v0 (array_like, shape (3,)) – Initial velocity. Default: [0.,0.,0.]

set_initial_velocity(v0)[source]

Sets the initial velocity

Parameters

v0 (array_like, shape (3,)) – Initial position. Default: [0.,0.,0.]

trapping_frequencies(axes, r=None, eps=0.01, **kwargs)[source]

Find the trapping frequency

Uses the find_equilibrium force() method to calculate the trapping frequency for the particular configuration.

Parameters
  • axes (array_like) – A list of axis indices to compute the trapping frequencies along. Here, \(\hat{x}\) is index 0, \(\hat{y}\) is index 1, and \(\hat{z}\) is index 2. For example, axes=[2] calculates the trapping frquency along \(\hat{z}\).

  • r (array_like, optional) – The position at which to calculate the damping coefficient. By default r=None, which defaults to calculating at the equilibrium position as found by the find_equilibrium_position() method. If this method has not been run, it defaults to the origin.

  • eps (float, optional) – The small numerical \(\epsilon\) parameter used for calculating the \(df/dr\) derivative. Default: 0.01

  • kwargs – Any additional keyword arguments to pass to find_equilibrium_force()

Returns

omega – The trapping frequencies along the selected axes.

Return type

list or float

class pylcp.heuristiceq(laserBeams, magField, a=array([0.0, 0.0, 0.0]), mass=100, gamma=1, k=1, r0=array([0.0, 0.0, 0.0]), v0=array([0.0, 0.0, 0.0]))[source]

Heuristic force equation

The heuristic equation governs the atom or molecule as if it has a single transition between an \(F=0\) ground state to an \(F'=1\) excited state.

Parameters
  • laserBeams (dictionary of pylcp.laserBeams, pylcp.laserBeams, or list of pylcp.laserBeam) –

    The laserBeams that will be used in constructing the optical Bloch equations. which transitions in the block diagonal hamiltonian. It can be any of the following:

    • A dictionary of pylcp.laserBeams: if this is the case, the keys of the dictionary should match available \(d^{nm}\) matrices in the pylcp.hamiltonian object. The key structure should be n->m. Here, it must be g->e.

    • pylcp.laserBeams: a single set of laser beams is assumed to address the transition g->e.

    • a list of pylcp.laserBeam: automatically promoted to a pylcp.laserBeams object assumed to address the transtion g->e.

  • magField (pylcp.magField or callable) – The function or object that defines the magnetic field.

  • hamiltonian (pylcp.hamiltonian) – The internal hamiltonian of the particle.

  • a (array_like, shape (3,), optional) – A default acceleraiton to apply to the particle’s motion, usually gravity. Default: [0., 0., 0.]

  • r0 (array_like, shape (3,), optional) – Initial position of the atom or molecule. Default: [0., 0., 0.]

  • v0 (array_like, shape (3,), optional) – Initial velocity of the atom or molecule. Default: [0., 0., 0.]

  • mass (float, optional) – Mass of the atom or molecule. Default: 100

  • gamma (float, optional) – Decay rate of the single transition in the atom or molecule. Default: 1

  • k (float, optional) – Magnitude of the k vector for the single transition in the atom or molecule. Default: 1

evolve_motion(t_span, freeze_axis=[False, False, False], random_recoil=False, random_force=False, max_scatter_probability=0.1, progress_bar=False, rng=Generator(PCG64) at 0x7F6787530A50, **kwargs)[source]

Evolve the motion of the atom in time.

Parameters
  • t_span (list or array_like) – A two element list or array that specify the initial and final time of integration.

  • freeze_axis (list of boolean) – Freeze atomic motion along the specified axis. Default: [False, False, False]

  • random_recoil (boolean) – Allow the atom to randomly recoil from scattering events. Default: False

  • random_force (boolean) – Rather than calculating the force using the heuristieq.force() method, use the calculated scattering rates from each of the laser beam to randomly add photon absorption events that cause the atom to recoil randomly from the laser beam(s). Default: False

  • max_scatter_probability (float) – When undergoing random recoils, this sets the maximum time step such that the maximum scattering probability is less than or equal to this number during the next time step. Default: 0.1

  • progress_bar (boolean) – If true, show a progress bar as the calculation proceeds. Default: False

  • rng (numpy.random.Generator()) – A properly-seeded random number generator. Default: calls numpy.random.default.rng()

  • **kwargs – Additional keyword arguments get passed to solve_ivp_random, which is what actually does the integration.

Returns

sol

Bunch object that contains the following fields:

  • t: integration times found by solve_ivp

  • v: atomic velocity

  • r: atomic position

It contains other important elements, which can be discerned from scipy’s solve_ivp documentation.

Return type

OdeSolution

find_equilibrium_force(return_details=False)[source]

Finds the equilibrium force at the initial position

Parameters

return_details (boolean, optional) – If True, returns the forces from each laser and the scattering rate matrix.

Returns

  • F (array_like) – total equilibrium force experienced by the atom

  • F_laser (array_like) – If return_details is True, the forces due to each laser.

  • R (array_like) – The scattering rate matrix.

force(r, v, t)[source]

Calculates the instantaneous force

Parameters
  • r (array_like) – Position at which to calculate the force

  • v (array_like) – Velocity at which to calculate the force

  • t (float) – Time at which to calculate the force

Returns

  • F (array_like) – total equilibrium force experienced by the atom

  • F_laser (dictionary of array_like) – If return_details is True, the forces due to each laser, indexed by the manifold the laser addresses. The dictionary is keyed by the transition driven, and individual lasers are in the same order as in the pylcp.laserBeams object used to create the governing equation.

generate_force_profile(R, V, name=None, progress_bar=False)[source]

Map out the equilibrium force vs. position and velocity

Parameters
  • R (array_like, shape(3, ...)) – Position vector. First dimension of the array must be length 3, and corresponds to \(x\), \(y\), and \(z\) components, repsectively.

  • V (array_like, shape(3, ...)) – Velocity vector. First dimension of the array must be length 3, and corresponds to \(v_x\), \(v_y\), and \(v_z\) components, repsectively.

  • name (str, optional) – Name for the profile. Stored in profile dictionary in this object. If None, uses the next integer, cast as a string, (i.e., ‘0’) as the name.

  • progress_bar (boolean, optional) – Displays a progress bar as the proceeds. Default: False

Returns

profile – Resulting force profile.

Return type

pylcp.common.base_force_profile

scattering_rate(r, v, t, return_kvecs=False)[source]

Calculates the scattering rate

Parameters
  • r (array_like) – Position at which to calculate the force

  • v (array_like) – Velocity at which to calculate the force

  • t (float) – Time at which to calculate the force

  • return_kvecs (bool) – If true, returns both the scattering rate and the k-vecotrs from the lasers.

Returns

  • R (array_like) – Array of scattering rates associated with the lasers driving the transition.

  • kvecs (array_like) – If return_kvecs is True, the k-vectors of each of the lasers. This is used in heuristiceq.force, where it calls this function to calculate the scattering rate first. By returning the k-vectors with the scattering rates, it prevents the need of having to recompute the k-vectors again.

class pylcp.common.base_force_profile(R, V, laserBeams, hamiltonian)[source]

Base force profile

The force profile object stores all of the calculated quantities created by the governingeq.generate_force_profile() method. It has the following attributes:

R

Positions at which the force profile was calculated.

Type

array_like, shape (3, …)

V

Velocities at which the force profile was calculated.

Type

array_like, shape (3, …)

F

Total equilibrium force at position R and velocity V.

Type

array_like, shape (3, …)

f_mag

Magnetic force at position R and velocity V.

Type

array_like, shape (3, …)

f

The forces due to each laser, indexed by the manifold the laser addresses. The dictionary is keyed by the transition driven, and individual lasers are in the same order as in the pylcp.laserBeams object used to create the governing equation.

Type

dictionary of array_like

Neq

Equilibrium population found.

Type

array_like

class pylcp.rateeq(laserBeams, magField, hamitlonian, a=array([0.0, 0.0, 0.0]), include_mag_forces=True, svd_eps=1e-10, r0=array([0.0, 0.0, 0.0]), v0=array([0.0, 0.0, 0.0]))[source]

The rate equations

This class constructs the rate equations from the given laser beams, magnetic field, and hamiltonian.

Parameters
  • laserBeams (dictionary of pylcp.laserBeams, pylcp.laserBeams, or list of pylcp.laserBeam) –

    The laserBeams that will be used in constructing the optical Bloch equations. which transitions in the block diagonal hamiltonian. It can be any of the following:

    • A dictionary of pylcp.laserBeams: if this is the case, the keys of the dictionary should match available \(d^{nm}\) matrices in the pylcp.hamiltonian object. The key structure should be n->m.

    • pylcp.laserBeams: a single set of laser beams is assumed to address the transition g->e.

    • a list of pylcp.laserBeam: automatically promoted to a pylcp.laserBeams object assumed to address the transtion g->e.

  • magField (pylcp.magField or callable) – The function or object that defines the magnetic field.

  • hamiltonian (pylcp.hamiltonian) – The internal hamiltonian of the particle.

  • a (array_like, shape (3,), optional) – A default acceleraiton to apply to the particle’s motion, usually gravity. Default: [0., 0., 0.]

  • include_mag_forces (boolean) – Optional flag to inculde magnetic forces in the force calculation. Default: True

  • r0 (array_like, shape (3,), optional) – Initial position. Default: [0., 0., 0.]

  • v0 (array_like, shape (3,), optional) – Initial velocity. Default: [0., 0., 0.]

construct_evolution_matrix(r, v, t=0.0, default_axis=array([0.0, 0.0, 1.0]))[source]

Constructs the evolution matrix at a given position and time.

Parameters
  • r (array_like, shape (3,)) – Position at which to calculate the equilibrium population

  • v (array_like, shape (3,)) – Velocity at which to calculate the equilibrium population

  • t (float) – Time at which to calculate the equilibrium population

equilibrium_populations(r, v, t, **kwargs)[source]

Returns the equilibrium population as determined by the rate equations

This method uses singular matrix decomposition to find the equilibrium state of the rate equations at a given position, velocity, and time.

Parameters
  • r (array_like, shape (3,)) – Position at which to calculate the equilibrium population

  • v (array_like, shape (3,)) – Velocity at which to calculate the equilibrium population

  • t (float) – Time at which to calculate the equilibrium population

  • return_details (boolean, optional) – In addition to the equilibrium populations, return the full population evolution matrix and the scattering rates for each of the lasers

Returns

  • Neq (array_like) – Equilibrium population vector

  • Rev (array_like) – If return details is True, the evolution matrix for the state populations.

  • Rijl (dictionary of array_like) – If return details is True, the scattering rates for each laser and each combination of states between the manifolds specified by the dictionary’s index.

evolve_motion(t_span, freeze_axis=[False, False, False], random_recoil=False, random_force=False, max_scatter_probability=0.1, progress_bar=False, record_force=False, rng=Generator(PCG64) at 0x7F6787530E50, **kwargs)[source]

Evolve the populations \(N\) and the motion of the atom in time.

This function evolves the rate equations, moving the atom through space, given the instantaneous force, for some period of time.

Parameters
  • t_span (list or array_like) – A two element list or array that specify the initial and final time of integration.

  • freeze_axis (list of boolean) – Freeze atomic motion along the specified axis. Default: [False, False, False]

  • random_recoil (boolean) – Allow the atom to randomly recoil from scattering events. Default: False

  • random_force (boolean) – Rather than calculating the force using the rateeq.force() method, use the calculated scattering rates from each of the laser beam (combined with the instantaneous populations) to randomly add photon absorption events that cause the atom to recoil randomly from the laser beam(s). Default: False

  • max_scatter_probability (float) – When undergoing random recoils and/or force, this sets the maximum time step such that the maximum scattering probability is less than or equal to this number during the next time step. Default: 0.1

  • progress_bar (boolean) – If true, show a progress bar as the calculation proceeds. Default: False

  • record_force (boolean) – If true, record the instantaneous force and store in the solution. Default: False

  • rng (numpy.random.Generator()) – A properly-seeded random number generator. Default: calls numpy.random.default.rng()

  • **kwargs – Additional keyword arguments get passed to solve_ivp_random, which is what actually does the integration.

Returns

sol

Bunch object that contains the following fields:

  • t: integration times found by solve_ivp

  • N: population vs. time

  • v: atomic velocity

  • r: atomic position

It contains other important elements, which can be discerned from scipy’s solve_ivp documentation.

Return type

OdeSolution

evolve_populations(t_span, **kwargs)[source]

Evolve the state population in time.

This function integrates the rate equations to determine how the populations evolve in time. Any initial velocity is kept constant. It is analogous to obe.evolve_density().

Parameters
  • t_span (list or array_like) – A two element list or array that specify the initial and final time of integration.

  • **kwargs – Additional keyword arguments get passed to solve_ivp, which is what actually does the integration.

Returns

sol

Bunch object that contains the following fields:

  • t: integration times found by solve_ivp

  • rho: density matrix

  • v: atomic velocity (constant)

  • r: atomic position

It contains other important elements, which can be discerned from scipy’s solve_ivp documentation.

Return type

OdeSolution

find_equilibrium_force(return_details=False, **kwargs)[source]

Finds the equilibrium force at the initial position

This method works by finding the equilibrium population through the rateeq.equilibrium_population() function, then calculating the resulting force.

Parameters
  • return_details (boolean, optional) – If true, returns the forces from each laser and the scattering rate matrix. Default: False

  • kwargs – Any additional keyword arguments to be passed to equilibrium_populations()

Returns

  • F (array_like) – total equilibrium force experienced by the atom

  • F_laser (array_like) – If return_details is True, the forces due to each laser.

  • Neq (array_like) – If return_details is True, the equilibrium populations.

  • Rijl (dictionary of array_like) – If return details is True, the scattering rates for each laser and each combination of states between the manifolds specified by the dictionary’s index.

  • F_mag (array_like) – If return_details is True, the forces due to the magnetic field.

  • ii (int) – Number of iterations needed to converge.

force(r, t, N, return_details=True)[source]

Calculates the instantaneous force

Parameters
  • r (array_like) – Position at which to calculate the force

  • t (float) – Time at which to calculate the force

  • N (array_like) – Relative state populations

  • return_details (boolean, optional) – If True, returns the forces from each laser and the magnetic forces.

Returns

  • F (array_like) – total equilibrium force experienced by the atom

  • F_laser (dictionary of array_like) – If return_details is True, the forces due to each laser, indexed by the manifold the laser addresses. The dictionary is keyed by the transition driven, and individual lasers are in the same order as in the pylcp.laserBeams object used to create the governing equation.

  • F_mag (array_like) – If return_details is True, the forces due to the magnetic field.

generate_force_profile(R, V, name=None, progress_bar=False, **kwargs)[source]

Map out the equilibrium force vs. position and velocity

Parameters
  • R (array_like, shape(3, ...)) – Position vector. First dimension of the array must be length 3, and corresponds to \(x\), \(y\), and \(z\) components, repsectively.

  • V (array_like, shape(3, ...)) – Velocity vector. First dimension of the array must be length 3, and corresponds to \(v_x\), \(v_y\), and \(v_z\) components, repsectively.

  • name (str, optional) – Name for the profile. Stored in profile dictionary in this object. If None, uses the next integer, cast as a string, (i.e., ‘0’) as the name.

  • progress_bar (boolean, optional) – Displays a progress bar as the proceeds. Default: False

  • kwargs – Any additional keyword arguments to be passed to rateeq.find_equilibrium_force()

Returns

profile – Resulting force profile.

Return type

pylcp.rateeq.force_profile

set_initial_pop(N0)[source]

Sets the initial populations

Parameters

N0 (array_like) – The initial state population vector \(N_0\). It must have \(n\) elements, where \(n\) is the total number of states in the system.

set_initial_pop_from_equilibrium()[source]

Sets the initial populations based on the equilibrium population at the initial position and velocity and time t=0.

class pylcp.rateeq.force_profile(R, V, laserBeams, hamiltonian)[source]

Rate equation force profile

The force profile object stores all of the calculated quantities created by the rateeq.generate_force_profile() method. It has following attributes:

R

Positions at which the force profile was calculated.

Type

array_like, shape (3, …)

V

Velocities at which the force profile was calculated.

Type

array_like, shape (3, …)

F

Total equilibrium force at position R and velocity V.

Type

array_like, shape (3, …)

f_mag

Magnetic force at position R and velocity V.

Type

array_like, shape (3, …)

f

The forces due to each laser, indexed by the manifold the laser addresses. The dictionary is keyed by the transition driven, and individual lasers are in the same order as in the pylcp.laserBeams object used to create the governing equation.

Type

dictionary of array_like

Neq

Equilibrium population found.

Type

array_like

Rijl

The pumping rates of each laser, indexed by the manifold the laser addresses. The dictionary is keyed by the transition driven, and individual lasers are in the same order as in the pylcp.laserBeams object used to create the governing equation.

Type

dictionary of array_like

class pylcp.obe(laserBeams, magField, hamitlonian, a=array([0.0, 0.0, 0.0]), transform_into_re_im=True, use_sparse_matrices=None, include_mag_forces=True, r0=array([0.0, 0.0, 0.0]), v0=array([0.0, 0.0, 0.0]))[source]

The optical Bloch equations

This class constructs the optical Bloch equations from the given laser beams, magnetic field, and hamiltonian.

Parameters
  • laserBeams (dictionary of pylcp.laserBeams, pylcp.laserBeams, or list of pylcp.laserBeam) –

    The laserBeams that will be used in constructing the optical Bloch equations. which transitions in the block diagonal hamiltonian. It can be any of the following:

    • A dictionary of pylcp.laserBeams: if this is the case, the keys of the dictionary should match available \(d^{nm}\) matrices in the pylcp.hamiltonian object. The key structure should be n->m.

    • pylcp.laserBeams: a single set of laser beams is assumed to address the transition g->e.

    • a list of pylcp.laserBeam: automatically promoted to a pylcp.laserBeams object assumed to address the transtion g->e.

  • magField (pylcp.magField or callable) – The function or object that defines the magnetic field.

  • hamiltonian (pylcp.hamiltonian) – The internal hamiltonian of the particle.

  • a (array_like, shape (3,), optional) – A default acceleraiton to apply to the particle’s motion, usually gravity. Default: [0., 0., 0.]

  • transform_into_re_im (boolean) – Optional flag to transform the optical Bloch equations into real and imaginary components. This helps to decrease computaiton time as it uses the symmetry \(\rho_{ji}=\rho_{ij}^*\) to cut the number of equations nearly in half. Default: True

  • use_sparse_matrices (boolean or None) – Optional flag to use sparse matrices. If none, it will use sparse matrices only if the number of internal states > 10, which would result in the evolution matrix for the density operators being a 100x100 matrix. At that size, there may be some speed up with sparse matrices. Default: None

  • include_mag_forces (boolean) – Optional flag to inculde magnetic forces in the force calculation. Default: True

  • r0 (array_like, shape (3,), optional) – Initial position. Default: [0., 0., 0.]

  • v0 (array_like, shape (3,), optional) – Initial velocity. Default: [0., 0., 0.]

evolve_density(t_span, progress_bar=False, **kwargs)[source]

Evolve the density operators \(\rho_{ij}\) in time.

This function integrates the optical Bloch equations to determine how the populations evolve in time. Any initial velocity is kept constant while the atom potentially moves through the light field. This function is therefore useful in determining average forces. Any constant acceleration set when the OBEs were generated is ignored. It is analogous to rateeq.evolve_populations().

Parameters
  • t_span (list or array_like) – A two element list or array that specify the initial and final time of integration.

  • progress_bar (boolean) – Show a progress bar as the calculation proceeds. Default:False

  • **kwargs – Additional keyword arguments get passed to solve_ivp, which is what actually does the integration.

Returns

sol

Bunch object that contains the following fields:

  • t: integration times found by solve_ivp

  • rho: density matrix

  • v: atomic velocity (constant)

  • r: atomic position

It contains other important elements, which can be discerned from scipy’s solve_ivp documentation.

Return type

OdeSolution

evolve_motion(t_span, freeze_axis=[False, False, False], random_recoil=False, max_scatter_probability=0.1, progress_bar=False, record_force=False, rng=Generator(PCG64) at 0x7F678622EE50, **kwargs)[source]

Evolve \(\rho_{ij}\) and the motion of the atom in time.

This function evolves the optical Bloch equations, moving the atom along given the instantaneous force, for some period of time.

Parameters
  • t_span (list or array_like) – A two element list or array that specify the initial and final time of integration.

  • freeze_axis (list of boolean) – Freeze atomic motion along the specified axis. Default: [False, False, False]

  • random_recoil (boolean) – Allow the atom to randomly recoil from scattering events. Default: False

  • max_scatter_probability (float) – When undergoing random recoils, this sets the maximum time step such that the maximum scattering probability is less than or equal to this number during the next time step. Default: 0.1

  • progress_bar (boolean) – If true, show a progress bar as the calculation proceeds. Default: False

  • record_force (boolean) – If true, record the instantaneous force and store in the solution. Default: False

  • rng (numpy.random.Generator()) – A properly-seeded random number generator. Default: calls numpy.random.default.rng()

  • **kwargs – Additional keyword arguments get passed to solve_ivp_random, which is what actually does the integration.

Returns

sol

Bunch object that contains the following fields:

  • t: integration times found by solve_ivp

  • rho: density matrix

  • v: atomic velocity

  • r: atomic position

It contains other important elements, which can be discerned from scipy’s solve_ivp documentation.

Return type

OdeSolution

find_equilibrium_force(deltat=500, itermax=100, Npts=5001, rel=1e-05, abs=1e-09, debug=False, initial_rho='rateeq', return_details=False, **kwargs)[source]

Finds the equilibrium force at the initial position

This method works by solving the OBEs in a chunk of time \(\Delta T\), calculating the force during that chunck, continuing the integration for another chunck, calculating the force during that subsequent chunck, and comparing the average of the forces of the two chunks to see if they have converged.

Parameters
  • deltat (float) – Chunk time \(\Delta T\). Default: 500

  • itermax (int, optional) – Maximum number of iterations. Default: 100

  • Npts (int, optional) – Number of points to divide the chunk into. Default: 5001

  • rel (float, optional) – Relative convergence parameter. Default: 1e-5

  • abs (float, optional) – Absolute convergence parameter. Default: 1e-9

  • debug (boolean, optional) – If true, pring out debug information as it goes.

  • initial_rho ('rateeq' or 'equally') – Determines how to set the initial rho at the start of the calculation.

  • return_details (boolean, optional) – If true, returns the forces from each laser and the scattering rate matrix.

Returns

  • F (array_like) – total equilibrium force experienced by the atom

  • F_laser (dictionary of array_like) – If return_details is True, the forces due to each laser, indexed by the manifold the laser addresses. The dictionary is keyed by the transition driven, and individual lasers are in the same order as in the pylcp.laserBeams object used to create the governing equation.

  • F_laser (dictionary of array_like) – If return_details is True, the forces due to each laser and its q component, indexed by the manifold the laser addresses. The dictionary is keyed by the transition driven, and individual lasers are in the same order as in the pylcp.laserBeams object used to create the governing equation.

  • F_mag (array_like) – If return_details is True, the forces due to the magnetic field.

  • Neq (array_like) – If return_details is True, the equilibrium populations.

  • ii (int) – Number of iterations needed to converge.

force(r, t, rho, return_details=False)[source]

Calculates the instantaneous force

Parameters
  • r (array_like) – Position at which to calculate the force

  • t (float) – Time at which to calculate the force

  • rho (array_like) – Density matrix with which to calculate the force

  • return_details (boolean, optional) – If true, returns the forces from each laser and the scattering rate matrix.

Returns

  • F (array_like) – total equilibrium force experienced by the atom

  • F_laser (array_like) – If return_details is True, the forces due to each laser.

  • F_laser_q (array_like) – If return_details is True, the forces due to each laser and it’s q component of the polarization.

  • F_mag (array_like) – If return_details is True, the forces due to the magnetic field.

full_OBE_ev(r, t)[source]

Calculate the evolution for the density matrix

This function calculates the OBE evolution matrix by assembling pre-stored versions of the component matries. This should be significantly faster than full_OBE_ev_scratch, but it may suffer bugs in the evolution that full_OBE_ev_scratch will not. If Bq is None, it will compute Bq based on r, t

Parameters
  • r (array_like, shape (3,)) – Position at which to calculate evolution matrix

  • t (float) – Time at which to calculate evolution matrix

Returns

ev_mat – Evolution matrix for the densities

Return type

array_like

full_OBE_ev_scratch(r, t)[source]

Calculate the evolution for the density matrix

This function calculates the OBE evolution matrix at position t and r from scratch, first computing the full Hamiltonian, then the OBE evolution matrix computed via commutators, then adding in the decay matrix evolution. If Bq is None, it will compute Bq.

Parameters
  • r (array_like, shape (3,)) – Position at which to calculate evolution matrix

  • t (float) – Time at which to calculate evolution matrix

Returns

ev_mat – Evolution matrix for the densities

Return type

array_like

generate_force_profile(R, V, name=None, progress_bar=False, **kwargs)[source]

Map out the equilibrium force vs. position and velocity

Parameters
  • R (array_like, shape(3, ...)) – Position vector. First dimension of the array must be length 3, and corresponds to \(x\), \(y\), and \(z\) components, repsectively.

  • V (array_like, shape(3, ...)) – Velocity vector. First dimension of the array must be length 3, and corresponds to \(v_x\), \(v_y\), and \(v_z\) components, repsectively.

  • name (str, optional) – Name for the profile. Stored in profile dictionary in this object. If None, uses the next integer, cast as a string, (i.e., ‘0’) as the name.

  • progress_bar (boolean, optional) – Displays a progress bar as the proceeds. Default: False

Returns

profile – Resulting force profile.

Return type

pylcp.obe.force_profile

observable(O, rho=None)[source]

Observable returns the obervable O given density matrix rho.

Parameters
  • O (array or array-like) – The matrix form of the observable operator. Can have any shape, representing scalar, vector, or tensor operators, but the last two axes must correspond to the matrix of the operator and have the same dimensions of the generating Hamiltonian. For example, a vector operator might have the shape (3, n, n), where n is the number of states and the first axis corresponds to x, y, and z.

  • rho ([optional] array or array-like) –

    The density matrix. The first two dimensions must have sizes (n, n), but there may be multiple instances of the density matrix tiled in the higher dimensions. For example, a rho with (n, n, m) could have m instances of the density matrix at different times.

    If not specified, will get rho from the current solution stored in memory.

Returns

observable – observable has shape (O[:-2])+(rho[2:])

Return type

float or array

set_initial_rho(rho0)[source]

Sets the initial \(\rho\) matrix

Parameters

rho0 (array_like) – The initial \(\rho\). It must have \(n^2\) elements, where \(n\) is the total number of states in the system. If a flat array, it will be reshaped.

set_initial_rho_equally()[source]

Sets the initial \(\rho\) matrix such that all states have the same population.

set_initial_rho_from_populations(Npop)[source]

Sets the diagonal elements of the initial \(\rho\) matrix

Parameters

Npop (array_like) – Array of the initial populations of the states in the system. The length must be \(n\), where \(n\) is the number of states.

set_initial_rho_from_rateeq()[source]

Sets the diagonal elements of the initial \(\rho\) matrix using the equilibrium populations as determined by pylcp.rateeq

class pylcp.obe.force_profile(R, V, laserBeams, hamiltonian)[source]

Optical Bloch equation force profile

The force profile object stores all of the calculated quantities created by the rateeq.generate_force_profile() method. It has the following attributes:

R

Positions at which the force profile was calculated.

Type

array_like, shape (3, …)

V

Velocities at which the force profile was calculated.

Type

array_like, shape (3, …)

F

Total equilibrium force at position R and velocity V.

Type

array_like, shape (3, …)

f_mag

Magnetic force at position R and velocity V.

Type

array_like, shape (3, …)

f

The forces due to each laser, indexed by the manifold the laser addresses. The dictionary is keyed by the transition driven, and individual lasers are in the same order as in the pylcp.laserBeams object used to create the governing equation.

Type

dictionary of array_like

f_q

The force due to each laser and its \(q\) component, indexed by the manifold the laser addresses. The dictionary is keyed by the transition driven, and individual lasers are in the same order as in the pylcp.laserBeams object used to create the governing equation.

Type

dictionary of array_like

Neq

Equilibrium population found.

Type

array_like