Input pulses#

Every propagation starts from an initial spectrum. solprop.pulses builds one from a physical description — a duration, a centre wavelength, an energy, a spectral phase — and hands back the array the solver expects. The keyword drivers (prop_gnlse, prop_capillary, prop_waveguide) build a Gaussian or sech pulse from fwhm, energy or peak_power and phases themselves; the builders on this page are what they use, and what to reach for when the input is anything else — a super-Gaussian, a measured spectrum, a continuous wave, several pulses at once, or a transverse beam.

This page covers the six temporal builders, the two ways of fixing a pulse’s size, the sign and phase-origin conventions that decide where in the window a pulse lands and which way it chirps, measured and continuous-wave inputs, multi-colour inputs, and the transverse beams of the spatially resolved geometries. The normative statements are in conventions, §2 and §4; this page derives them and shows them at work.

One method, two grids#

A pulse in solprop is a specification, not an array. Each builder is a frozen dataclass with a single method:

from solprop.grid import EnvelopeGrid, FieldGrid
from solprop.pulses import GaussianPulse

pulse = GaussianPulse(wavelength=800e-9, fwhm=30e-15, energy=1e-6)

field = FieldGrid.from_limits(1.0, 800e-9, (400e-9, 2e-6), 2e-12)
envelope = EnvelopeGrid.from_limits(1.0, 800e-9, (400e-9, 2e-6), 2e-12)

spectrum_field = pulse.build(field)        # (1, 1, n_omega), real field's spectrum
spectrum_envelope = pulse.build(envelope)  # (1, 1, n_omega), envelope's spectrum

The same specification builds against either grid and gives the same physical pulse. On a field grid the time samples are \(\sqrt{I(t)}\cos(\omega_p t)\), whose analytic signal is the pulse; on an envelope grid they are \(\sqrt{I(t)}\,e^{i(\omega_p-\omega_0)t}\), which is that analytic signal with the grid’s carrier removed. Because conventions, §4, gives the two grids one energy functional, the spectral energy density that comes out is the same function of \(\omega\) either way, and the test suite checks that it is.

The returned array is always (n_modes, n_pol, n_omega) with both leading axes equal to one, which is the canonical kernel layout of conventions, §6. Nothing squeezes it; a multimode input is a stack of these.

The six builders are GaussianPulse, SechPulse, CustomPulse, DataPulse, TimeDataPulse and CWPulse. They all satisfy the Pulse protocol, which is what downstream code types against. The first three share ShapedPulse, which holds the parameters they have in common and the whole of build; subclass it and supply an intensity(t) to add a shape of your own. A shape function that returns the wrong length, a negative value, a non-finite value or nothing but zeros is rejected with a message that names the problem, so a bad envelope never reaches the initial state as a NaN.

builder	intensity \(I(t)\)	duration given as
`GaussianPulse`	\(\exp\!\left[-\tfrac12 (t/\sigma)^{2m}\right]\)	`fwhm` (power FWHM), order `m`
`SechPulse`	\(\operatorname{sech}^2(t/\tau_w)\)	`fwhm` or the natural width `width`
`CustomPulse`	any callable `I(t)`	implicit in the callable
`DataPulse`	measured spectrum, interpolated	implicit in the data
`TimeDataPulse`	measured field, transformed exactly	implicit in the data
`CWPulse`	monochromatic or \(\operatorname{sech}^2\) in \(\omega\)	`bandwidth`, a wavelength FWHM

The FWHM conversions are exact for every order. A super-Gaussian of order \(m\) uses \(\sigma = \tau/\left[2(2\ln 2)^{1/2m}\right]\), so that \(I(\tau/2) = \tfrac12\) whatever \(m\) is; a \(\operatorname{sech}^2\) pulse uses \(\tau_w = \tau/[2\ln(1+\sqrt2)]\).

Energy or peak power, never both#

A pulse is scaled by its energy or by its peak_power. Giving both, or neither, is an error rather than a precedence rule. Both are realised exactly on the grid: the spectrum is built, the phase is applied, the realised quantity is measured with the grid’s own functional, and one multiplication by a square root fixes it. Nothing is left to a quadrature error, and the round trip closes to a part in \(10^{12}\) or better.

The order matters: the peak power is the peak power after the spectral phase. For a transform-limited input, broadening it by chirp while holding its peak power fixed therefore raises its energy. Holding the energy fixed is a different experimental constraint. Scaling by energy is insensitive to the order, because a spectral phase moves energy about in time but does not create or destroy it.

Warning

A pulse whose spectrum falls outside the grid’s window, or whose duration is longer than the time window, is not silently truncated into something plausible: it is scaled so that whatever the grid does hold carries the requested energy. Check the duration against trange and the bandwidth against the wavelength limits before reading a result.

Signs: where a pulse lands and which way it chirps#

A phase convention determines the sign of a delay or chirp and the time origin of a measured spectrum. Use the same convention when specifying the input and interpreting the output.

The interface is the optics convention. A spectral phase given to a builder (its phase field, the phases keyword of the drivers) is \(\tilde E(\omega) = |\tilde E|\,e^{i\phi(\omega)}\) with

\[ \phi(\omega) = \omega t_0 + \sum_{n} \frac{\phi_n}{n!}\,(\omega - \omega_p)^n, \]

\(\omega_p\) being the pulse’s own centre frequency, not necessarily the grid’s reference. The group delay is \(T_g(\omega) = \mathrm{d}\phi/\mathrm{d}\omega\), so

a positive t0 makes the pulse arrive later (it must lie strictly inside the window: a delay is a spectral phase, and one past the half-window would wrap the pulse round to the other end, so it is refused);
phase[1], the first-order coefficient, is a group delay and does exactly what t0 does;
a positive phase[2] (GDD) gives \(T_g = \phi_2(\omega-\omega_p)\), so the blue arrives after the red: an upchirp, the sign of propagation through normally dispersive glass;
phase[0] is the carrier-envelope phase at the pulse’s own peak, on either grid and independent of the window length.

The stored array is the FFT convention. The conversion happens once per build, where the phase factor is formed, through to_internal(). The pulse builders are one of only three places where the optics convention meets the internal one (see the design decisions); nothing else in solprop.pulses conjugates anything.

The phase origin, and why `(−1)^k` exists#

A stored spectrum is the raw forward DFT of the sample array, so its phase origin is the first sample, at \(t = -n_t\,\delta t/2\), not the centre of the window. Referenced to \(t = 0\) the stored bin \(k\) therefore carries an extra \(e^{+i\omega_k n_t \delta t/2}\), and since \(\omega_k = k\,\delta\omega\) with \(\delta\omega = 2\pi/(n_t\delta t)\),

\[ \omega_k \frac{n_t \delta t}{2} = \pi k \qquad\Longrightarrow\qquad e^{+i\omega_k n_t \delta t / 2} = (-1)^k . \]

On an envelope grid the transform’s own variable is the relative frequency, whose bin \(k\) gives \(\pi k\) as well because \(n_t\) is even, so the factor is \((-1)^k\) on the storage index for both layouts. Putting the two halves together,

\[ \tilde E_{\text{optics},k} = (-1)^k\,\overline{\tilde E_k}, \]

which is what spectrum_to_optics() and spectrum_to_internal() do. Both are involutive and each is the other’s inverse.

The practical consequence: a transform-limited pulse sitting at \(t = 0\) reads back a flat optics phase, while its stored phase alternates by \(\pi\) from bin to bin. That alternation is not a bug and not a chirp; it is where the window’s centre is.

A phase increment, such as a Taylor phase or an extra delay, multiplies an existing spectrum. Its phasor needs the sign conversion, while the existing spectrum already carries its time origin. A complete measured spectrum specified about \(t=0\) needs both the conjugation and the origin shift in the spectrum helpers. DataPulse handles that conversion. This distinction concerns whether an array is a complete spectrum or a phase factor; multiplication of spectral phase factors itself is commutative.

The delay uses the absolute \(\omega\) on both grids, so it moves the whole field rather than only its envelope: the carrier-envelope phase slips by \(\omega_0 t_0\) and the field and envelope representations stay in step.

Measured pulses#

DataPulse stores a measured complex spectrum, because that is what a spectrometer and a phase retrieval give. TimeDataPulse stores a measured field in time, because that is what a field measurement gives.

from solprop.pulses import DataPulse, TimeDataPulse

DataPulse(omega=omega, spectrum=complex_spectrum, energy=1e-6)
DataPulse.from_spectral_density(omega, intensity, spectral_phase, energy=1e-6)
TimeDataPulse(t=t, field=measured_field, energy=1e-6)

A DataPulse interpolates the magnitude and the phase separately with a cubic B-spline. Interpolating the complex spectrum itself would be wrong: its real and imaginary parts oscillate at the pulse’s arrival time, and no smooth interpolant follows them. from_spectral_density keeps the phase it is given, in the pulse’s spectral_phase, and that phase is interpolated as it stands. A complex spectrum holds its phase only modulo \(2\pi\), so it is unwrapped first, and unwrapping cannot follow a change of more than \(\pi\) between neighbouring samples: 20 000 fs² of group-delay dispersion on a 1.6 nm axis at 800 nm is built quite wrongly, and with no warning. Give a measured phase through from_spectral_density. Outside the data’s own frequency range the built spectrum is exactly zero, the grid’s spectral collar is applied so that nothing is left sitting on the window edge, and a warning is logged when the grid’s active band reaches past the data.

A spline through a magnitude and a phase cannot follow a spectral fringe to its null, so a DataPulse suits one compact pulse and not several far apart. If the measurement is a field in time, use a TimeDataPulse, which interpolates nothing.

A TimeDataPulse takes out the record’s mean \(\bar{E}\) and evaluates the transform at the grid’s own frequencies,

\[ \tilde{E}(\omega_k) = \sum_n \bigl(E_n - \bar{E}\bigr)\, e^{+i\omega_k t_n}, \]

with the chirp-z transform, whose cost is that of a few FFTs of the record and the grid together. This is exact whatever the record holds: any number of pulses, at any separation, with any fringe contrast, and each lands where it sat wherever the record starts, because \(t_n\) is the record’s own axis. It is what the record says rather than what the pulse was:

The spectrum is zero at frequencies above the record’s Nyquist frequency \(\pi/\Delta t\), and a warning is logged when the grid’s active band reaches past it.
A baseline is an offset of the instrument and not part of a field that propagates, so the mean is taken out: left in, it would be a step at each end of the record, and 1 % of the peak would change a 20 fs Gaussian by 0.15 %. A background that varies across the record stays, and is such a step where it does not vanish at the ends. A pulse with a net area, such as a half-cycle pulse, makes a mean of its own, and taking it out is that step: a unipolar 10 fs Gaussian in an 800 fs record changes by 3.3 % on a 1 to 50 µm band.
A pulse that reaches the ends of the record is cut there.
A record longer than the grid’s window folds back into it, as any spectrum sampled at the grid’s spacing must, and a warning is logged when that takes more than a millionth of the energy.

A real array is the electric field; a complex one is the analytic signal in the optics convention. The run record holds the record itself, which costs the same as the spectrum at its own bins: two real numbers for each real sample.

Measured spectra from a file#

load_spectrum() reads a spectrum and its retrieved phase from a text file and returns a DataPulse. The default layout is the one a spectrometer and a phase retrieval usually write, and the one Luna.jl’s DataField documents: three columns of frequency in Hz, spectral power density in arbitrary units and spectral phase in radians.

# frequency (Hz)   power (a.u.)   phase (rad)
3.4e14             0.0012         4.1
3.5e14             0.0310         2.7
...

from solprop.pulses import load_spectrum

pulse = load_spectrum("retrieved.txt", energy=1e-6)
pulse = load_spectrum(
    "spectrum.csv", energy=1e-6, axis="wavelength", axis_unit="nm", columns=(0, 1)
)
pulse = load_spectrum("for_luna.dat", energy=1e-6, phase_convention="fft")

The file. Columns are separated by commas, or by whitespace on a line with no comma. A # starts a comment, blank lines are skipped, and the first line that is not a comment may be a header of column names. The rows may run either way along the axis: they are sorted. A two-column file has a flat phase; columns picks other columns, as zero-based (axis, power) or (axis, power, phase), and further columns are ignored. path may also be an open text stream such as an io.StringIO. A file that cannot be read as a spectrum is refused with a message that names the problem: no such file, a line that is not a row of numbers, too few columns or rows, a value that is not finite, an axis value that is not positive or repeats, or a power column with no positive value. Negative powers, which a background subtraction leaves where there is no signal, are read as zero; one more negative than a tenth of the peak is refused as the wrong column.

The axis. axis_unit defaults to the axis’s own SI unit, Hz or m. A frequency may be in "Hz", "GHz", "THz" or "PHz", or as an angular frequency in "rad/s" or "rad/fs"; a wavelength in "m", "um" (or "µm") or "nm". A density per unit wavelength is converted to one per unit angular frequency, which is what a DataPulse holds:

\[ P_\omega = P_\lambda \left|\frac{\mathrm{d}\lambda}{\mathrm{d}\omega}\right| = P_\lambda\,\frac{\lambda^2}{2\pi c}, \]

so that each band holds the same energy on either axis. The constant goes with the arbitrary units, since energy sets the scale.

The sign. The phase is read in the optics convention of this page: a delay \(\tau\) is \(+\omega\tau\) and puts the pulse at \(+\tau\). A file written for Luna.jl holds the FFT convention, in which the same delay is \(-\omega\tau\); phase_convention="fft" reads it with the opposite sign, so that the same file puts the pulse at \(-\tau\).

The phase. The pulse is built through DataPulse.from_spectral_density, so the phase is interpolated as read, by a cubic spline through the file’s unwrapped values, and is never re-derived from a complex spectrum. Wrapping a phase and unwrapping it again cannot recover a change of more than \(\pi\) between samples, so a large chirp is best given unwrapped. The column’s own range decides how it is read:

A column whose values all lie within \(2\pi\) (and a milliradian, for the rounding of a file’s last digit) is taken as wrapped and unwrapped; any other is taken as it stands. wrapped=True or wrapped=False says which where that guess would be wrong.
A phase in degrees needs phase_unit="deg": wrapped, it spans 360, and as radians it would be taken as unwrapped.
A wrapped column is unwrapped only over the samples whose power is above floor times the peak of the power column (\(10^{-4}\) by default). Between them the phase is interpolated linearly across any gap, and beyond the first and last it is held at their value. The random phase of the noise floor, or of a sample that drops into it within the band, therefore cannot add a turn to the phase of the pulse. An unwrapped column is taken as it stands, noise floor and all.
A wrapped phase that is undersampled cannot be recovered, and a SolpropWarning says so. The sign of it is a step between neighbouring samples above the floor within a tenth of \(\pi\) of \(\pm\pi\): a chirp whose step per sample grows past \(\pi\) passes through there in every turn, and every sample beyond the innermost such step, outward from the peak, may carry a wrong turn. The warning is given when those samples hold more than a thousandth of the power. A phase that steps by well over \(\pi\) everywhere, such as a delay sampled too coarsely, aliases to a smaller step and cannot be detected at all.

Continuous wave#

CWPulse writes its spectrum straight onto the grid’s own bins, so the field is periodic in the time window by construction and needs no apodisation. bandwidth = 0 puts everything on the single nearest bin; a positive bandwidth, given as a wavelength FWHM and converted with \(\Delta\omega = 2\pi c\, \Delta\lambda/\lambda^2\), gives a \(\operatorname{sech}^2\) power spectrum with a uniform random spectral phase. The conversion is first order at the centre wavelength: the spectrum’s frequency FWHM is \(\Delta\omega\), the frequency equivalent of the requested width there, and its wavelength FWHM equals bandwidth only in the narrowband limit (a wide spectrum is asymmetric in wavelength). The phase is drawn from a generator of its own, never the global one: give seed and the realisation is reproducible (the seed is a field of the dataclass and appears in its description); leave it None and each build draws a fresh phase. A CWPulse builds on either grid.

The average_power is realised in one pass, over the samples where the temporal collar is exactly one: the spectrum is written onto the bins and scaled once so that the mean power over that flat part of the window is the requested one. There is no iteration between a fixed spectral magnitude and a windowed field; the single pass trades a slightly perturbed spectral shape for an exact power and a bounded cost.

Several pulses at once#

Two-colour and multi-colour experiments — a pump and a seed for four-wave mixing, a strong pump and a weak probe for cross-phase modulation, a pulse and a delayed replica — launch the sum of independently specified pulses. PulseSum holds such a set: each member is any pulse above, with its own centre wavelength, duration, energy or peak power, spectral phase and delay t0, and is built and scaled on the grid exactly as it would be alone; the input is the sum of the members’ fields. The energies of members that are separated in time or in frequency therefore add, while members that overlap interfere, as fields do: two identical pulses in phase carry four times the energy of one.

The keyword interfaces (prop_gnlse, prop_capillary, prop_waveguide and their Simulation constructors) take the members as pulses=[...], in place of the single-pulse keywords (fwhm, width, shape, energy, peak_power, phases, which must then be left unset); wavelength still sets the grid’s reference wavelength, and the grid’s wavelength_limits must cover every member:

from solprop.driver import prop_capillary
from solprop.pulses import GaussianPulse

pump = GaussianPulse(wavelength=267e-9, fwhm=35e-15, energy=95e-6)
seed = GaussianPulse(wavelength=800e-9, fwhm=35e-15, energy=45e-6, t0=20e-15)
run = prop_capillary(
    75e-6,
    1.5,
    "He",
    0.48,
    wavelength=267e-9,
    pulses=[pump, seed],
    wavelength_limits=(120e-9, 3e-6),
    trange=400e-15,
)

A positive t0 makes a member arrive later, so here the seed trails the pump by 20 fs. The members are recorded in the run’s configuration, and a Simulation built with Simulation.from_capillary(..., pulses=...) or from any initial=PulseSum(...) behaves the same way. The worked examples on four-wave mixing to the vacuum ultraviolet and on cross-phase modulation in the deep ultraviolet use this.

Launching into chosen modes#

A pulse given by keywords enters the first mode of a multimode run, and the others fill through the nonlinearity. launch= puts it into chosen modes instead, on prop_capillary, prop_waveguide and their Simulation constructors (SEP 41): a mapping from a mode to its share, or to a share and a phase in radians.

import warnings

from solprop.driver import Simulation
from solprop.exceptions import MaterialRangeWarning

with warnings.catch_warnings():
    warnings.simplefilter("ignore", MaterialRangeWarning)
    simulation = Simulation.from_capillary(
        125e-6,
        0.1,
        "He",
        1.0,
        wavelength=800e-9,
        wavelength_limits=(300e-9, 3e-6),
        trange=100e-15,
        fwhm=10e-15,
        energy=5e-6,
        modes=3,
        plasma=False,
        n_r=16,
        n_save=2,
        launch={"HE1,1": 0.9, "HE1,2": {"share": 0.1, "phase": 1.57}},
    )
simulation.configuration["modal_launch"]
# {'split': 'energy', 'pulses': [{'labels': ['HE1,1:even', 'HE1,2:even'],
#   'shares': [0.9, 0.1], 'phases': [0.0, 1.57]}]}

Modes are named as channel= names a general guide’s channels: by display label ("HE1,1:even") or by index. A label may leave out a parity that the set carries once, so "HE1,1" is "HE1,1:even" in a capillary; in a set that carries both polarisations it is refused, naming both. A mode the set does not have is refused, naming those it has.
Shares are energy fractions, normalised to one, so {"HE1,1": 9, "HE1,2": 1} is the launch above, and each mode carries its share of the pulse’s energy.
The phase is in the optics convention (conventions): a mode’s spectrum is the pulse’s times \(e^{i\phi}\), so the phase adds to the pulse’s carrier-envelope phase in that mode. It sets the relative phase of the modes where they overlap, which decides how they beat.
Several pulses. One mapping launches every pulse of pulses= the same way; a list of mappings, one per pulse, launches each into its own modes, so a pump can enter the fundamental and a seed a higher mode: launch=[{"HE1,1": 1.0}, {"HE1,2": 1.0}].
Beside the other launches. launch= generalises channel=, so the two are not given together; nor is launch= given beside incident= and coupling=, which compute the launch from a beam. An array is launched with Simulation.with_input, which also drops the launch record, since the array replaces it.

On a general modal basis. A guide in its provider’s own coordinates is orthonormal in true power, and the shares are energy fractions there too. A guide that declares a change of basis (general modal guides) carries a flux metric \(W\) whose cross terms make the energy \(a^\dagger W a\) of the coefficients depend on their relative phases. There the shares are relative amplitudes squared: the coefficients point where the shares and phases say, and are scaled so that their flux-metric energy is the pulse’s. The per-mode split is then a coordinate one, not an energy split, and the record says so with "split": "coordinate".

The record. The run records each pulse’s labels, normalised shares and phases as configuration["modal_launch"], which a resume compares (resuming), so a file written with one launch is not continued with another. A run with no launch records no entry and reads as before.

Transverse beams#

The builders above describe a pulse in time. For the spatially resolved geometries a temporal pulse is wrapped in a beam: SpatialPulse gives it a (super-)Gaussian transverse profile with a 1/e² intensity radius waist and a super-Gaussian order, a signed launch distance from the waist (propz, negative to launch before it) and a Jones vector (jones, see polarisation), for the free-space radial, Cartesian and line grids; PlanarPulse projects an astigmatic Gaussian beam onto the slab modes of the hybrid planar grid (see spatial decompositions). In both, the temporal pulse’s energy or peak_power is the whole beam’s physical energy or peak axial power through the launch plane, in the medium there — per metre of y (J/m, W/m) on a line grid, which is uniform in y. A beam built directly, SpatialPulse.build(grid, space, medium=..., z=..., side=...), names that medium and plane; there is no vacuum default. A medium of index \(n_0\) therefore carries a field \(\sqrt{n_0}\) lower than vacuum for the same energy (energy and flux in a medium).

A beam launched into a guide’s modes, or the output field of another run used as an input, is a coupling rather than a pulse: the beam keeps its energy over its whole incident plane and the modes receive only what they capture (coupling simulations).

Where the quantum floor sits#

Building a pulse says nothing about the noise it sits on. The shot-noise background is a separate object with its own page, noise; the two are added by the driver, never by a builder. photons_per_bin() is the bridge: it reports how far above one photon per bin a built pulse sits, bin by bin.

Where solprop’s choices differ from Luna.jl’s, the differences are listed in Noteworthy differences from Luna.