StrongFieldDynamics.jl

Appearance

StrongFieldDynamics.StrongFieldDynamics Module
julia
StrongFieldDynamics

A Julia package for strong field atomic physics calculations using the Strong Field Approximation (SFA).

This package provides tools for:

  • Computing electron wavefunctions in strong laser fields

  • Calculating photoionization and high harmonic generation amplitudes

  • Analyzing pulse shapes and field integrals

  • Unit conversions for atomic physics calculations

  • Visualization of strong field dynamics

Main Components

  • Amplitude calculations for various strong field processes

  • Electron wavefunction computations in intense laser fields

  • Pulse shape analysis and field integration utilities

  • Atomic units and unit conversion functions

  • Plotting utilities for visualization

Example

julia
using StrongFieldDynamics
source
StrongFieldDynamics.AngularDistribution Type
julia
AngularDistribution

Results from angle-resolved photoelectron distribution at fixed energy and theta.

Fields

  • energy::Float64: Fixed photoelectron energy (a.u.)

  • θ::Float64: Fixed polar angle (radians)

  • ϕ::Vector{Float64}: Azimuthal angle grid (radians)

  • distribution::Vector{Float64}: Angular distribution P(ϕ) at fixed θ

  • pulse::Pulse: Laser pulse

source
StrongFieldDynamics.AtomicElectron Type
julia
AtomicElectron

Represents a bound atomic electron with quantum numbers and radial wavefunction.

Fields

  • Z::Int64: Atomic number

  • n::Int64: Principal quantum number

  • l::Int64: Orbital angular momentum quantum number

  • j::Rational{Int64}: Total angular momentum quantum number

  • ε::Float64: Binding energy (positive value in a.u.)

  • r::Vector{Float64}: Radial grid points (a.u.)

  • P::Vector{Float64}: Radial wavefunction P(r) = r·R(r)

source
StrongFieldDynamics.ContinuumElectron Type
julia
ContinuumElectron

Defines a continuum (photo) electron for a particular energy.

Fields

  • ε::Float64

  • p::Float64

  • solution::ContinuumSolution

source
StrongFieldDynamics.ContinuumSolution Type
julia
ContinuumSolution

Enumeration of different approximations for the continuum electron wavefunction in strong-field ionization calculations.

Variants

  • Bessel

  • Distorted

Bessel

Volkov states - Free electron wavefunction in an oscillating electric field.

Distorted

Distorted waves - Numerical solution including both Coulomb and laser field effects.

Usage

julia
# Create continuum electrons with different wavefunctions
bessel_electron = ContinuumElectron(1.0, sqrt(2.0), Bessel)       # SFA calculation
distorted_electron = ContinuumElectron(1.0, sqrt(2.0), Distorted) # Full calculation
source
StrongFieldDynamics.EnergyDistribution Type
julia
EnergyDistribution

Results from energy-resolved photoelectron distribution calculation at fixed angles.

Fields

  • θ::Float64: Polar angle (radians)

  • ϕ::Float64: Azimuthal angle (radians)

  • energies::Vector{Float64}: Energy grid (a.u.)

  • distribution::Vector{Float64}: Differential ionization probability d²P/dΩdE

  • pulse::Pulse: Laser pulse

source
StrongFieldDynamics.MomentumDistribution Type
julia
MomentumDistribution

Results from momentum-resolved photoelectron distribution calculation in spherical coordinates.

Fields

  • p::Vector{Float64}: Momentum magnitude grid (a.u.)

  • θ::Vector{Float64}: Polar angle grid (radians, 0 to π)

  • φ::Vector{Float64}: Azimuthal angle grid (radians, 0 to 2π)

  • distribution::Array{Float64,3}: 3D momentum distribution P(p,θ,φ)

  • pulse::Pulse: Laser pulse

source
StrongFieldDynamics.PartialWave Type
julia
PartialWave

Represents a partial wave component of the continuum electron wavefunction.

Fields

  • ε::Float64: Kinetic energy (a.u.)

  • l::Int64: Orbital angular momentum quantum number

  • j::Rational{Int64}: Total angular momentum quantum number

  • P::Vector{Float64}: Radial wavefunction component

  • δ::Float64: Phase shift (radians)

source
StrongFieldDynamics.Pulse Type
julia
Pulse

Represents a laser pulse with all necessary parameters for strong-field calculations. All internal values are stored in atomic units.

Fields

  • I₀::Float64: Peak intensity (W/cm²) - stored in atomic units

  • A₀::Float64: Peak vector potential amplitude (a.u.)

  • λ::Float64: Wavelength (nm) - stored in atomic units

  • ω::Float64: Angular frequency (a.u.)

  • cycles::Int64: Number of optical cycles

  • Tp::Float64: Pulse duration (a.u.)

  • Up::Float64: Ponderomotive energy (a.u.)

  • f::Function: Envelope function f(t)

  • cep::Float64: Carrier-envelope phase (radians)

  • helicity::Int64: Helicity (+1 or -1)

  • ϵ::Float64: Ellipticity parameter (0=linear, 1=circular)

  • u::OffsetVector{Float64, Vector{Float64}}: Polarization unit vector

  • Sv::Function: Volkov-phase function

Required Arguments

  • I₀: Peak intensity (can include units, e.g., 1e14u"W/cm^2", or numeric value assuming W/cm²)

  • λ: Wavelength (can include units, e.g., 800u"nm", or numeric value assuming nm)

  • cycles: Number of optical cycles (positive integer)

  • envelope: Envelope type (:sin2, :gauss, :flat) or custom function f(t)

  • cep: Carrier-envelope phase in radians (default: 0.0)

  • helicity: Helicity (+1 or -1, default: +1)

  • ϵ: Ellipticity parameter (0=linear, 1=circular, default: 0.0)

Optional Arguments

  • Sv: Custom Volkov-phase function (default: t->0.0)
source
StrongFieldDynamics.ClebschGordan Method
julia
ClebschGordan(ja, ma, jb, mb, Jab, Mab)

Computes the Clebsch-Gordan coefficients

source
StrongFieldDynamics.F1_integral_jac Method
julia
F1_integral_jac(pulse::Pulse, a_electron::AtomicElectron, p_electron::ContinuumElectron, thetap::Float64, phip::Float64 ; sign=1)

Returns the integration result of type ComplexF64.

source
StrongFieldDynamics.F1_integral_levin_approxfun Method
julia
F1_integral_levin_approxfun(pulse::Pulse, a_electron::AtomicElectron, p_electron::ContinuumElectron, θ::Float64, ϕ::Float64 ; sign=1)

Computes the pulse shape integral of the form: Returns the integration result of type ComplexF64.

source
StrongFieldDynamics.F1_integral_quadgk Method
julia
F1_integral_quadgk(pulse::Pulse, a_electron::AtomicElectron, p_electron::ContinuumElectron, θ::Float64 ; sign=1)
source
StrongFieldDynamics.F2_integral_jac Method
julia
F2_integral_jac(pulse::Pulse, a_electron::AtomicElectron, p_electron::ContinuumElectron, thetap::Float64, phip::Float64)

Returns the integration result of type ComplexF64.

source
StrongFieldDynamics.F2_integral_levin_approxfun Method
julia
F2_integral_levin_approxfun(pulse::Pulse, a_electron::AtomicElectron, p_electron::ContinuumElectron, θ::Float64, ϕ::Float64)

Computes the pulse shape integral of the form: Returns the integration result of type ComplexF64.

source
StrongFieldDynamics.F2_integral_quadgk Method
julia
F2_integral_quadgk(pulse::Pulse, a_electron::AtomicElectron, p_electron::ContinuumElectron, θ::Float64, ϕ::Float64)
source
StrongFieldDynamics.T0 Method
julia
T0(pulse::Pulse, a_electron::AtomicElectron, p_electron::ContinuumElectron, mj::Rational{Int64}, msp::Rational{Int64}, θ::Float64, ϕ::Float64)

Computes the direct scattering amplitude

source
StrongFieldDynamics.T0_uncoupled Method
julia
T0_uncoupled(pulse::Pulse, a_electron::AtomicElectron, p_electron::ContinuumElectron, m::Rational{Int64}, θ::Float64, ϕ::Float64)

Computes the direct scattering amplitude

source
StrongFieldDynamics.Ylm Method
julia
Ylm(l, m, θ, ϕ)

Computes the Spherical hamonics using the Associated Legender Polynomials as Ylm(l, m, θ, ϕ) = Plm(l, m, cos(θ)) * exp(im_m_ϕ)

source
StrongFieldDynamics.bessel_electron Method
julia
bessel_electron(ε::Float64, l::Int64, r::Vector{Float64}) → (Vector{Float64}, Vector{Float64}, Float64)

Computes the free electron continuum wavefunction using spherical Bessel functions.

Arguments

  • ε::Float64: Kinetic energy of the free electron (in atomic units)

  • l::Int64: Orbital angular momentum quantum number

  • r::Vector{Float64}: Radial grid points (in bohr)

Returns

  • r::Vector{Float64}: Input radial grid

  • P::Vector{Float64}: Radial wavefunction P(r) = r * jₗ(pr)

  • δ::Float64: Phase shift (always 0.0 for free electrons)

Physics

For a free electron in the absence of any potential, the radial wavefunction is:

P(r) = r * jₗ(pr)

where jₗ is the spherical Bessel function of order l, and p = √(2ε) is the momentum.

This represents the exact solution to the Schrödinger equation with V(r) = 0.

Example

julia
r_grid = range(0.1, 50.0, 500)
energy = 1.0  # 1 Hartree kinetic energy
l = 1  # p-wave

r, P, phase = bessel_electron(energy, l, r_grid)
source
StrongFieldDynamics.compute_angular_distribution Method
julia
compute_angular_distribution(electron::ContinuumElectron, bound::AtomicElectron,
                            pulse::Pulse; energy::Float64=1.0,
                            θ::Float64=0.0,
                            ϕ_range::Tuple{Float64,Float64}=(0.0, 2π),
                            n_ϕ::Int=100,
                            coupled::Bool=true) -> AngularDistribution

Computes the angular distribution of photoelectrons at fixed energy and theta using SFA.

Arguments

  • a_electron::AtomicElectron: Initial bound state

  • pulse::Pulse: Laser pulse parameters

  • energy::Float64=1.0: Fixed photoelectron energy (a.u.)

  • θ::Real=pi/2: Fixed polar angle (radians)

  • ϕ_range::Tuple{Float64,Float64}=(0.0, 2π): Azimuthal angle range (radians)

  • n_ϕ::Int=200: Number of azimuthal angle points

source
StrongFieldDynamics.compute_atomic_electron Method
julia
compute_atomic_electron(Z::Int64, n::Int64, l::Int64; ip::Float64=0.0) → AtomicElectron

Computes the atomic electron wavefunction for specified quantum numbers.

Arguments

  • Z::Int64: Atomic number (nuclear charge)

  • n::Int64: Principal quantum number

  • l::Int64: Orbital angular momentum quantum number

  • ip::Float64=0.0: Ionization potential in atomic units (Hartree). If 0.0, uses default values.

Returns

  • AtomicElectron: Structure containing the radial wavefunction and quantum numbers

Supported Atoms and Orbitals

  • Lithium (Z=3): 1s, 2s orbitals

  • Neon (Z=10): 1s, 2p orbitals

  • Argon (Z=18): 1s, 3p orbitals

  • Krypton (Z=36): 1s, 4p orbitals

  • Xenon (Z=54): 1s, 5p orbitals

Physics

The function loads pre-computed Hartree-Fock or Dirac-Fock wavefunctions from data files. These represent accurate many-electron atomic wavefunctions that account for electron correlation and relativistic effects where applicable.

Examples

julia
# Lithium 2s electron
li_2s = compute_atomic_electron(3, 2, 0)

# Neon 2p electron with custom ionization potential
ne_2p = compute_atomic_electron(10, 2, 1; ip=0.8)  # 0.8 Hartree

# Argon 3p electron
ar_3p = compute_atomic_electron(18, 3, 1)

Notes

  • Ionization potentials are converted from eV to atomic units (÷27.21138)

  • Data files contain radial positions and corresponding wavefunction values

  • Interpolation is performed to ensure smooth wavefunction representation

source
StrongFieldDynamics.compute_energy_distribution Method
julia
compute_energy_distribution(a_electron::AtomicElectron, pulse::Pulse; θ::Real=pi/2, ϕ::Real=0.0,
                           energy_range::Tuple{Float64,Float64}=(1e-6, 0.5), n_points::Int64=200, coupled::Bool=true) -> EnergyDistribution

Computes the energy-resolved photoelectron spectrum at specified angles using SFA.

Arguments

  • a_electron::AtomicElectron: Initial bound state

  • pulse::Pulse: Laser pulse parameters

  • θ::Real=π/2: Polar detection angle (radians)

  • ϕ::Real=0.0: Azimuthal detection angle (radians)

  • energy_range::Tuple{Float64,Float64}=(0.0, 0.5): Energy range (a.u.)

  • n_points::Int64=200: Number of energy points

source
StrongFieldDynamics.compute_momentum_distribution Method
julia
compute_momentum_distribution(electron::ContinuumElectron, bound::AtomicElectron,
                             pulse::Pulse; p_max::Float64=2.0,
                             n_p::Int=50, n_θ::Int=25, n_φ::Int=50) -> MomentumDistribution

Computes the 3D momentum distribution of photoelectrons using SFA in spherical coordinates.

Arguments

  • a_electron::AtomicElectron: Initial bound state

  • pulse::Pulse: Laser pulse parameters

  • energy_range::Tuple{Float64,Float64}=(0.0, 0.5): Energy range (a.u.)

  • n_p::Int=50: Number of momentum magnitude points

  • n_θ::Int=1: Number of polar angle points

  • n_φ::Int=50: Number of azimuthal angle points

source
StrongFieldDynamics.compute_partial_wave Method
julia
compute_partial_wave(l::Int64, j::Rational{Int64}, p_electron::ContinuumElectron, a_electron::AtomicElectron) → PartialWave

Computes a partial wave for photoionization calculations.

Arguments

  • l::Int64: Orbital angular momentum quantum number

  • j::Rational{Int64}: Total angular momentum quantum number (j = l ± 1/2)

  • p_electron::ContinuumElectron: Continuum electron specification (energy, solution type)

  • a_electron::AtomicElectron: Bound atomic electron data

Returns

  • PartialWave: Structure containing energy, quantum numbers, wavefunction, and phase shift

Physics

Computes the continuum electron partial wave corresponding to the photoionization transition from the bound atomic orbital. The choice of solution method (Bessel vs Distorted) affects the accuracy:

  • Bessel: Free electron approximation, ignores atomic potential

  • Distorted: Includes scattering from atomic potential (more accurate)

Example

julia
# Set up bound and continuum electrons
atom = compute_atomic_electron(18, 3, 1)  # Ar 3p
continuum = ContinuumElectron(1.0, :Distorted)  # 1 Hartree, distorted wave

# Compute p-wave (l=1) partial wave
partial_wave = compute_partial_wave(1, 3//2, continuum, atom)

Notes

  • The radial grid from the atomic electron is used for the continuum calculation

  • Phase shifts are important for interference effects in photoionization

  • Distorted waves require the atomic potential (currently needs rV to be defined)

source
StrongFieldDynamics.distorted_electron Method
julia
distorted_electron(ε::Float64, l::Int64, r::Vector{Float64}, rV::Vector{Float64}) → (Vector{Float64}, Vector{Float64}, Float64)

Computes the distorted wave continuum electron wavefunction in an atomic potential.

Arguments

  • ε::Float64: Kinetic energy of the electron (in atomic units)

  • l::Int64: Orbital angular momentum quantum number

  • r::Vector{Float64}: Radial grid points (in bohr)

  • rV::Vector{Float64}: Potential energy array V(r) at grid points (in atomic units)

Returns

  • r::Vector{Float64}: Input radial grid

  • P::Vector{Float64}: Normalized radial wavefunction P(r)

  • δ::Float64: Total phase shift (inner + Coulomb contributions)

Physics

Solves the radial Schrödinger equation with the given potential:

d²P/dr² + [2ε - 2V(r) - l(l+1)/r²]P = 0

The wavefunction asymptotically behaves as:

P(r) → sin(pr - lπ/2 + δₗ) / p    as r → ∞

This accounts for scattering from the atomic potential, providing more accurate wavefunctions for photoionization calculations than free electron approximations.

Implementation

Uses a Fortran library (mod_sfree.so) for numerical integration of the radial Schrödinger equation with proper boundary conditions.

Example

julia
r_grid = range(0.1, 50.0, 500)
V_coulomb = -2.0 ./ r_grid  # Hydrogen-like potential
energy = 0.5  # 0.5 Hartree kinetic energy

r, P, delta = distorted_electron(energy, 0, r_grid, V_coulomb)
println("s-wave phase shift: ", delta, " radians")
source
StrongFieldDynamics.inner_product Method
julia
inner_product(p_partialwave::PartialWave, a_electron::AtomicElectron, r::Vector{Float64})

Calculates the overlap integral between continuum and bound electron radial wavefunctions.

This function computes the inner product ⟨ψ_continuum|ψ_bound⟩ of the radial parts of the electron wavefunctions, which appears in the third term of the scattering amplitude calculation and represents the direct overlap between initial and final states.

Arguments

  • p_partialwave::PartialWave: Continuum electron partial wave containing radial wavefunction P

  • a_electron::AtomicElectron: Atomic electron containing bound state radial wavefunction P

  • r::Vector{Float64}: Radial grid points used for numerical integration

Returns

  • ComplexF64: The overlap integral multiplied by the phase factor (-i)^l

Implementation Details

  • Uses cubic spline interpolation (Dierckx.jl) for smooth integration

  • Applies numerical quadrature (QuadGK.jl) for accurate integral evaluation

  • Includes the angular momentum phase factor (-i)^l_p for proper quantum mechanical normalization

source
StrongFieldDynamics.levin_integrate_approxfun Method
julia
levin_integrate_approxfun(f::Function, gp::Function, a::Float64, b::Float64, ga::Float64, gb::Float64)

Compute the highly oscillatory integral ∫[a,b] f(x) exp(ig(x)) dx using Levin's method with ApproxFun.

This implementation uses ApproxFun's automatic differentiation and adaptive function approximation to solve the differential equation ψ'(x) + ig'(x)ψ(x) = f(x) with ψ(a) = 0.

Arguments:

  • f: amplitude function

  • gp: derivative of phase function g'(x)

  • a, b: integration bounds

  • ga, gb: g(a) and g(b) at the integration limits

Returns:

  • ComplexF64: Levin integral approximation ψ(b)e^{ig(b)} - ψ(a)e^
source
StrongFieldDynamics.plot_angular_distribution Method
julia
plot_angular_distribution(ad::AngularDistribution; title::String="", save_path::String="", 
                         plot_type::Symbol=:polar)

Plot the angular distribution of photoelectrons at fixed energy and theta.

Arguments

  • ad::AngularDistribution: Angular distribution data to plot

  • title::String="": Plot title (auto-generated if empty)

  • save_path::String="": Path to save the plot (optional)

  • plot_type::Symbol=:polar: Plot type (:polar for polar plot, :cartesian for line plot)

Returns

  • Figure: CairoMakie figure object
source
StrongFieldDynamics.plot_energy_distribution Method
julia
plot_energy_distribution(ed::EnergyDistribution; title::String="", xlabel::String="Energy (a.u.)", 
                        ylabel::String="Probability", save_path::String="")

Plot the energy distribution of photoelectrons with proper labels.

Arguments

  • ed::EnergyDistribution: Energy distribution data to plot

  • title::String="": Plot title (auto-generated if empty)

  • xlabel::String="Energy (a.u.)": X-axis label

  • ylabel::String="Probability": Y-axis label

  • save_path::String="": Path to save the plot (optional)

Returns

  • Figure: CairoMakie figure object
source
StrongFieldDynamics.plot_momentum_distribution Method
julia
plot_momentum_distribution(md::MomentumDistribution; title::String="", save_path::String="", 
                           slice_type::Symbol=:p_slice, slice_value::Float64=1.0, 
                           colormap=nothing)

Plot the momentum distribution of photoelectrons with different visualization options.

Arguments

  • md::MomentumDistribution: Momentum distribution data to plot

  • title::String="": Plot title (auto-generated if empty)

  • save_path::String="": Path to save the plot (optional)

  • slice_type::Symbol=:p_slice: Type of slice (:p_slice for fixed p, :theta_slice for fixed θ, :phi_slice for fixed φ)

  • slice_value::Float64=1.0: Value at which to take the slice

  • colormap=:viridis: Colormap for heatmaps

Returns

  • Figure: CairoMakie figure object
source
StrongFieldDynamics.probability Method
julia
probality(pulse::Pulse, a_electron::AtomicElectron, p_electron::ContinuumElectron)

Calculates the ionization probability

source
StrongFieldDynamics.probability_uncoupled Method
julia
probality_uncoupled(pulse::Pulse, a_electron::AtomicElectron, p_electron::ContinuumElectron)

Calculates the ionization probability

source
StrongFieldDynamics.reduced_matrix_element Method
julia
reduced_matrix_element(p_partialwave::PartialWave, a_electron::AtomicElectron, r::Vector{Float64})

Computes the reduced matrix element <εp lp jp || p || n l j> for the momentum operator between a continuum electron partial wave and an atomic (bound) electron state.

Arguments

  • p_partialwave::PartialWave: Continuum electron partial wave with quantum numbers lp, jp

  • a_electron::AtomicElectron: Atomic (bound) electron with quantum numbers l, j

  • r::Vector{Float64}: Radial grid points for numerical integration

Returns

  • ComplexF64: The reduced matrix element value

References

  • Equation (16) from PRA 2023 paper for the radial integral formulation
source
StrongFieldDynamics.reduced_matrix_element_uncoupled Method
julia
reduced_matrix_element_uncoupled(p_partialwave::PartialWave, a_electron::AtomicElectron, r::Vector{Float64})

Computes the reduced matrix element <εp lp jp || p || n l j> for the momentum operator between a continuum electron partial wave and an atomic (bound) electron state.

Arguments

  • p_partialwave::PartialWave: Continuum electron partial wave with quantum numbers lp, jp

  • a_electron::AtomicElectron: Atomic (bound) electron with quantum numbers l, j

  • r::Vector{Float64}: Radial grid points for numerical integration

Returns

  • ComplexF64: The reduced matrix element value

References

  • Equation (16) from PRA 2023 paper for the radial integral formulation
source
StrongFieldDynamics.sin2Sv Method
julia
sin2Sv(t::Float64, θ::Float64, pulse::Pulse, p_electron::ContinuumElectron)

Defines the Volkov phase integral over time for a sin² envelope. Only for circular polarization from the Birger PRA 2020.

source
StrongFieldDynamics.sin2Sv_danish Method
julia
sin2Sv_danish(t::Float64, theta::Float64, phi::Float64, pulse::Pulse, p_electron::ContinuumElectron)
source
StrongFieldDynamics.sin2Sv_general Method
julia
sin2Sv_general(t::Float64, θ::Float64, ϕ::Float64, pulse::Pulse, p_electron::ContinuumElectron)

Computes the Volkov phase for a sin² pulse envelope at time t.

source
StrongFieldDynamics.sin2Sv_prime Method

sin2Sv_prime(t::Float64, θ::Float64, ϕ::Float64, pulse::Pulse, p_electron::ContinuumElectron)

source
StrongFieldDynamics.sin2Sv_prime_b Method
julia
sin2Sv_prime_b(t::Float64, θ::Float64, ϕ::Float64, pulse::Pulse, p_electron::ContinuumElectron)

Only for circular polarizaton, Takes the analytical derivative from the expression in Appendix B, Birger 2020 for sin2 pulse

source
StrongFieldDynamics.sin2Sv_prime_general Method
julia
sin2Sv_prime_general(t::Float64, θ::Float64, ϕ::Float64, pulse::Pulse, p_electron::ContinuumElectron)

Computes the time derivative of the Volkov phase for a sin² pulse envelope at time t.

source
StrongFieldDynamics.sin2Sv_quadgk Method
julia
sin2Sv_quadgk(t::Float64, θ::Float64, ϕ::Float64, pulse::Pulse, p_electron::ContinuumElectron)
source