stochkin.mfpt

Mean first-passage time and transition-rate estimation.

stochkin.mfpt

Mean first-passage time (MFPT) and transition-rate estimation.

This module provides tools for computing MFPTs between metastable basins on a free-energy landscape:

Two-basin (pairwise) MFPT

• compute_mfpt() / compute_mfpt_1d() – trajectory-shooting MFPT between two basins (Langevin or overdamped BD).

• compute_bidirectional_mfpt() / compute_bidirectional_mfpt_1d() – convenience wrappers that compute A→B and B→A in one call.

Multi-basin MFPT network

• compute_mfpt_network() – first-exit overdamped BD simulations across all basins of a BasinNetwork.

• compute_mfpt_network_fpe() – grid-based backward Kolmogorov solve for all-pairs MFPTs (no trajectories needed).

CTMC rate matrix estimation

• estimate_rate_matrix() – build a continuous-time Markov chain generator \(K\) from first-exit statistics or inverse-MFPT fallback.

Diagnostics

• estimate_transition_rates() – simple two-state rate and equilibrium-population summary.

• plot_mfpt_statistics() – histogram of passage-time distributions.

stochkin.mfpt.single_passage_time(args)

Run a single trajectory and return the first-passage time.

This is a top-level, picklable worker designed for use with multiprocessing.Pool.map. It integrates the dynamics from a starting position inside start_basin until target_basin is reached (returns the elapsed time) or max_time is exceeded (returns None).

Parameters:

args (tuple) –

Packed positional arguments (for multiprocessing compatibility):

0. potential : callable – (U, F) = potential(x)

1. start_basin : callable – bool = start_basin(x)

2. target_basin : callable – bool = target_basin(x)

3. x0 : array_like – initial position

4. v0 : array_like – initial velocity

5. max_time : float – maximum simulation time

6. dt : float – integration time-step

7. gamma : float – friction

8. kT : float – thermal energy

9. m : float – mass

10. regime : {‘underdamped’, ‘overdamped’}

11. diffusion : scalar, callable, or None

Optional (appended for newer code):

12. bounds : sequence of (lo, hi) or None

13. boundary : {‘reflect’, ‘clip’} or None

14. seed : int or None – per-trajectory RNG seed

15. overdamped_eps : float – FD step for ∇D

Returns:

First-passage time, or None if max_time was exceeded.

Return type:

float or None

stochkin.mfpt.generate_basin_position(basin_func, bounds, max_attempts=1000)

Sample a random 2D position inside a basin (rejection sampling).

Parameters:

• basin_func (callable) – basin_func(pos) -> bool indicating membership.

• bounds (array_like, shape (2, 2)) – [[xmin, xmax], [ymin, ymax]] bounding box.

• max_attempts (int) – Maximum rejection-sampling iterations.

Returns:

A point inside the basin, or None if sampling failed.

Return type:

ndarray of shape (2,) or None

stochkin.mfpt.generate_basin_position_1d(basin_func, bounds, max_attempts=1000)

1D variant of generate_basin_position for x in [xmin, xmax].

Parameters:

• basin_func (callable) – Returns True if position pos = np.array([x]) is in the basin.

• bounds ((xmin, xmax)) – Interval from which to sample candidate x.

• max_attempts (int) – Maximum attempts to find a valid position.

Return type:

np.ndarray of shape (1,) or None if no valid position found.

stochkin.mfpt.compute_mfpt(potential, start_basin, target_basin, n_trials=1000, max_time=1000.0, dt=0.01, gamma=10.0, kT=0.05, m=1.0, start_bounds=None, initial_velocity=(0.0, 0.0), processes=None, verbose=True, regime='underdamped', diffusion=None, overdamped_eps=1e-06, bounds=None, boundary=None, base_seed=None)

Compute the mean first-passage time from start_basin to target_basin.

Launches n_trials independent trajectories (Langevin or overdamped BD) from random positions inside start_basin and records the time of first entry into target_basin.

Parameters:

• potential (callable) – potential(x) -> (U, F).

• start_basin (callable) – basin(x) -> bool.

• target_basin (callable) – basin(x) -> bool.

• n_trials (int) – Number of independent trajectories.

• max_time (float) – Maximum simulation time per trajectory.

• dt (float) – Integration time step.

• gamma (float) – Friction coefficient.

• kT (float) – Thermal energy.

• m (float) – Mass.

• start_bounds (sequence of (lo, hi) or None) – Bounding box for initial-position sampling.

• initial_velocity (array_like or 'thermal') – Fixed initial velocity vector, or 'thermal' to draw from the Maxwell–Boltzmann distribution.

• processes (int or None) – Worker processes (None = single process).

• verbose (bool) – Print progress.

• regime ({'underdamped', 'overdamped'}) – Dynamics type.

• diffusion (scalar, callable, or None) – Diffusion coefficient for overdamped mode.

• overdamped_eps (float) – FD step for ∇D.

• bounds (sequence of (lo, hi) or None) – Domain bounds (essential for cropped-FES workflows).

• boundary ({'reflect', 'clip'} or None) – Bound enforcement method.

• base_seed (int or None) – Base seed for reproducible per-trajectory seeds.

Returns:

Keys: 'mean', 'std', 'successful_trials', 'passage_times' (list), 'success_rate', 'total_trials'.

Return type:

dict

stochkin.mfpt.compute_mfpt_1d(potential, start_basin, target_basin, n_trials=1000, max_time=1000.0, dt=0.01, gamma=10.0, kT=0.05, m=1.0, start_bounds=None, initial_velocity=(0.0,), processes=None, verbose=True, regime='underdamped', diffusion=None, overdamped_eps=1e-06, bounds=None, boundary=None, base_seed=None)

1D analogue of compute_mfpt().

Identical interface to compute_mfpt() but initial positions are 1D scalars and start_bounds is a simple (xmin, xmax) interval. Supports domain bounds for cropped 1D free-energy surfaces.

Returns:

Same keys as compute_mfpt().

Return type:

dict

stochkin.mfpt.compute_bidirectional_mfpt(potential, basinA, basinB, n_trials=500, max_time=1000.0, dt=0.01, gamma=10.0, kT=0.05, m=1.0, boundsA=None, boundsB=None, initial_velocity=(0.0, 0.0), processes=None, verbose=True, regime='underdamped', diffusion=None)

Compute MFPT in both directions: A→B and B→A.

Convenience wrapper that calls compute_mfpt() twice.

Parameters:

• potential (callable) – potential(x) -> (U, F).

• basinA (callable) – Basin membership tests.

• basinB (callable) – Basin membership tests.

• n_trials (int) – Trajectories per direction.

• max_time (float) – Dynamics parameters.

• dt (float) – Dynamics parameters.

• gamma (float) – Dynamics parameters.

• kT (float) – Dynamics parameters.

• m (float) – Dynamics parameters.

• boundsA (sequence of (lo, hi) or None) – Sampling boxes for each basin.

• boundsB (sequence of (lo, hi) or None) – Sampling boxes for each basin.

• initial_velocity (array_like or 'thermal') – Shared initial-velocity setting.

• processes (int or None) – Worker processes.

• verbose (bool) – Print progress.

• regime ({'underdamped', 'overdamped'}) – Dynamics type.

• diffusion (scalar, callable, or None) – Diffusion for overdamped mode.

Returns:

{'A_to_B': <mfpt_dict>, 'B_to_A': <mfpt_dict>}.

Return type:

dict

stochkin.mfpt.compute_bidirectional_mfpt_1d(potential, basinA, basinB, n_trials=1000, max_time=1000.0, dt=0.01, gamma=10.0, kT=0.05, m=1.0, boundsA=None, boundsB=None, initial_velocity=0.0, processes=None, verbose=True, regime='underdamped', diffusion=None)

1D analogue of compute_bidirectional_mfpt().

Computes MFPT in both directions (A→B and B→A) using compute_mfpt_1d().

Returns:

results_AB, results_BA – Two MFPT result dicts as returned by compute_mfpt_1d().

Return type:

dict, dict

stochkin.mfpt.plot_mfpt_statistics(mfpt_results, title='MFPT Results')

Plot histograms of passage-time distributions.

Handles both single-direction results (one panel) and bidirectional results (two panels side by side).

Parameters:

• mfpt_results (dict) – Output of compute_mfpt(), compute_mfpt_1d(), or compute_bidirectional_mfpt(). If the dict contains 'A_to_B' and 'B_to_A' keys, two histograms are plotted.

• title (str) – Figure super-title.

stochkin.mfpt.estimate_transition_rates(mfpt_results, verbose=True)

Estimate transition rates from MFPT results (simple two-state model).

For a single direction, \(k = 1/\tau\). For bidirectional results, also computes the equilibrium constant \(K_{eq} = k_{AB}/k_{BA}\) and the steady-state populations \(p_A, p_B\).

Parameters:

• mfpt_results (dict) – Output of compute_mfpt() or compute_bidirectional_mfpt().

• verbose (bool) – Print summary.

Returns:

For bidirectional: {'k_AB', 'k_BA', 'K_eq', 'p_A_eq', 'p_B_eq'}. For single direction: {'rate'}. None if MFPTs are NaN.

Return type:

dict or None

stochkin.mfpt.basinA_mfpt(x, rA=0.3)

Basin A for MFPT: central well around the origin.

stochkin.mfpt.basinB_mfpt(x, r_inner=0.9, r_outer=1.2)

Basin B for MFPT: ring around the origin.

stochkin.mfpt.dw1d_basin_left_mfpt(x, x_cut=0.0)

Left basin for 1D double-well in MFPT calculations: x < x_cut (default: negative well).

stochkin.mfpt.dw1d_basin_right_mfpt(x, x_cut=0.0)

Right basin for 1D double-well in MFPT calculations: x > x_cut (default: positive well).

stochkin.mfpt.compute_mfpt_network(potential, basin_network, dt=0.01, max_time=100.0, D=None, beta=None, bounds=None, boundary='reflect', trials_per_basin=None, n_procs=None, seed=None, batch_size=None, mp_start_method=None, **legacy_kwargs)

Estimate a multi-basin MFPT network via overdamped BD first-exit simulations.

For each basin in basin_network, trials_per_basin overdamped Brownian dynamics trajectories are launched from the basin minimum. The first basin hit (if any) within max_time is recorded. The returned dictionary contains MFPT matrices, exit-count tables, and first-exit statistics suitable for CTMC rate-matrix estimation via estimate_rate_matrix().

Parameters:

• potential (callable) – potential(x) -> (U, F).

• basin_network (BasinNetwork) – Multi-basin partition on a 2D grid.

• dt (float) – Integration time step.

• max_time (float) – Maximum simulation time per trajectory.

• D (float or ndarray or None) – Diffusion coefficient (default 1.0).

• beta (float or None) – Inverse temperature \(1/(k_BT)\) (default 1.0).

• bounds (sequence of (lo, hi) or None) – Domain bounds.

• boundary (str) – Bound enforcement mode (default 'reflect').

• trials_per_basin (int or None) – Number of trajectories per basin (default 1000).

• n_procs (int or None) – Worker processes (default 1 = serial).

• seed (int or None) – Master seed for reproducibility.

• batch_size (int or None) – Trajectories per pool task; increase (20–100) to reduce IPC overhead.

• mp_start_method (str or None) – Force multiprocessing context (e.g. 'fork', 'spawn').

Returns:

• dict – Keys include 'mfpt', 'mfpt_matrix', 'n_samples', 'exit_to_counts', 'censored_counts', 'transition_times', 'first_exit_times', 'attempts_per_basin', 'params', 'method'.

• Performance Notes

• —————–

• When *n_procs > 1, a multiprocessing pool is created with an*

• initializer so that large objects (potential, basin_network) are

• transmitted only once per worker. Trials are grouped into small

• batches to reduce inter-process communication overhead.

stochkin.mfpt.compute_mfpt_network_fpe(basin_network, D, beta, initial_distribution='boltzmann', sparse=True, verbose=True)

Compute MFPTs between all basins in a BasinNetwork by solving the discrete Fokker–Planck / backward Kolmogorov equation on the 2D grid.

Parameters:

• basin_network (BasinNetwork) – Basin partition built on top of a 2D FES grid (xs, ys, U, labels).

• D (float or 2D array (nx, ny)) – Scalar diffusion coefficient D(x,y) on the same grid as U.

• beta (float) – Inverse temperature 1 / (k_B T). Must be consistent with the units of U.

• initial_distribution ({"boltzmann", "uniform"}, optional) – How to initialise within each basin when averaging MFPTs. “boltzmann” (default) uses exp(-beta U) restricted to the basin, “uniform” uses a flat distribution over grid points in the basin.

• sparse (bool, optional) – Passed through to build_fp_generator_from_fes.

• verbose (bool, optional) – If True, prints a small table of MFPTs for each target basin.

Returns:

results – Compatible with estimate_rate_matrix, with keys: - “n_basins” - “mfpt_matrix” - “std_matrix” (zeros, placeholder) - “generator” (the FP generator L) - “beta” - “D” - “method” = “fpe_grid”

Return type:

dict

stochkin.mfpt.estimate_rate_matrix(mfpt_network_results, verbose=True, method='auto')

Estimate a continuous-time Markov chain (CTMC) generator K between basins.

Preferred (and default) behavior:

• If first-exit CTMC statistics are available (from trajectory MFPT-network runs), return the consistent generator constructed as:

  k_out(i) = 1 / E[t_exit | exit] p_ij = P(next basin = j | exit from i) K_ij = k_out(i) * p_ij (j != i) K_ii = -sum_{j != i} K_ij

  This is the proper CTMC model implied by first-exit trajectories.

Fallback:

• If those stats are not present (or method=”inverse_mfpt”), estimate pairwise rates as k_ij ≈ 1 / MFPT(i->j). This is generally NOT a consistent generator for multi-basin systems, but is kept for backward compatibility.

Parameters:

• mfpt_network_results (dict) – Output of compute_mfpt_network / compute_mfpt_network_fpe / etc.

• verbose (bool) – Print diagnostics.

• method ({"auto", "ctmc", "inverse_mfpt"}) – Selection of estimator.

Returns:

K – CTMC generator (rows sum to zero).

Return type:

(n_basins, n_basins) ndarray