excursion_set.h#

#include <sif/model/excursion_set.h>
SIF_MODEL_EXCURSION_SET_H#

The excursion-set first-crossing problem: how often a random walk in density first drops below a barrier.

Two paths to the same multiplicity function. sif_ep_multiplicity_function() runs correlated walks by Monte Carlo, which needs the full packed covariance and costs quadratic work in the radius count. sif_ep_multiplicity_function_emu() evaluates a trained network over the same features, needs only the diagonal, and is linear – about a million times faster, at a fitted accuracy quoted with the emulator.

The upcrossing form – Verza et al. (2024) eq. (3.15), the Musso-Sheth rate with the exact scale-dependent <(d delta / dS)^2> – is the analytic limit both are measured against.

sif_real *sif_svdw_multiplicity_function(const sif_real *sigma, uint32_t n, sif_real delta_v, sif_real delta_c)#

Excursion-set void multiplicity function f_ln(sigma), for a constant barrier and uncorrelated steps.

Sheth & van de Weygaert (2004) as given by Jennings, Li & Hu (2013) eq. (8).

Note

The mode series is summed to convergence rather than truncated at the four terms of the reference, which keeps it valid as D = |delta_v| / (delta_c + |delta_v|) approaches 1.

Parameters:

  • sigma – R.m.s. density contrast per entry, strictly positive

  • n – Length of sigma

  • delta_v – Linear void barrier, strictly negative

  • delta_c – Collapse barrier, strictly positive

Returns:

Newly allocated array of n values, released with sif_free_aligned, or NULL on invalid input.

sif_real *sif_ep_barrier_smt(const sif_real *sigma, uint32_t n, sif_real alpha, sif_real beta, sif_real gamma)#

The Sheth-Mo-Tormen moving barrier, B(sigma) = alpha [1 + (beta / sigma)^gamma].

gamma = 0 gives the constant 2 alpha; large sigma tends to alpha.

Note

Positive, and compared as an upcrossing by the first-crossing routines below. For a void that means alpha carries the magnitude of the underdensity, not its sign, and the walk is that of -delta.

Parameters:

  • sigma – R.m.s. density contrast per smoothing scale, strictly positive. Take it from sif_delta_covariance_pk, so the barrier and the walk share one variance.

  • n – Length of sigma

  • alpha – Barrier amplitude, strictly positive

  • beta – Barrier scale in units of sigma, strictly positive

  • gamma – Barrier slope

Returns:

Newly allocated array of n barrier heights, released with sif_free_aligned, or NULL on invalid input.

uint64_t *sif_ep_first_crossing_counts(const sif_real *radii, uint32_t n_radii, const double *cov, const sif_real *barrier, uint64_t n_paths, uint64_t seed, sif_option opt)#

Raw first-crossing counts of a correlated random walk against a moving barrier, one count per smoothing radius.

Each path draws n_radii standard normals xi and forms the walk delta(R_j) = sum_{m <= j} L_jm xi_m, with L L^T = cov, stepping from the largest radius down and recording the first radius at which delta reaches the barrier. Paths that never reach it are not recorded, so the counts sum to at most n_paths.

Note

The barrier is an upcrossing condition, delta >= B, which makes the sign convention of the barrier the caller’s responsibility. Reversing it produces a plausible-looking curve rather than an obvious failure.

Note

The result is a deterministic function of (seed, n_paths) alone: the counts do not depend on the thread count or on how the loop was scheduled.

Parameters:

  • radii – Smoothing radii, strictly positive and strictly increasing, at least 2 and at most SIF_COV_MAX_RADII

  • n_radii – Number of radii

  • cov – Covariance of the smoothed field between every pair of radii, packed lower triangle in the ascending order of radii: S(i, j) at SIF_COV_INDEX(i, j) for j <= i. Must be symmetric positive semi-definite with a strictly positive diagonal. Double rather than sif_real: a single-precision factorization of a realistic radius grid reaches a non-positive pivot and fails.

  • barrier – n_radii barrier heights, in the ascending order of radii

  • n_paths – Number of walks, strictly positive

  • seed – Seed for the walk ensemble

  • opt – Reserved; pass SIF_DEFAULT

Returns:

Newly allocated array of n_radii counts in the ascending order of radii, released with sif_free_aligned, or NULL on invalid input or a failed factorization.

sif_real *sif_ep_multiplicity_function(const sif_real *radii, uint32_t n_radii, const double *cov, const sif_real *barrier, uint64_t n_paths, uint64_t seed, uint64_t *counts, sif_option opt)#

Lagrangian void multiplicity function from the first crossing of a moving barrier by a correlated random walk.

f[i] = (walks first crossing at radii[i])

/ (n_paths * (radii[i+1] - radii[i]))

Note

The result has n_radii - 1 entries, on the bin centres 0.5 (radii[i] + radii[i+1]), not n_radii point values on the radii themselves; the largest radius is the walk’s first step and has no bin above it. This is a different shape from sif_svdw_multiplicity_function, which is point-evaluated.

Note

The Monte Carlo error on a bin is sqrt(n_i) / (n_paths dr_i), so a bin holding a fraction p of the walks is known to a relative 1 / sqrt(p n_paths). Take counts and check it before trusting the tails; re-running the walk through sif_ep_first_crossing_counts to get the same numbers doubles the cost of the most expensive call in the library.

Parameters:

  • counts – Optional output, n_radii entries: the raw first-crossing counts this multiplicity was built from, exactly as sif_ep_first_crossing_counts would have returned them. Pass NULL to skip.

Remaining parameters, and the barrier and reproducibility conventions, are as sif_ep_first_crossing_counts.

Returns:

Newly allocated array of n_radii - 1 values, released with sif_free_aligned, or NULL on invalid input or a failed factorization.

struct sif_emu_domain_t#

Report on where an emulated call sat relative to what the emulator was trained over. Optional; pass NULL to sif_ep_multiplicity_function_emu if the answer is all that is wanted.

Leaving the trained region is not an error and does not fail the call. It costs accuracy – around 0.33% against the usual 0.13% at the edge that was measured – and a sampler should be able to notice that without parsing the log or aborting a chain over one proposal.

sif_real *sif_ep_multiplicity_function_emu(const sif_real *radii, uint32_t n_radii, const sif_real *sigma, const sif_real *barrier, const double *deriv_variance, sif_emu_domain_t *domain, sif_option opt)#

The same multiplicity function as sif_ep_multiplicity_function, emulated: no random walks, no paths, well under a millisecond.

The up-crossing rate of Verza et al. (2024) eq. (3.15) corrected by a small trained network. The correction multiplies a hazard rather than the multiplicity itself, and the result is rebuilt through the survival recursion, so it is non-negative and integrates to at most one whatever the network predicts.

Note

Takes sigma, not the packed covariance that sif_ep_multiplicity_function needs. That is deliberate rather than an oversight: the emulator reads only the diagonal, so it needs n numbers where the Monte Carlo needs n(n+1)/2, and its cost is linear rather than quadratic in the radius count.

Note

Returns the answer the Monte Carlo CONVERGES to, not the answer it gives on the caller’s grid. Those differ: a first-crossing walk sampled at 50 radii sits about 1% from its own continuum limit. Comparing the two on a coarse grid shows a disagreement larger than either method’s error, and the emulator is the one to trust. The emulated result moves by under 0.25% between 64 and 256 radii.

Note

Trained over CDM-like spectra and Sheth-Mo-Tormen barriers. Other spectral families are outside the contract, and the domain report does not reliably catch them: a per-feature range check passes power-law spectra that are several per cent wrong, because what distinguishes them is the shape of the whole trajectory rather than any pointwise value.

Parameters:

  • radii – Smoothing radii, strictly positive and strictly increasing, at least 3. Extend the grid outward far enough that few walks start above the barrier; domain reports whether that was achieved.

  • n_radii – Number of radii

  • sigma – R.m.s. density contrast per radius, strictly positive. Take it from sif_delta_covariance_pk so the barrier and the walk share one variance.

  • barrier – n_radii barrier heights, in the ascending order of radii, strictly positive – stricter than the Monte Carlo, which takes any finite value. The features here are built from log(B / sigma), and the correction was trained over nu from 0.23 upward, so the convention is a positive barrier crossed from below. That is what sif_ep_barrier_smt() produces; a void threshold enters through its magnitude, as it does in the BBKS routines.

  • deriv_variance – n_radii entries: <(d delta / dS)^2> from sif_delta_covariance_pk. REQUIRED, unlike the Monte Carlo path, which has no use for it: differencing it off a covariance converges only at first order, which would make the answer depend on how finely radii was sampled – the one property this entry point exists to avoid.

  • domain – Optional report; pass NULL to skip

  • opt – Reserved; pass SIF_DEFAULT

Returns:

Newly allocated array of n_radii - 1 values on the same bin centres as sif_ep_multiplicity_function, released with sif_free_aligned, or NULL on invalid input. An input outside the trained region is NOT invalid: the answer is returned, domain records it, and the log names the quantity responsible.