JumpProcess v10 Attempt 2

.github/workflows/Downgrade.yml

@@ -13,7 +13,11 @@ on:
 jobs:
   downgrade:
     name: "Downgrade"
-    # Runs strict (allow_reresolve=false); expected RED until the upstream fix
-    # ForwardDiff/LogExpFunctions — red until ForwardDiff 1.3.4 (#811). Auto-greens on upstream release.
+    # The downgrade action discards URL sources. Restore this job after the
+    # StochasticDiffEqCore branch for JumpProcesses 10 is released and the
+    # temporary pin is removed.
+    if: ${{ false }}
     uses: "SciML/.github/.github/workflows/downgrade.yml@v1"
+    with:
+      group: All
     secrets: "inherit"

.github/workflows/ThreadSafety.yml

@@ -0,0 +1,67 @@
+name: "Thread Safety Tests"
+
+on:
+  push:
+    branches:
+      - master
+    paths-ignore:
+      - 'docs/**'
+  pull_request:
+    branches:
+      - master
+    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:
+          # Restore "lts" when the temporary StochasticDiffEq pins are removed: Julia 1.10
+          # ignores [sources], and registered StochasticDiffEq packages do not admit
+          # JumpProcesses 10, so the test environment cannot resolve.
+          - "1"
+          - "pre"
+    runs-on: ubuntu-latest
+    steps:
+      - uses: actions/checkout@v7
+
+      - name: "Setup Julia ${{ matrix.version }}"
+        uses: julia-actions/setup-julia@v3
+        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@v7
+        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,178 @@
 # Breaking updates and feature summaries across releases
 
-## JumpProcesses unreleased (master branch)
+## 10.0 (Breaking)
+
+JumpProcesses 10.0 is a transitional release. StochasticDiffEq, ModelingToolkit, and
+Catalyst do not yet support it, so environments that include them keep JumpProcesses
+9.x until their updates, and JumpProcesses 10.0.1, are released. Jump problems coupled
+to delay differential equations (DDEs) are not supported in JumpProcesses 10.0 at this
+time.
+
+  - **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))
+    ```
+  - **Breaking**: `JumpProcesses.DEFAULT_RNG` has been removed. Custom jump
+    affects should draw from `SciMLBase.get_rng(integrator)`. Add SciMLBase as
+    a direct project dependency and load it with `using SciMLBase` to use this
+    interface. Pass `rng` or `seed` to `solve`/`init` to control the solver's RNG.
+  - Jump callbacks now sample from the integrator's RNG. `JumpProblem`s,
+    aggregators, and callbacks no longer store their own RNGs. Aggregators and
+    callbacks still contain mutable working state, so concurrent solves must
+    isolate that state; sharing an arbitrary `JumpProblem` across threads or
+    tasks is not made safe by this change.
+  - `SSAStepper` and the OrdinaryDiffEq ODE/DAE pathways accept `rng` and `seed`
+    through `solve`/`init`; all five CPU tau-leaping algorithms accept them
+    through `solve`. These pathways use `rng` first, then `Xoshiro(seed)` if a
+    seed is supplied, then `Random.default_rng()`.
+  - StochasticDiffEq SDE/RODE `solve` and `init` follow the backend's RNG policy.
+    An explicit RNG other than `TaskLocalRNG` takes priority over seeds. With
+    no RNG or with `TaskLocalRNG`, a nonzero solve/init seed takes priority,
+    followed by the underlying problem's nonzero seed, then a randomly seeded
+    `Xoshiro`. `TaskLocalRNG` is converted to `Xoshiro`; `seed = 0` means no seed
+    override on this pathway, whereas the JumpProcesses RNG resolver uses
+    `Xoshiro(0)`.
+  - **Breaking**: Solvers owned by JumpProcesses (`SSAStepper` and the
+    OrdinaryDiffEq ODE/DAE pathways) no longer copy a `JumpProblem`'s jump state.
+    On every thread they reuse the problem's aggregator and jump callbacks,
+    re-initializing them at each `init`; previously, solves on threads other than
+    thread 1 copied this state. A `JumpProblem` therefore supports one active solve
+    or integrator at a time. Concurrent solves from threads or tasks, and
+    integrators alive at the same time, must each use an independent copy, such as
+    `deepcopy(jprob)`. Problems created with `remake` share the original's jump
+    state, so they are not independent either. `EnsembleProblem` provides the
+    needed copies: `EnsembleThreads` copies the problem once per spawned task, and
+    `safetycopy = true` copies it for every trajectory. See the new tutorial on
+    ensembles and problem reuse for patterns that balance performance and safety.
+  - **Breaking**: The `alias_jump` keyword argument has been removed. Passing it to
+    `solve`, `init`, `JumpProblem`, or a multi-trajectory `EnsembleGPUKernel`
+    solve, or storing it on the problem wrapped by a `JumpProblem`, raises an
+    error; remove it, and solve a copy of the problem when
+    independent jump state is needed. StochasticDiffEq SDE/RODE solvers ignore the
+    keyword and keep their own policy: pass
+    `alias = SciMLBase.SDEAliasSpecifier(; alias_jumps = false)` (or
+    `SciMLBase.RODEAliasSpecifier` for RODEs) to copy jump state, or
+    `alias_jumps = true` to reuse it. When unspecified, those solvers reuse jump
+    state on thread 1 and copy it on other threads.
+  - **Breaking**: Callbacks stored on the problem wrapped by a `JumpProblem` are not
+    supported. The `JumpProblem` constructor, and `remake` with a new wrapped problem,
+    now raise an `ArgumentError` if that problem stores a `callback` (other than
+    `nothing` or an empty `CallbackSet`). Previously, OrdinaryDiffEq solvers
+    (including `FunctionMap`) ran such callbacks, while `SSAStepper` and
+    StochasticDiffEq SDE/RODE solvers silently ignored them. Pass callbacks to the
+    `JumpProblem` constructor or to `solve`/`init` instead:
+    ```julia
+    # Before (no longer works):
+    oprob = ODEProblem(f, u0, tspan, p; callback = cb)
+    jprob = JumpProblem(oprob, Direct(), jump)
+    # After:
+    oprob = ODEProblem(f, u0, tspan, p)
+    jprob = JumpProblem(oprob, Direct(), jump; callback = cb)
+    ```
+    `remake(prob; callback = nothing)` removes callbacks from an existing problem.
+  - `SSAStepper` now supports the common `alias` keyword argument for its `u0`, `p`,
+    and `tstops` inputs, given as a `Bool`, a `SciMLBase.DiscreteAliasSpecifier`, or
+    an `ODEAliasSpecifier` (which also controls `tstops`). The defaults match
+    OrdinaryDiffEq: `u0` is copied, `p` is reused, and a caller's `tstops` array is
+    never modified. `alias` never controls jump state.
+  - `SSAStepper` now saves states, and evaluates `integrator(t)`, with
+    `recursivecopy` as OrdinaryDiffEq does, so nested state types are saved as
+    independent snapshots. `tstops` given as a number, a tuple, or an
+    `AbstractVector` other than a `Vector` of the time type (such as a range or
+    view) is copied into a `Vector` at `init`, so `add_tstop!` works with it.
+  - **Breaking**: Evaluating an `SSAStepper` integrator, as `integrator(t)` or
+    `integrator(out, t)`, at a time other than `integrator.t` now raises an error,
+    unless the problem is solved with `SSAStepper(; save_uprev = true)`, which also
+    allows times in the current step, `[integrator.tprev, integrator.t]`.
+    Previously the current state was returned for any `t`, which was silently wrong
+    for times before the most recent event, for example in a `SavingCallback` or
+    `FunctionCallingCallback` with save times. Alternatively, pass such times as
+    `tstops`, so that the integrator stops exactly at them. DiffEqCallbacks'
+    integrating callbacks are not yet supported with `SSAStepper`.
+  - New `SSAStepper(; save_uprev = true)` option: the integrator keeps the state at
+    the start of each step in `integrator.uprev`, copying the state on every step,
+    and provides `get_tmp_cache` for callbacks that evaluate it in place. The
+    default, `SSAStepper()`, keeps no copy. `SSAStepper` is now a parametric type;
+    `SSAStepper()` is unchanged and has type `SSAStepper{false}`. After `solve!`,
+    `integrator.tprev` is the time of the last event or tstop.
+  - Fixed `SSAStepper` saving pending `saveat` times before the end time after its
+    final callback pass. A callback that changed the state at the end time altered the
+    values saved at those earlier times, and the saved times came out of order. They
+    are now saved before the final callback pass.
+  - Fixed `SSAStepper` erroring with a `BoundsError` when saving only the final state:
+    `save_end = true` with nothing saved earlier, for example with `save_start = false`,
+    `save_positions = (false, false)`, and no `saveat`. The final state is now saved.
+  - Fixed `init` on ODE, DAE, and `FunctionMap` integrators for `JumpProblem`s:
+    callbacks stored in the `JumpProblem` (`JumpProblem(...; callback)`) were merged
+    twice, so they ran twice per step, and were reintroduced even with
+    `merge_callbacks = false`. `init` now combines stored and call-level callbacks
+    exactly as `solve` does.
+  - Stochastic `init` of a `JumpProblem` now reaches StochasticDiffEqCore's own
+    per-algorithm `JumpProblem` initializers, so `init` followed by `solve!` matches
+    `solve`. This requires a StochasticDiffEqCore release that supports JumpProcesses 10.
+  - The state types `SSAStepper` supports are now documented: vectors and
+    `SVector`s of integers or floats, species × sites matrices for spatial
+    problems, and scalars for models without mass-action jumps. See "Supported
+    state types" in the solver documentation.
+  - **Breaking**: Multi-trajectory `EnsembleGPUKernel` solves reject `alias` values
+    other than `nothing` and the `alias_jumps` keyword, whether passed to `solve` or
+    stored on the problem, since kernels always build device-owned state. They
+    also reject `SSAStepper(; save_uprev = true)`, since kernels never expose an
+    integrator.
+  - The minimum supported SciMLBase version is now 3.34.
+  - `SSAIntegrator` now supports the `SciMLBase` RNG interface (`has_rng`,
+    `get_rng`, `set_rng!`).
+  - **Breaking**: Passing `seed` to a multi-trajectory `EnsembleGPUKernel`
+    solve with a GPU device backend now raises `ArgumentError`, instead of
+    accepting a seed that does not control the device RNG. This applies to
+    `SSAStepper`, `SimpleTauLeaping`, and `SimpleExplicitTauLeaping`. Seed the
+    device through the backend's API before solving (for example,
+    `CUDA.seed!`). Multi-trajectory kernel solves also reject host `rng`
+    objects and `rng_func`. The `CPU()` backend still supports `seed` using
+    Julia's task-local RNG; repeatability requires the same thread count,
+    workload/trajectory count, CPU backend settings, and Julia/package versions.
+    Identical trajectories across thread counts are not guaranteed. With `trajectories = 1`,
+    `EnsembleGPUKernel` retains its serial CPU fallback and that pathway's
+    usual RNG/seed rules.
+  - **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`) retain `scaled_rates === nothing`.
+    Solvers and aggregators own the working rate buffers and fill them from
+    current parameters at initialization and reset. `remake` no longer mutates
+    shared rate coefficients in the jump definition. Fixed-rate and symbolic
+    definitions retain their supported `scaled_rates` representation. 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.
+      + The `update_jump_params` keyword has been removed from
+        `reset_aggregated_jumps!`; supplying either `true` or `false` raises an
+        explanatory error. Rate refresh is automatic during reset.
+      + Custom parameter mappers (e.g. ModelingToolkitBase's
+        `JumpSysMajParamMapper`) must replace the old extraction/mutation
+        callables with `(mapper)(dest::AbstractVector, maj::MassActionJump, params)`
+        and fill `dest` with the scaled rates. The built-in mapper honors
+        `maj.rescale_rates_on_update`; custom mappers own their scaling policy
+        and must avoid scaling already-scaled symbolic expressions twice.
+        See the [custom mapper migration guide](docs/src/migration.md#custom-parameter-mappers)
+        for an example.
+  - 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.33.1
 
@@ -28,7 +200,6 @@
 
     ```julia
     using JumpProcesses
-
     rate(u, p, t) = p[1] * u[1] / (1 + u[2])
     function affect!(integrator)
         integrator.u[1] -= 1
@@ -41,10 +212,8 @@
             urate = p[1] * uhigh[1] / (1 + ulow[2])
         )
     end
-
     jump = ConstantRateJump(rate, affect!; bounds)
     prob = DiscreteProblem([10, 0], (0.0, 10.0), [1.0])
-
     # Each species affects the rate of jump 1; jump 1 changes both species.
     vartojumps_map = [[1], [1]]
     jumptovars_map = [[1, 2]]
@@ -67,7 +236,7 @@
   - Added a default aggregator selection algorithm based on the number of passed
     in jumps. i.e. the following now auto-selects an aggregator (`Direct` in this
     case):
-    
+
     ```julia
     using JumpProcesses
     rate(u, p, t) = u[1]
@@ -82,10 +251,11 @@
     `ConstantRateJump`s, and bounded `VariableRateJump`s, one no longer needs to
     specify `SSAStepper()` when calling `solve`, i.e. the following now works for
     the previous example and is equivalent to manually passing `SSAStepper()`:
-    
+
     ```julia
     sol = solve(jprob)
     ```
+
   - Plotting a solution generated with `save_positions = (false, false)` now uses
     piecewise linear plots between any saved time points specified via `saveat`
     instead (previously the plots appeared piecewise constant even though each

LICENSE.md

@@ -1,17 +1,17 @@
 The JumpDiffEq.jl package is licensed under the MIT "Expat" License:
 
 > Copyright (c) 2017: ChrisRackauckas.
-> 
+>
 > Permission is hereby granted, free of charge, to any person obtaining a copy
 > of this software and associated documentation files (the "Software"), to deal
 > in the Software without restriction, including without limitation the rights
 > to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
 > copies of the Software, and to permit persons to whom the Software is
 > furnished to do so, subject to the following conditions:
-> 
+>
 > The above copyright notice and this permission notice shall be included in all
 > copies or substantial portions of the Software.
-> 
+>
 > THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
 > IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
 > FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE

Project.toml

@@ -1,7 +1,7 @@
 name = "JumpProcesses"
 uuid = "ccbc3e58-028d-4f4c-8cd5-9ae44345cda5"
 authors = ["Chris Rackauckas <accounts@chrisrackauckas.com>"]
-version = "9.33.1"
+version = "10.0.0"
 
 [deps]
 ADTypes = "47edcb42-4c32-4615-8424-f2b9edc5f35b"
@@ -52,22 +52,35 @@ OrdinaryDiffEq = "7.0.0, 7"
 OrdinaryDiffEqCore = "4.3.0, 4"
 OrdinaryDiffEqFunctionMap = "2.0.0, 2"
 Pkg = "1"
-PoissonRandom = "0.4"
+PoissonRandom = "0.4.8"
 PrecompileTools = "1.2.1"
 Random = "1"
 RecursiveArrayTools = "4.3.0, 4"
 SafeTestsets = "0.1"
-SciMLBase = "3.18.0, 3.1"
+SciMLBase = "3.34"
 SciMLTesting = "2.8"
 SimpleNonlinearSolve = "2, 3"
 StableRNGs = "1"
 StaticArrays = "1.9.18"
 Statistics = "1"
 StochasticDiffEq = "7.0.0, 7"
+StochasticDiffEqCore = "2.2.3"
+StochasticDiffEqImplicit = "2.2.1"
+StochasticDiffEqLeaping = "2.1.0"
 SymbolicIndexingInterface = "0.3.48"
 Test = "1"
 julia = "1.10"
 
+[sources]
+# Temporary: tests use a commit of the StochasticDiffEq branch for JumpProcesses 10
+# (OrdinaryDiffEq.jl#4836), which removes the need for a JumpProcesses StochasticDiffEqCore
+# extension. Implicit and Leaping are pinned too, since their registered versions do not
+# admit JumpProcesses 10. Remove when those releases are registered. Julia 1.10 ignores
+# [sources]; see test/test_groups.toml.
+StochasticDiffEqCore = {url = "https://github.com/isaacsas/OrdinaryDiffEq.jl.git", rev = "1de403a5e59474fab9a43bf60d376391734eb65b", subdir = "lib/StochasticDiffEqCore"}
+StochasticDiffEqImplicit = {url = "https://github.com/isaacsas/OrdinaryDiffEq.jl.git", rev = "1de403a5e59474fab9a43bf60d376391734eb65b", subdir = "lib/StochasticDiffEqImplicit"}
+StochasticDiffEqLeaping = {url = "https://github.com/isaacsas/OrdinaryDiffEq.jl.git", rev = "1de403a5e59474fab9a43bf60d376391734eb65b", subdir = "lib/StochasticDiffEqLeaping"}
+
 [extras]
 Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e"
 FastBroadcast = "7034ab61-46d4-4ed7-9d0f-46aef9175898"
@@ -83,7 +96,10 @@ SciMLTesting = "09d9d899-5365-40a9-917a-5f67fddea283"
 StableRNGs = "860ef19b-820b-49d6-a774-d7a799459cd3"
 Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"
 StochasticDiffEq = "789caeaf-c7a9-5a7d-9973-96adeb23e2a0"
+StochasticDiffEqCore = "19c5a474-6cd1-4a5f-be79-46dc34e54d7f"
+StochasticDiffEqImplicit = "5080b986-4c76-4669-b5dc-373a41579d5b"
+StochasticDiffEqLeaping = "aefaaa88-39f2-4e89-b162-d61b7e8cc81b"
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
 
 [targets]
-test = ["Adapt", "FastBroadcast", "ForwardDiff", "KernelAbstractions", "LinearAlgebra", "LinearSolve", "OrdinaryDiffEq", "OrdinaryDiffEqFunctionMap", "Pkg", "SafeTestsets", "SciMLTesting", "StableRNGs", "Statistics", "StochasticDiffEq", "Test"]
+test = ["Adapt", "FastBroadcast", "ForwardDiff", "KernelAbstractions", "LinearAlgebra", "LinearSolve", "OrdinaryDiffEq", "OrdinaryDiffEqFunctionMap", "Pkg", "SafeTestsets", "SciMLTesting", "StableRNGs", "Statistics", "StochasticDiffEq", "StochasticDiffEqCore", "StochasticDiffEqImplicit", "StochasticDiffEqLeaping", "Test"]