"""bngsim.NfsimSession — Stateful session API for network-free simulation.
Wraps the C++ NfsimSimulator to provide a safe, Pythonic interface for
multi-action workflows: initialize, simulate a segment, modify molecules
or parameters, simulate again, and clean up.
Supports the context manager protocol for guaranteed resource cleanup.
Example
-------
>>> with bngsim.NfsimSession("model.xml") as session:
... session.set_param("kf", 0.5)
... session.initialize(seed=42)
... r1 = session.simulate(0, 100, n_points=101)
... session.add_species("L(r)", 500)
... r2 = session.simulate(100, 200, n_points=101)
"""
from __future__ import annotations
import contextlib
import logging
from pathlib import Path
from bngsim._exceptions import ParameterError, SimulationError, SimulationTimeout
from bngsim._result import Result
from bngsim._rounding import round_half_up
from bngsim._sample_times import normalize_sample_times
logger = logging.getLogger("bngsim")
[docs]
class NfsimSession:
"""Stateful NFsim simulation session with per-step control.
Use this class when you need fine-grained control over an NFsim
session: setting parameters before initialization, running multiple
simulation segments, and mutating molecule or exact species counts between
segments.
For simple one-shot NFsim simulations, ``bngsim.Simulator(model,
method="nf", xml_path=...)`` is easier.
Parameters
----------
xml_path : str or Path
Path to a BNG XML file (generated by BNG2.pl ``writeXML``).
molecule_limit : int, optional
Global molecule limit (NFsim ``-gml`` flag).
connectivity : bool, optional
Enable NFsim's reaction-connectivity optimization.
v1143_compat : bool, optional
Enable NFsim v1.14.3 selector compatibility mode for
same-seed trajectory reproducibility with the standalone CLI.
block_same_complex_binding : bool, optional
NFsim ``-bscb``: when True, block bimolecular reactions whose two
reactant patterns match molecules in the same complex. Default:
``True`` in bngsim — the NFsim CLI default is off, but bngsim defaults
it on for correctness on BLBR/aggregation models that polymerize within
a single complex otherwise. Pass ``False`` to allow same-complex
binding (matching BNG2.pl's ``complex=>1``). Note: this flag governs
only the binding policy; complex bookkeeping needed to count
``Species``-typed observables is enabled automatically whenever the
model declares one, regardless of this setting.
traversal_limit : int, optional
NFsim ``-utl N``: universal traversal limit. ``None`` (default)
lets NFsim auto-compute a suggested limit from the XML.
Raises
------
RuntimeError
If NFsim support is not available in this bngsim build.
Examples
--------
Basic one-shot simulation::
with bngsim.NfsimSession("model.xml") as nf:
nf.initialize(seed=42)
result = nf.simulate(0, 100, n_points=101)
Multi-step protocol::
with bngsim.NfsimSession("model.xml", molecule_limit=100000) as nf:
nf.set_param("kf", 0.5)
nf.initialize(seed=42)
r1 = nf.simulate(0, 50, n_points=51)
nf.add_species("L(r)", 500)
r2 = nf.simulate(50, 100, n_points=51)
Parameter scan (new session per point)::
for kf_val in [0.1, 1.0, 10.0]:
with bngsim.NfsimSession("model.xml") as nf:
nf.set_param("kf", kf_val)
nf.initialize(seed=42)
result = nf.simulate(0, 100, n_points=101)
"""
__slots__ = (
"_core",
"_xml_path",
"_initialized",
"_destroyed",
"_seed",
)
def __init__(
self,
xml_path: str | Path,
*,
molecule_limit: int | None = None,
connectivity: bool | None = None,
v1143_compat: bool = False,
block_same_complex_binding: bool = True,
traversal_limit: int | None = None,
) -> None:
from bngsim._bngsim_core import HAS_NFSIM
if not HAS_NFSIM:
raise RuntimeError(
"NFsim support is not present in this bngsim install. "
"The vendored NFsim backend at third_party/nfsim/ is built "
"by default; this install was either configured with "
"-DBNGSIM_BUILD_NFSIM=OFF or installed from a wheel that "
"excludes NFsim."
)
from bngsim._bngsim_core import NfsimSimulator
self._xml_path = str(xml_path)
self._initialized = False
self._destroyed = False
self._core = NfsimSimulator(self._xml_path)
if molecule_limit is not None:
self._core.set_molecule_limit(int(molecule_limit))
if connectivity is not None:
self._core.set_connectivity(bool(connectivity))
if v1143_compat:
self._core.set_nfsim_v1143_compat(True)
# The C++ default is also True; we always pass through so explicit
# False propagates correctly.
self._core.set_block_same_complex_binding(bool(block_same_complex_binding))
if traversal_limit is not None:
self._core.set_traversal_limit(int(traversal_limit))
self._seed: int | None = None
logger.debug("NfsimSession created: %s", self._xml_path)
# ── Context manager ──────────────────────────────────────────
def __enter__(self) -> NfsimSession:
return self
def __exit__(self, exc_type, exc_val, exc_tb) -> None:
self.destroy()
# ── Session lifecycle ────────────────────────────────────────
[docs]
def set_molecule_limit(self, limit: int) -> None:
"""Set the global molecule limit.
Can be called before or after ``initialize()`` to change
the NFsim ``-gml`` cap.
Parameters
----------
limit : int
Maximum number of molecules per molecule type.
"""
self._require_alive()
self._core.set_molecule_limit(int(limit))
[docs]
def set_param(self, name: str, value: float) -> None:
"""Set a parameter value and propagate to dependents.
May be called before *or* after ``initialize()``. Pre-init writes
are stored and applied on session start; post-init writes update
the live NFsim System in place. In both cases bngsim re-evaluates
every parameter whose BNG XML ``expr=`` transitively depends on
``name``, so derived values like ``LT = LT_conc_M*NA*V_sim`` and
``_rateLaw*`` follow automatically.
Pre-init writes are baked into the XML that NFsim parses, so any
``<Species concentration="X">`` seed-species expression that
references ``name`` (directly or through a derived parameter) is
evaluated against the override-resolved namespace when the agent
population is created. This matches BNG2.pl's behavior when
``setParameter`` is invoked between ``generate_network`` and a
``method=>"nf"`` simulate.
Post-init writes do *not* alter the live agent population — they
only update reaction rates from the new parameter namespace. To
change live counts after ``initialize()``, use
:meth:`add_species` / :meth:`remove_species`.
Parameters
----------
name : str
Parameter name.
value : float
Parameter value.
Raises
------
ParameterError
If the parameter name is not recognized or if expression
propagation fails for some dependent parameter.
"""
self._require_alive()
try:
self._core.set_param(name, float(value))
except (RuntimeError, ValueError) as e:
raise ParameterError(f"Failed to set NFsim parameter '{name}': {e}") from e
[docs]
def clear_param_overrides(self) -> None:
"""Clear all parameter overrides set via ``set_param``.
When called on a live session, parameters revert to their
XML-time values and dependent expressions are re-evaluated so
the System matches the original BNG XML namespace.
"""
self._require_alive()
self._core.clear_param_overrides()
[docs]
def evaluate(self, expr: str, overrides: dict[str, float] | None = None) -> float:
"""Evaluate a BNG expression in the current parameter namespace.
Uses the same ExprTk evaluator that powers ``set_param``
propagation, so expressions can reference any BNG parameter,
BNG built-in (``_pi``, ``_NA``, ...), or BNG-style construct
(``if(cond, then, else)``, ``&&``, ``||``, ``mratio``).
Parameters
----------
expr : str
Expression string. May reference any BNG parameter name.
overrides : dict[str, float], optional
Extra ``name → value`` bindings layered on top of the
simulator's persistent ``set_param`` overrides for this
single evaluation. Useful for "what if X were Y?" probes
without mutating session state.
Returns
-------
float
Numeric value of the expression.
Raises
------
ParameterError
If the expression fails to compile.
"""
self._require_alive()
try:
return self._core.evaluate_expression(expr, overrides or {})
except (RuntimeError, ValueError) as e:
raise ParameterError(f"Failed to evaluate expression '{expr}': {e}") from e
[docs]
def initialize(self, seed: int | None = None) -> None:
"""Initialize the NFsim session.
Parses the XML, seeds the RNG, and prepares the simulation
system. Must be called before ``simulate()``.
Parameters
----------
seed : int, optional
Random number generator seed. When omitted (or ``None``),
bngsim draws a fresh seed from system entropy so consecutive
session lifetimes produce independent trajectories. Pass an
explicit integer for reproducibility. The actual seed used
is exposed via the :attr:`seed` property and stamped onto
every ``Result`` returned by :meth:`simulate`.
Raises
------
SimulationError
If initialization fails.
"""
from bngsim._seed import _resolve_seed
self._require_alive()
used_seed = _resolve_seed(seed)
try:
self._core.initialize(used_seed)
except RuntimeError as e:
raise SimulationError(f"NFsim initialization failed: {e}") from e
self._initialized = True
self._seed = used_seed
logger.info("NfsimSession initialized (seed=%d)", used_seed)
[docs]
def destroy(self) -> None:
"""Destroy the session and release C++ resources.
Safe to call multiple times. Called automatically by
``__exit__`` when using the context manager.
"""
if self._destroyed:
return
with contextlib.suppress(Exception):
self._core.destroy_session()
self._destroyed = True
self._initialized = False
logger.debug("NfsimSession destroyed")
# ── Simulation ───────────────────────────────────────────────
[docs]
def simulate(
self,
t_start: float | None = None,
t_end: float | None = None,
n_points: int | None = None,
*,
timeout: float | None = None,
relative_time: bool = False,
sample_times: list[float] | None = None,
) -> Result:
"""Run a simulation segment and return a Result.
Session state is preserved across calls, enabling multi-step
protocols (simulate, mutate, simulate again).
Parameters
----------
t_start : float, optional
Segment start time. Required unless ``sample_times`` is given.
t_end : float, optional
Segment end time. Required unless ``sample_times`` is given.
n_points : int, optional
Number of output time points (including endpoints). Required
unless ``sample_times`` is given.
sample_times : list[float], optional
Explicit output instants (GH #184). When provided, overrides
``t_start`` / ``t_end`` / ``n_points``: the returned
:attr:`Result.time` array equals these times exactly (sorted
ascending, treated as absolute), instead of a uniform grid.
The live NFsim clock advances from the current session time by
the span of the list, so this works mid-protocol (after
``setConcentration`` / a prior segment), continuing the same
trajectory. Must contain at least two finite points. Mutually
exclusive with ``relative_time`` (sample_times are absolute).
timeout : float, optional
Wall-clock budget in seconds. ``None`` or ``0`` disables the
budget. Positive values arm a check between each ``stepTo()``
output point; on overrun this method raises
:class:`bngsim.SimulationTimeout` and the session must be
destroyed (or the context manager exited) before reuse — the
live NFsim System has advanced an indeterminate distance into
the segment and the session clock is not updated.
relative_time : bool, default False
When ``True``, the returned :attr:`Result.time` array is
offset so it begins at ``0.0`` and ends at ``t_end -
t_start`` — the convention used by BNG2.pl's
``simulate({method=>"nf", ...})`` action, which prints
"NFsim timepoints are reported as time elapsed since
``t_start=$t_start``" at runtime. The internal NFsim clock
and ``session_logical_time`` are unaffected, so multi-segment
threading (simulate, mutate, simulate again) keeps the same
absolute physical time it would have under the default.
Recommended for hosts that thread NF protocols across
multiple bngsim ``simulate()`` calls and want their stitched
time axis to match BNG2.pl's. Default ``False`` preserves
the absolute-time labelling used by every other bngsim
backend.
Returns
-------
Result
Simulation result for this segment.
Raises
------
SimulationError
If the session is not initialized or simulation fails.
SimulationTimeout
If ``timeout`` is set and the wall-clock budget is exceeded.
"""
self._require_initialized()
# Mirror Simulator.run: None or non-positive disables the budget.
timeout_seconds: float = 0.0
if timeout is not None:
timeout_seconds = float(timeout)
if timeout_seconds < 0.0:
raise ValueError(f"timeout must be non-negative or None, got {timeout!r}")
# Resolve the output schedule. Explicit sample_times (GH #184) override
# t_start/t_end/n_points; the C++ core ignores those in that branch but
# still takes them as positional args, so pass the resolved endpoints.
explicit_times: list[float] = []
if sample_times is not None:
if relative_time:
raise ValueError(
"sample_times and relative_time are mutually exclusive: "
"sample_times are treated as absolute output instants."
)
explicit_times = normalize_sample_times(sample_times)
t_start = explicit_times[0]
t_end = explicit_times[-1]
n_points = len(explicit_times)
else:
if t_start is None or t_end is None or n_points is None:
raise ValueError(
"t_start, t_end and n_points are required unless sample_times is given."
)
try:
core_result = self._core.simulate(
float(t_start),
float(t_end),
int(n_points),
timeout_seconds,
bool(relative_time),
explicit_times,
)
except SimulationTimeout:
# Already a typed bngsim exception (raised via the C++ translator).
# Pass through unchanged so callers can classify wall-clock
# terminations distinctly from solver errors.
raise
except RuntimeError as e:
raise SimulationError(f"NFsim simulation failed: {e}") from e
result = Result(core_result)
result._seed = self._seed
return result
[docs]
def step_to(self, time: float) -> None:
"""Advance the simulation to *time* without recording output.
Useful for fast-forwarding before a protocol action.
Parameters
----------
time : float
Target time to advance to.
Raises
------
SimulationError
If the session is not initialized.
"""
self._require_initialized()
try:
self._core.step_to(float(time))
except RuntimeError as e:
raise SimulationError(f"NFsim step_to failed: {e}") from e
# ── Session queries & mutations ──────────────────────────────
[docs]
def get_molecule_count(self, mol_type: str) -> int:
"""Get the current count of a molecule type.
Parameters
----------
mol_type : str
Molecule type name (e.g. ``"L"``, ``"R"``).
Returns
-------
int
Current molecule count.
"""
self._require_initialized()
return self._core.get_molecule_count(mol_type)
[docs]
def add_molecules(self, mol_type: str, count: int) -> None:
"""Add molecules of a given type in their default (unbound) state.
Parameters
----------
mol_type : str
Molecule type name.
count : int
Number of molecules to add. Must be > 0. A fractional value is
rounded to the nearest integer (round-half-up, GH #51) so warm /
dosing / scan re-seeds match bngsim SSA and the cold-start seed.
Raises
------
ValueError
If *count* is not positive.
"""
self._require_initialized()
count = round_half_up(count)
if count <= 0:
raise ValueError(
f"count must be positive, got {count}. "
"Use remove_species() for exact species removal."
)
self._core.add_molecules(mol_type, count)
[docs]
def get_species_count(self, pattern: str) -> int:
"""Get the current count of an exact BNGL species pattern.
The current implementation supports exact single-molecule species
patterns such as ``"X(p~0,y)"`` or ``"TNF()"``. Mutation-oriented
species patterns must list every component and specify every stateful
component state so that the live-session operation is unambiguous.
Multi-molecule complex patterns fail clearly.
Symmetric components (e.g. ``L(r~u~c~g, r~u~c~g)`` declares two sites
named ``r``) are matched by sorted multiset of state/bond constraints
across class members, so heterogeneous patterns like
``L(r~u,r~c)`` and ``L(r~c,r~u)`` resolve to the same species count.
Parameters
----------
pattern : str
BNGL species pattern.
Returns
-------
int
Current exact species count.
"""
self._require_initialized()
try:
return self._core.get_species_count(pattern)
except RuntimeError as e:
raise SimulationError(f"NFsim get_species_count failed: {e}") from e
[docs]
def add_species(self, pattern: str, count: int) -> None:
"""Add exact BNGL species instances to the live NFsim session.
Parameters
----------
pattern : str
Exact single-molecule BNGL species pattern.
count : int
Number of instances to add. Must be > 0. A fractional value is
rounded to the nearest integer (round-half-up, GH #51).
"""
self._require_initialized()
count = round_half_up(count)
if count <= 0:
raise ValueError(f"count must be positive, got {count}")
try:
self._core.add_species(pattern, count)
except RuntimeError as e:
raise SimulationError(f"NFsim add_species failed: {e}") from e
[docs]
def remove_species(self, pattern: str, count: int) -> None:
"""Remove exact BNGL species instances from the live NFsim session.
Parameters
----------
pattern : str
Exact single-molecule BNGL species pattern.
count : int
Number of instances to remove. Must be > 0. A fractional value is
rounded to the nearest integer (round-half-up, GH #51).
"""
self._require_initialized()
count = round_half_up(count)
if count <= 0:
raise ValueError(f"count must be positive, got {count}")
try:
self._core.remove_species(pattern, count)
except RuntimeError as e:
raise SimulationError(f"NFsim remove_species failed: {e}") from e
[docs]
def set_species_count(self, pattern: str, count: int) -> None:
"""Set the live count for an exact BNGL species pattern.
This is implemented as a diff against ``get_species_count(pattern)``
and may add or remove exact species instances. A fractional *count* is
rounded to the nearest integer (round-half-up, GH #51).
"""
self._require_initialized()
count = round_half_up(count)
if count < 0:
raise ValueError(f"count must be nonnegative, got {count}")
try:
self._core.set_species_count(pattern, count)
except RuntimeError as e:
raise SimulationError(f"NFsim set_species_count failed: {e}") from e
[docs]
def get_parameter(self, name: str) -> float:
"""Evaluate and return a parameter value from the live session.
Can evaluate expressions, not just named parameters.
Parameters
----------
name : str
Parameter or expression name.
Returns
-------
float
Evaluated value.
"""
self._require_initialized()
return self._core.get_parameter(name)
[docs]
def get_observable_names(self) -> list[str]:
"""Return observable names from the live session.
Returns
-------
list[str]
Observable names defined in the BNG XML.
"""
self._require_initialized()
return self._core.get_observable_names()
[docs]
def get_observable_values(self) -> list[float]:
"""Return current observable values from the live session.
Returns
-------
list[float]
Current observable values.
"""
self._require_initialized()
return self._core.get_observable_values()
[docs]
def save_species(self, path: str | Path) -> None:
"""Write the live System's molecular species to a BNG ``.species`` file.
Mirrors the artifact BNG2.pl emits when running
``simulate({method=>"nf", get_final_state=>1, ...})``. The file
is a plain-text listing of species patterns and counts,
formatted exactly as NFsim's standalone CLI writes it. A host
like PyBioNetGen that drives multi-segment NF protocols can read
this file to thread state across segments without needing
in-process snapshot support (see issue #52 for the in-process
path).
Parameters
----------
path : str or Path
Output file path. Overwritten if it exists.
Raises
------
SimulationError
If the session is not initialized or the file cannot be
written.
"""
self._require_initialized()
try:
self._core.save_species(str(path))
except RuntimeError as e:
raise SimulationError(f"NFsim save_species failed: {e}") from e
[docs]
def save_concentrations(self) -> None:
"""Snapshot the live session's full molecular state for later restore.
Captures every molecule's count, component states, and bonds into an
in-memory snapshot held by the session. A later
:meth:`restore_concentrations` call rewinds the session to exactly
this state. Overwrites any previous snapshot. Mirrors the BNG
``saveConcentrations()`` action.
Unlike :meth:`save_species` (which writes a ``.species`` file for
out-of-process state threading), this keeps the state in process —
useful for equilibrate → snapshot → perturb → restore workflows
without touching disk.
Raises
------
SimulationError
If the session is not initialized.
"""
self._require_initialized()
try:
self._core.save_concentrations()
except RuntimeError as e:
raise SimulationError(f"NFsim save_concentrations failed: {e}") from e
[docs]
def restore_concentrations(self) -> None:
"""Restore the molecular state captured by :meth:`save_concentrations`.
Rewinds the live session to the most recently snapshotted state
(counts, component states, and bonds), discarding any simulation
progress since the snapshot. Mirrors the BNG ``resetConcentrations()``
action.
Raises
------
SimulationError
If the session is not initialized, or no snapshot has been saved
with :meth:`save_concentrations`.
"""
self._require_initialized()
if not self._core.has_saved_concentrations():
raise SimulationError(
"No saved concentrations to restore. "
"Call save_concentrations() before restore_concentrations()."
)
try:
self._core.restore_concentrations()
except RuntimeError as e:
raise SimulationError(f"NFsim restore_concentrations failed: {e}") from e
[docs]
def has_saved_concentrations(self) -> bool:
"""Whether a :meth:`save_concentrations` snapshot is available to restore."""
self._require_alive()
if not self._initialized:
return False
return self._core.has_saved_concentrations()
# ── Properties ───────────────────────────────────────────────
@property
def initialized(self) -> bool:
"""Whether ``initialize()`` has been called."""
return self._initialized
@property
def seed(self) -> int | None:
"""The integer RNG seed used by ``initialize()``, or ``None``
if the session has not been initialized yet."""
return self._seed
@property
def destroyed(self) -> bool:
"""Whether ``destroy()`` has been called."""
return self._destroyed
@property
def xml_path(self) -> str:
"""Path to the BNG XML file."""
return self._xml_path
# ── Internal helpers ─────────────────────────────────────────
def _require_alive(self) -> None:
if self._destroyed:
raise SimulationError(
"NfsimSession has been destroyed. Create a new session to continue."
)
def _require_initialized(self) -> None:
self._require_alive()
if not self._initialized:
raise SimulationError("NfsimSession is not initialized. Call initialize(seed) first.")
def __del__(self) -> None:
if not self._destroyed:
self.destroy()
def __repr__(self) -> str:
state = (
"destroyed" if self._destroyed else "initialized" if self._initialized else "created"
)
return f"NfsimSession(xml_path={self._xml_path!r}, state={state})"