JumpProcesses 10

.github/workflows/Downgrade.yml

@@ -3,6 +3,7 @@ on:
   pull_request:
     branches:
       - master
+      - v10_working
     paths-ignore:
       - 'docs/**'
   push:

.github/workflows/Tests.yml

@@ -8,6 +8,7 @@ on:
   pull_request:
     branches:
       - master
+      - v10_working
     paths-ignore:
       - 'docs/**'
 
@@ -29,6 +30,7 @@ jobs:
         group:
           - InterfaceI
           - InterfaceII
+          - ThreadSafety
           - QA
         exclude:
           - version: "pre"

.github/workflows/ThreadSafety.yml

@@ -0,0 +1,66 @@
+name: "Thread Safety Tests"
+
+on:
+  push:
+    branches:
+      - master
+    paths-ignore:
+      - 'docs/**'
+  pull_request:
+    branches:
+      - master
+      - v10_working
+    paths-ignore:
+      - 'docs/**'
+
+concurrency:
+  group: ${{ github.workflow }}-${{ github.ref }}
+  cancel-in-progress: ${{ github.ref_name != github.event.repository.default_branch || github.ref != 'refs/tags/v*' }}
+
+jobs:
+  thread-safety:
+    name: "Thread Safety"
+    strategy:
+      fail-fast: false
+      matrix:
+        version:
+          - "1"
+          - "lts"
+          - "pre"
+    runs-on: ubuntu-latest
+    steps:
+      - uses: actions/checkout@v4
+
+      - name: "Setup Julia ${{ matrix.version }}"
+        uses: julia-actions/setup-julia@v2
+        with:
+          version: "${{ matrix.version }}"
+
+      - name: "Cache Julia packages"
+        uses: julia-actions/cache@v2
+        with:
+          token: "${{ secrets.GITHUB_TOKEN }}"
+
+      - name: "Build package"
+        uses: julia-actions/julia-buildpkg@v1
+
+      - name: "Run thread safety tests (4 threads)"
+        run: |
+          julia --threads=4 --code-coverage=user --check-bounds=yes --compiled-modules=yes \
+                --project=@. --color=yes -e '
+            using Pkg
+            Pkg.test()
+          '
+        env:
+          GROUP: ThreadSafety
+
+      - name: "Process Coverage"
+        uses: julia-actions/julia-processcoverage@v1
+
+      - name: "Report Coverage"
+        uses: codecov/codecov-action@v5
+        if: ${{ github.event.pull_request.head.repo.full_name == github.repository }}
+        with:
+          files: lcov.info
+          token: "${{ secrets.CODECOV_TOKEN }}"
+          fail_ci_if_error: false

HISTORY.md

@@ -1,6 +1,51 @@
 # Breaking updates and feature summaries across releases
 
-## JumpProcesses unreleased (master branch)
+## 10.0 (Breaking)
+
+  - **Breaking**: The `rng` keyword argument has been removed from
+    `JumpProblem`. Pass `rng` to `solve` or `init` instead:
+    ```julia
+    # Before (no longer works):
+    jprob = JumpProblem(dprob, Direct(), jump; rng = Xoshiro(1234))
+    sol = solve(jprob, SSAStepper())
+
+    # After:
+    jprob = JumpProblem(dprob, Direct(), jump)
+    sol = solve(jprob, SSAStepper(); rng = Xoshiro(1234))
+    ```
+  - RNG state is now owned by the integrator, not the aggregator. This
+    eliminates data races when sharing a `JumpProblem` across threads and
+    ensures a single, consistent RNG priority across all solver pathways:
+    `rng` > `seed` > `Random.default_rng()`.
+  - `rng` and `seed` kwargs are fully supported on `solve`/`init` for all
+    solver pathways (SSAStepper, ODE, SDE, tau-leaping).
+  - `SSAIntegrator` now supports the `SciMLBase` RNG interface (`has_rng`,
+    `get_rng`, `set_rng!`).
+  - **Breaking**: The `scale_rates` and `useiszero` keyword arguments have been
+    removed from `JumpProblem`. Set them on the `MassActionJump` directly:
+    ```julia
+    # Before (no longer works):
+    jprob = JumpProblem(dprob, Direct(), maj; scale_rates = false)
+
+    # After:
+    maj = MassActionJump(rates, reactant_stoch, net_stoch; scale_rates = false)
+    jprob = JumpProblem(dprob, Direct(), maj)
+    ```
+  - **Breaking**: Parameterized `MassActionJump`s (those constructed with
+    `param_idxs` or a custom `param_mapper`) are now immutable — rates are
+    computed from parameters at aggregator initialization rather than being
+    materialized into the jump at `JumpProblem` construction time. This means:
+      - `update_parameters!` has been removed. Mass action rates are now
+        automatically recomputed from the current parameter values whenever the
+        aggregator reinitializes. After modifying parameters (e.g. in a
+        callback), call `reset_aggregated_jumps!(integrator)` to trigger
+        reinitialization with the updated parameter values.
+      - Custom parameter mappers (e.g. ModelingToolkitBase's
+        `JumpSysMajParamMapper`) must implement the 3-arg callable API:
+        `(mapper)(dest::AbstractVector, maj::MassActionJump, params)`.
+        See [`MassActionJumpParamMapper`](@ref) for details.
+  - Scalar `param_idxs` (e.g. `param_idxs = 1`) is now internally converted to
+    a one-element vector. The scalar form continues to work as before.
 
 ## 9.14
 

Project.toml

@@ -54,15 +54,18 @@ Random = "1"
 RecursiveArrayTools = "3.35"
 Reexport = "1.2"
 SafeTestsets = "0.1"
-SciMLBase = "2.115"
+SciMLBase = "2.147"
 StableRNGs = "1"
 StaticArrays = "1.9.8"
 Statistics = "1"
-StochasticDiffEq = "6.82"
+StochasticDiffEq = "6.95"
 SymbolicIndexingInterface = "0.3.36"
 Test = "1"
 julia = "1.10"
 
+[sources]
+SciMLBase = {url = "https://github.com/isaacsas/SciMLBase.jl.git", rev = "ensemble_rng_redesign"}
+
 [extras]
 ADTypes = "47edcb42-4c32-4615-8424-f2b9edc5f35b"
 Aqua = "4c88cf16-eb10-579e-8560-4a9242c79595"

docs/src/api.md

@@ -40,9 +40,91 @@ RSSACR
 SortingDirect
 ```
 
-# Private API Functions
+## Random Number Generator Control
+
+JumpProcesses supports controlling the random number generator (RNG) used for
+jump sampling via the `rng` and `seed` keyword arguments to `solve` or `init`.
+
+### `rng` keyword argument
+
+Pass any `AbstractRNG` to `solve` or `init`:
+
+```julia
+using Random, StableRNGs
+
+# Using a StableRNG for cross-version reproducibility
+sol = solve(jprob, SSAStepper(); rng = StableRNG(1234))
+
+# Using Julia's built-in Xoshiro
+sol = solve(jprob, Tsit5(); rng = Xoshiro(42))
+```
+
+### `seed` keyword argument
+
+As a shorthand, pass an integer `seed` to create a `Xoshiro` generator:
+
+```julia
+sol = solve(jprob, SSAStepper(); seed = 1234)
+# equivalent to: solve(jprob, SSAStepper(); rng = Xoshiro(1234))
+```
+
+### Resolution priority
+
+When both `rng` and `seed` are passed to the same `solve`/`init` call, `rng`
+takes priority:
+
+| User provides | Result |
+|---|---|
+| `rng` via `solve`/`init` | Uses that `rng` |
+| `seed` via `solve`/`init` | Creates `Xoshiro(seed)` |
+| Nothing | Uses `Random.default_rng()` (SSAStepper, ODE, tau-leaping) or a randomly-seeded `Xoshiro` (SDE) |
+
+### Behavior by solver pathway
+
+| Solver | Default RNG (nothing passed) | `rng` / `seed` support |
+|---|---|---|
+| `SSAStepper` | `Random.default_rng()` | Full support via `solve`/`init` kwargs |
+| ODE solvers (e.g., `Tsit5`) | `Random.default_rng()` | Full support via `solve`/`init` kwargs |
+| SDE solvers (e.g., `SRIW1`) | Randomly-seeded `Xoshiro` | Full support; `TaskLocalRNG` is auto-converted to `Xoshiro` |
+| `SimpleTauLeaping` | `Random.default_rng()` | Full support via `solve` kwargs |
+
+!!! note
+    For reproducible simulations, always pass an explicit `rng` or `seed`.
+    The default RNG is shared global state and may produce different results
+    depending on prior usage.
+
+# Private / Developer API
 
 ```@docs
 ExtendedJumpArray
 SSAIntegrator
 ```
+
+## Internal Dispatch Pathways
+
+The following table documents which code handles `solve`/`init` for each solver
+type. This is relevant for developers working on JumpProcesses or its solver
+backends.
+
+| Solver type | `__solve` handled by | `__init` handled by | Uses `__jump_init`? |
+|---|---|---|---|
+| `SSAStepper` | JumpProcesses (`solve.jl`) | JumpProcesses (`SSA_stepper.jl`) | No |
+| ODE (e.g., `Tsit5`) | JumpProcesses (`solve.jl`) | JumpProcesses (`solve.jl`) → OrdinaryDiffEq | Yes |
+| SDE (e.g., `SRIW1`) | StochasticDiffEq | StochasticDiffEq | No |
+| `SimpleTauLeaping` | JumpProcesses (`simple_regular_solve.jl`, custom `DiffEqBase.solve`) | N/A | No |
+
+For **SSAStepper**, `rng` is resolved via `resolve_rng` in `SSA_stepper.jl`'s
+`__init` and stored on the [`SSAIntegrator`](@ref).
+
+For **ODE solvers**, `rng` is resolved via `resolve_rng` in `__jump_init`
+(`solve.jl`) and forwarded to OrdinaryDiffEq's `init`, which stores it on the
+`ODEIntegrator`.
+
+For **SDE solvers**, StochasticDiffEq handles the full solve/init pathway
+directly (JumpProcesses' ambiguity-fix `__solve` method is never dispatched to).
+StochasticDiffEq has its own `_resolve_rng` that additionally handles
+`TaskLocalRNG` conversion and the problem's stored seed.
+
+For **tau-leaping**, JumpProcesses defines a custom `DiffEqBase.solve` that
+bypasses the standard `__solve`/`__init` pathway. It calls `resolve_rng`
+directly with the `rng` and `seed` kwargs from the `solve` call.

docs/src/applications/advanced_point_process.md

@@ -387,10 +387,10 @@ function Base.rand(rng::AbstractRNG,
     out = Array{History, 1}(undef, n)
     p = params(pp)
     dprob = DiscreteProblem([0], tspan, p)
-    jprob = JumpProblem(dprob, Coevolve(), jumps...; dep_graph = pp.g, save_positions, rng)
+    jprob = JumpProblem(dprob, Coevolve(), jumps...; dep_graph = pp.g, save_positions)
     for i in 1:n
         params!(pp, p)
-        solve(jprob, SSAStepper())
+        solve(jprob, SSAStepper(); rng)
         out[i] = deepcopy(p.h)
     end
     return out

docs/src/faq.md

@@ -55,21 +55,27 @@ jset = JumpSet(; constant_jumps = cjvec, variable_jumps = vjtuple,
 
 ## How can I set the random number generator used in the jump process sampling algorithms (SSAs)?
 
-Random number generators can be passed to `JumpProblem` via the `rng` keyword
+Random number generators can be passed to `solve` or `init` via the `rng` keyword
 argument. Continuing the previous example:
 
 ```julia
-#] add RandomNumbers
-using RandomNumbers
-jprob = JumpProblem(dprob, Direct(), maj,
-    rng = Xorshifts.Xoroshiro128Star(rand(UInt64)))
+using Random
+jprob = JumpProblem(dprob, Direct(), maj)
+sol = solve(jprob, SSAStepper(); rng = Xoshiro(1234))
+```
+
+Any `AbstractRNG` can be used. For example, to use a generator from
+[StableRNGs.jl](https://github.com/JuliaRandom/StableRNGs.jl):
+
+```julia
+using StableRNGs
+sol = solve(jprob, SSAStepper(); rng = StableRNG(1234))
 ```
 
-uses the `Xoroshiro128Star` generator from
-[RandomNumbers.jl](https://github.com/JuliaRandom/RandomNumbers.jl).
+A `seed` keyword argument is also supported as a shorthand for creating a `Xoshiro`
+generator: `solve(jprob, SSAStepper(); seed = 1234)`.
 
-On version 1.7 and up, JumpProcesses uses Julia's built-in random number generator by
-default. On versions below 1.7 it uses `Xoroshiro128Star`.
+By default, JumpProcesses uses Julia's built-in `Random.default_rng()`.
 
 ## What are these aggregators and aggregations in JumpProcesses?
 

docs/src/tutorials/simple_poisson_process.md

@@ -371,11 +371,10 @@ with ``N(t)`` a Poisson counting process with constant transition rate
 ``\lambda``, and the ``C_i`` independent and identical samples from a uniform
 distribution over ``\{-1,1\}``. We can simulate such a process as follows.
 
-We first ensure that we use the same random number generator as JumpProcesses. We
-can either pass one as an input to [`JumpProblem`](@ref) via the `rng` keyword
-argument, and make sure it is the same one we use in our `affect!` function, or
-we can just use the default generator chosen by JumpProcesses if one is not
-specified, `JumpProcesses.DEFAULT_RNG`. Let's do the latter
+We first ensure that we use the same random number generator as JumpProcesses.
+Custom RNGs can be passed to `solve` or `init` via the `rng` keyword argument.
+If no RNG is specified, JumpProcesses uses `Random.default_rng()`, which is also
+available as `JumpProcesses.DEFAULT_RNG`. Let's use the default
 
 ```@example tut1
 rng = JumpProcesses.DEFAULT_RNG

src/JumpProcesses.jl

@@ -5,7 +5,7 @@ using Reexport: Reexport, @reexport
 
 # Explicit imports from standard libraries
 using LinearAlgebra: LinearAlgebra, mul!
-using Random: Random, randexp, seed!
+using Random: Random, randexp
 
 # Explicit imports from external packages
 using DocStringExtensions: DocStringExtensions, FIELDS, TYPEDEF
@@ -40,7 +40,7 @@ using DiffEqBase: DiffEqBase, CallbackSet, ContinuousCallback, DAEFunction,
                   ODESolution, ReturnCode, SDEFunction, SDEProblem, add_tstop!,
                   deleteat!, isinplace, remake, savevalues!, step!,
                   u_modified!
-using SciMLBase: SciMLBase, DEIntegrator
+using SciMLBase: SciMLBase, DEIntegrator, has_rng, get_rng, set_rng!
 
 abstract type AbstractJump end
 abstract type AbstractMassActionJump <: AbstractJump end