Solvers

GNLSESolver

Abstract base type for all GNLSE propagation solvers.

ERK4IP(; rtol=1e-6, atol=1e-8, dz_init=nothing, dz_min=0.0)

Embedded Runge-Kutta 4(3) solver in the Interaction Picture.

Fields

• rtol::Float64: Relative error tolerance for adaptive stepping.
• atol::Float64: Absolute error tolerance for adaptive stepping.
• dz_init::Union{Float64, Nothing}: Optional initial step size [m].
• dz_min::Float64: Minimum allowed step size [m]. If the adaptive controller shrinks below this floor (e.g. because the field diverged to NaN/Inf, or the local error never converges), propagation stops with a clear error instead of looping forever. 0.0 (default) auto-selects max(1e-15, fiber_length * 1e-12).

SSFM(dz)

Fixed-step Symmetric Split-Step Fourier Method (SSFM) solver configuration.

AdaptiveSSFM(; phi_max=1e-3, dz_init=nothing, dz_min=1e-12, dz_max=1.0)

Phase-controlled adaptive Split-Step Fourier Method (SSFM) solver. Controls step size via maximum nonlinear phase shift per step: Δzopt = phimax / (γphys * Pmax).

SimParams

Simulation parameters for GNLSE propagation.

Fields

• medium::Medium: Fiber medium parameters
• z_saves::Int: Number of snapshots to save along fiber
• raman_model::Union{RamanModel,Nothing}: Raman scattering model (Nothing to disable)
• self_steepening::Bool: Enable self-steepening/shock effect
• solver::GNLSESolver: Numerical solver (default: ERK4IP()). Tolerances are set on the solver, e.g. ERK4IP(; rtol=1e-6, atol=1e-8).

Notes

Following gnlse-python's GNLSESetup structure:

• z_saves determines memory allocation for result arrays
• Raman disabled by setting raman_model = nothing or fr = 0
• Self-steepening modifies W frequency grid

Solution{T<:Complex}

Solution to GNLSE propagation.

Fields

• t::Vector{Float64}: Time domain grid [s]
• W::Vector{Float64}: Absolute angular frequency grid [rad/s]
• omega0::Float64: Central angular frequency [rad/s]
• Z::Vector{Float64}: Propagation distances [m]
• At::Matrix{T}: Time domain solution (N × z_saves)
• AW::Matrix{T}: Frequency domain solution (N × z_saves)

Notes

Following gnlse-python's Solution structure:

• At[i, j] = pulse at time t[i], distance Z[j]
• AW[i, j] = frequency component at W[i], distance Z[j]

VectorialSolution{T<:Complex}

Solution to coupled vector GNLSE propagation.

Fields

• t::Vector{Float64}: Time domain grid [s]
• W::Vector{Float64}: Absolute angular frequency grid [rad/s]
• omega0::Float64: Central angular frequency [rad/s]
• Z::Vector{Float64}: Propagation distances [m]
• At::Array{T, 3}: Time domain solution array of size N × 2 × z_saves
• AW::Array{T, 3}: Frequency domain solution array of size N × 2 × z_saves

solve(pulse::AbstractPulse, params::SimParams; progress::Bool=true, rng::AbstractRNG=default_rng())

Solve GNLSE following gnlse-python conventions using adaptive ERK4IP method.

rng seeds the Amplified Spontaneous Emission (ASE) noise for AmplifyingMedium propagation (ignored for all other media).

solve(pulse::AbstractPulse, stages::AbstractVector; progress::Bool=true)

Propagate an optical pulse through a sequence of stages (fibers, amplifiers, filters, etc.).

Returns a vector of results corresponding to each stage (either a Solution for propagation stages, or a new Pulse for lumped elements / functional stages).

solve(pulse::VectorialPulse, params::SimParams; progress=true)

Propagate a VectorialPulse through a BirefringentMedium according to SimParams.

solve_sweep(setup_fn::Function, param_space::AbstractArray; progress::Bool=true)
solve_sweep(pulses, params_list; progress::Bool=true)

Perform a multi-threaded parallel parameter sweep over GNLSE simulations.

Concurrently executes simulations across available Julia worker threads (Threads.nthreads()), ensuring thread-safe, isolated workspace memory per trial.

Arguments

• setup_fn::Function: Callable param -> (pulse, params) that maps a parameter combination to its pulse and simulation parameters.
• param_space::AbstractArray: Collection or grid of parameters.
• progress::Bool: Display a progress bar tracking trial completions (default: true).

Returns

• Array of Solution (or VectorialSolution) objects matching the shape of param_space.

Example

# Sweep peak power P0 from 1 kW to 10 kW across 10 trials
powers = range(1e3, 10e3; length=10)
sols = solve_sweep(powers) do P0
    p = sech_pulse(grid, P0, 100e-15)
    params = SimParams(; medium=medium, z_saves=20)
    return (p, params)
end