sigway.integrators¶

The integrator performs the double integral over the two internal momenta \(s\in[0,1]\) and \(t\in[0,\infty)\).

Base class¶

sigway.integrators.Integrator ¶

Abstract base class for SIGW double-integral strategies.

An integrator performs the numerical double integral over the two rescaled internal momenta \(s \in [0, 1]\) and \(t \in [0, \infty)\) that defines the scalar-induced gravitational-wave tensor power spectrum. Given a Kernel (which encodes the radiation-transfer physics) and a primordial scalar power spectrum ScalarPerturbations, it returns the un-normalised tensor power spectrum on a wavenumber grid.

The global normalisation prefactor (kernel.norm), the conversion between frequency and wavenumber, and any output upsampling are applied by OmegaGW after calling this method.

Subclasses must override integrate.

Methods:

Name	Description
`integrate`	Evaluate the un-normalised tensor power spectrum on `kvec`. Must be implemented by every concrete subclass.

Simpson quadrature¶

sigway.integrators.SimpsonIntegrator ¶

SimpsonIntegrator(s, t)

Bases: Integrator

2-D Simpson quadrature over the \((s, t)\) momentum plane (Gaussian spectrum).

This integrator performs the double integral that enters the scalar-induced gravitational-wave energy density. It assumes a Gaussian primordial curvature perturbation, so the integrand contains the product of two copies of the dimensionless scalar power spectrum \(\mathcal{P}_\zeta\):

\[ \overline{I^2}(k) = \int_0^1 ds \int_0^\infty dt \; \overline{I^2}(t, s, k) \cdot J(t, s) \cdot \mathcal{P}_\zeta(ku) \cdot \mathcal{P}_\zeta(kv), \]

where the dimensionless internal momenta are defined by

\[ u = \frac{t + s + 1}{2}, \qquad v = \frac{t - s + 1}{2}, \]

and \(J(t, s)\) is the Jacobian of the \((u, v) \to (s, t)\) change of variables (see polynomial).

The integral is discretised on user-supplied grids for \(s\) and \(t\) and evaluated with composite Simpson quadrature. Grid choice matters physically: the integrand peaks near \(t \approx 1\), so the \(t\)-grid should be linear below 1 (to resolve the peak) and geometric above 1 (to cover the slowly-decaying tail efficiently). A few hundred points total is typically sufficient.

If the Kernel declares a resonant feature at a specific value \(t_\mathrm{res}\) (kernel.resonant_t is set), an additional 1-D Simpson integral over \(s\) at that fixed \(t\) slice is computed and added to the 2-D result.

Parameters:

Name	Type	Description	Default
`s`	`array_like or callable`	Quadrature nodes for the momentum variable \(s \in [0, 1]\). May be a fixed 1-D array or a callable `s(kvec, theta)` returning a 1-D array, where `theta` collects all physical parameters `(theta_pz, *theta_k)`. A callable is useful when the required \(s\)-resolution depends on the spectral peak position.	required
`t`	`array_like or callable`	Quadrature nodes for the momentum variable \(t \in [0, \infty)\). May be a fixed 1-D array, a 2-D array of shape `(nt, nk)` providing a separate node sequence per wavenumber, or a callable `t(kvec, *theta)` returning either shape. When a 1-D array is provided it is broadcast internally to `(nt, nk)` so that the non-uniform Simpson rule over \(t\) applies consistently at every \(k\). A \(k\)-dependent grid is useful for tracking a spectral feature whose relevant \(t\)-range shifts with \(k\).	required

Attributes:

Name	Type	Description
`s`	`array_like or callable`	The \(s\)-grid as supplied at construction time.
`t`	`array_like or callable`	The \(t\)-grid as supplied at construction time.

Methods:

Name	Description
`integrate`	Evaluate the un-normalised tensor power spectrum on `kvec`.

Examples:

Build a SimpsonIntegrator with a hybrid \(t\)-grid that is linear near the integrand peak (\(t \lesssim 1\)) and geometric in the tail (\(t \gtrsim 1\)):

>>> import jax.numpy as jnp
>>> from sigway.integrators import SimpsonIntegrator
>>> s = jnp.linspace(0.0, 1.0, 10)
>>> t = jnp.concatenate([jnp.linspace(1e-5, 0.999, 200),
...                      jnp.geomspace(1.0, 1e3, 800)])
>>> integ = SimpsonIntegrator(s, t)   # or pass s=, t= straight to OmegaGW

integrate ¶

integrate(kernel, pzeta, kvec, theta_pz, theta_k)

Evaluate the un-normalised tensor power spectrum \(\overline{I^2}(k)\).

Resolves the \(s\) and \(t\) grids (evaluating them as callables if needed), then performs the composite Simpson quadrature over both momentum variables and returns the result on the requested wavenumber grid.

For analytic primordial spectra (pzeta.jittable = True) the quadrature is compiled at first call and reused for subsequent calls with different physical parameters at the same array shapes, making parameter scans fast. Spectra computed by the Mukhanov-Sasaki solver run the same quadrature without compilation.

Parameters:

Name	Type	Description	Default
`kernel`	`Kernel`	A Kernel instance providing the radiation-transfer function \(\overline{I^2}(t, s, k)\). If the kernel has a resonant feature it additionally provides `overline_Isq_resonant` and `resonant_t`.	required
`pzeta`	`ScalarPerturbations`	The primordial scalar power spectrum \(\mathcal{P}_\zeta(k)\). Its `prepare` method returns a callable that evaluates \(\mathcal{P}_\zeta\) on an arbitrary wavenumber array.	required
`kvec`	`(array_like, shape(nk))`	Wavenumber grid at which \(\overline{I^2}(k)\) is evaluated. Units must be consistent with those used by `kernel` and `pzeta`.	required
`theta_pz`	`tuple`	Physical parameter values for the scalar power spectrum (e.g. amplitude, tilt, peak scale).	required
`theta_k`	`tuple`	Physical parameter values for the kernel (e.g. equation-of-state parameter \(w\)).	required

Returns:

Type	Description
`(Array, shape(nk))`	Un-normalised tensor power spectrum values \(\overline{I^2}(k)\) at each wavenumber in `kvec`. OmegaGW multiplies by `kernel.norm` and the standard prefactors to produce \(\Omega_\mathrm{GW}(f)\).