Skip to content

ENH: reproducible Monte Carlo via per-simulation-index seeding#1054

Open
thc1006 wants to merge 3 commits into
RocketPy-Team:developfrom
thc1006:enh/reproducible-montecarlo-seeding
Open

ENH: reproducible Monte Carlo via per-simulation-index seeding#1054
thc1006 wants to merge 3 commits into
RocketPy-Team:developfrom
thc1006:enh/reproducible-montecarlo-seeding

Conversation

@thc1006

@thc1006 thc1006 commented Jul 7, 2026

Copy link
Copy Markdown

Pull request type

  • Code changes (bugfix, features)

Current behavior

MonteCarlo.simulate() seeds the stochastic models per worker in parallel mode (from a fresh, unseeded np.random.SeedSequence().spawn(n_workers)) and once at construction in serial mode. So the sampled inputs depend on the execution mode and the number of workers, and parallel runs aren't reproducible run to run. This is #1053.

New behavior

Adds a keyword-only random_seed to simulate(). From that root I spawn one child seed per simulation index and reseed the models from child_seeds[i] before simulation i. SeedSequence.spawn is prefix-stable, so index i maps to the same seed regardless of which worker runs it. The sampled inputs come out identical across serial, parallel(2) and parallel(N), and reproducible from the seed.

random_seed is a seed, not a live RNG: it takes an int, a sequence of ints, or a SeedSequence, with None = fresh entropy so existing behavior is unchanged unless you pass a seed. A supplied SeedSequence is copied from its full state before spawning, so it is never mutated and repeated calls with the same object reproduce the same run. This is informed by SPEC 7 (a per-call keyword-only seed) and NumPy's parallel idiom (SeedSequence(root).spawn(n)), but it keeps seed-snapshot semantics rather than SPEC 7's stateful rng: a Generator/BitGenerator is not accepted, because reducing it to its underlying SeedSequence would ignore how far it has been consumed. Pass rng.bit_generator.seed_seq to seed from an existing generator.

Tests

tests/unit/simulation/test_monte_carlo_determinism.py unit-tests the seed handling (the accepted seed types, the SeedSequence copy leaving the caller untouched, the three-way env/rocket/flight split) and the parallel index claim (a deterministic barrier-based test that fails if the claim's lock is dropped). tests/integration/simulation/test_monte_carlo_determinism.py covers end-to-end reproducibility (serial, and serial == parallel(2) == parallel(4)) with a stubbed Flight.

Notes from review

  • Parallel workers claimed the next index with an unlocked keep_simulating() + increment(). Near the end of a run two workers could both pass count < n and then claim sim_idx == n, and the per-index child_seeds lookup turned that into an IndexError. The claim now holds the shared mutex across the check and the increment, so each index is handed out once.
  • A supplied SeedSequence was returned as-is, and spawn() advances its child counter, so passing the same object twice was not reproducible. It is now copied from its full state, and Generator/BitGenerator are no longer accepted (see above).

Known limitation

From #1053: list-valued stochastic attributes are sampled with the stdlib random.choice (an unseeded global), so a model with a multi-element thrust_source isn't fully seed-controlled yet. A full fix would seed the stdlib random per index too.

Breaking change

  • Yes

The exact numbers a run produces change (per-index seeding, the env/rocket/flight decorrelation, and the serial index now counting from 0 to match parallel), so external code that pinned exact Monte Carlo samples would need to re-baseline. The in-repo Monte Carlo tests don't pin exact values (test_monte_carlo_simulate checks apogee and impact velocity within a tolerance), so none of them change. There's no API break for users who don't pin exact samples, and random_seed is opt-in.

One more heads-up: Python 3.14 flipped the Linux multiprocessing default from fork to forkserver, so workers no longer inherit memory and the seed objects crossing the process boundary must be picklable. SeedSequence is, and the tests pass under fork.

Closes #1053

@thc1006 thc1006 marked this pull request as ready for review July 8, 2026 19:35
@thc1006 thc1006 requested a review from a team as a code owner July 8, 2026 19:35
Copilot AI review requested due to automatic review settings July 8, 2026 19:35

Copilot AI left a comment

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Copilot was unable to review this pull request because the user who requested the review has reached their quota limit.

@codecov

codecov Bot commented Jul 8, 2026

Copy link
Copy Markdown

Codecov Report

❌ Patch coverage is 85.18519% with 4 lines in your changes missing coverage. Please review.
✅ Project coverage is 82.30%. Comparing base (4a67a41) to head (761c092).
⚠️ Report is 4 commits behind head on develop.

Files with missing lines Patch % Lines
rocketpy/simulation/monte_carlo.py 85.18% 4 Missing ⚠️
Additional details and impacted files
@@             Coverage Diff             @@
##           develop    #1054      +/-   ##
===========================================
+ Coverage    81.74%   82.30%   +0.56%     
===========================================
  Files          118      118              
  Lines        15248    15262      +14     
===========================================
+ Hits         12464    12561      +97     
+ Misses        2784     2701      -83     

☔ View full report in Codecov by Harness.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:
  • ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

@thc1006

thc1006 commented Jul 8, 2026

Copy link
Copy Markdown
Author

Pushed 0d37ed6 to cover the seeding code that the earlier tests only reached through a full (slow) Monte Carlo run.

  • New unit tests in tests/unit/simulation/test_monte_carlo_determinism.py drive __root_seed_sequence and __seed_simulation directly: every supported random_seed type normalizes to the same root stream, None draws fresh entropy, an existing SeedSequence/Generator/BitGenerator is reused rather than copied, and each child seed splits three ways so environment, rocket and flight get independent streams.
  • The end-to-end simulate reproducibility tests moved to tests/integration/ next to test_monte_carlo_simulate. The serial run is in the normal suite now; only the fork-based worker-invariance test stays slow, and it imports multiprocess lazily the way the library does.

That takes the patch from 21 uncovered lines down to 4, all in the parallel branch. I looked into whether it is worth covering those too, and I do not think it is:

  • Three of them (__run_in_parallel's dispatch, the child-seed spawn, and the worker loop) run in the main process. The slow worker-invariance test already exercises them; the coverage jobs just do not pass --runslow.
  • The fourth, the __seed_simulation call inside __sim_producer, only runs inside a forked worker, and coverage.py does not instrument subprocesses unless the project sets concurrency = multiprocessing. So no test can count that line without a repo-wide coverage change.

So un-slowing the parallel test would recover at most the three main-process lines and still not the subprocess one, at the cost of forking workers in the default suite. I kept it slow to match the other Monte Carlo multiprocessing tests. Glad to set up multiprocessing coverage separately if you want the parallel path counted, but that felt out of scope for this PR.

@Gui-FernandesBR Gui-FernandesBR force-pushed the enh/reproducible-montecarlo-seeding branch from 0d37ed6 to 761c092 Compare July 9, 2026 21:08

@phmbressan phmbressan left a comment

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The implementation is very clear and throughout, nice work.

The explanation on the concepts behind per index seeding (both in the issue and PR description) were rather helpful. I agree having reproducible results was an issue with the parallel per worker seeding.

Regarding the decisions on parameter naming, I agree with most of the decisions taken here. Moreover, the rng attribute is well docstringed, so it shouldn't be a matter of confusion to the user.

@MateusStano could you give your two cents on the changes here before we proceed with a merge?

rocket=rocket,
flight=flight,
)
montecarlo.simulate(**simulate_kwargs)

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Could you double check if after running this test, no files are left dangling on the test file system. If that is the case, using a try-finally block to remove those (as done in a few other MonteCarlo tests) is an option.

Copy link
Copy Markdown
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks for the review. On the dangling files: I checked, and nothing gets left behind. The test passes filename=tmp_path / tag into MonteCarlo, so the .inputs.txt, .outputs.txt and .errors.txt all land under pytest's per-test tmp_path and get cleaned up with it. I reran it and confirmed the repo and working directory stay clean, so there's nothing here for a try/finally to remove.

The other MonteCarlo tests need that cleanup because their fixture uses a fixed filename="monte_carlo_test", which writes into the cwd. This one keeps everything inside tmp_path instead. I'm happy to match the try/finally style if you'd prefer consistency across the file, but the tmp_path route seemed cleaner to me. Your call.

Comment thread rocketpy/simulation/monte_carlo.py Outdated
Comment on lines 467 to 474
while sim_monitor.keep_simulating():
sim_idx = sim_monitor.increment() - 1
inputs_json, outputs_json = "", ""

self.__seed_simulation(child_seeds[sim_idx])
flight = self.__run_single_simulation()
inputs_json = self.__evaluate_flight_inputs(sim_idx)
outputs_json = self.__evaluate_flight_outputs(flight, sim_idx)

@MateusStano MateusStano Jul 10, 2026

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I am not sure about this so @thc1006 could you verify, please?

But I believe adding this child_seeds[sim_idx] can break at the very end of the simulation. If you are running 8 simulations, and we are at step 7. We could have one worker check 7<8 -> yes, and a second worker also check 7<8 -> yes, because keep_simulating() and increment() are two separate steps and nothing locks them together. Then both workers increment: the first gets sim_idx = 7 (ok), but the second gets sim_idx = 8. Since child_seeds only has 8 item, child_seeds[8] does not exist

I think this race was already there before, but it didnt crash the simulation. Could you check and fix, please?

Copy link
Copy Markdown
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You're right, this is a real race. keep_simulating() and increment() are separate manager calls, so near the end multiple workers can pass the count < n check before any of them increments and then claim sim_idx == n, which the new child_seeds[sim_idx] lookup turned into an IndexError (before, it just wrote an extra record). Fixed in 9c020b6: the claim now goes through a small _claim_next_index helper that holds the shared mutex across the check and the increment, so each index is handed out once and the counter never runs past number_of_simulations. There's a deterministic unit test for it (a barrier plus a widened check-to-increment window) that over-claims and fails if the lock is removed. Thanks for catching it.

Comment on lines +298 to +310
def __seed_simulation(self, child_seed):
"""Reseed the stochastic models for a single simulation index.

The per-index child seed is split three ways so the environment,
rocket and flight draw from independent streams instead of sharing
one. Seeding per simulation index (not per worker) is what makes the
sampled inputs invariant to the execution mode and to the number of
workers.
"""
env_seed, rocket_seed, flight_seed = child_seed.spawn(3)
self.environment._set_stochastic(env_seed)
self.rocket._set_stochastic(rocket_seed)
self.flight._set_stochastic(flight_seed)

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This .spawn() call seems to make the case of using SeedSequence/Generator not repeatable

From what I understand if a user builds a SeedSequence (or Generator) once and passes the same object to simulate() twice, the first will spawns a set of children and the second run will spawn another set

Copy link
Copy Markdown
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Agreed, the current code mutated the caller's seed lineage. SeedSequence.spawn() advances n_children_spawned, and the helper returned the original object, so a second simulate() with the same object spawned a different set. Fixed in 9c020b6 with two changes:

  1. A supplied SeedSequence is now copied from its full state (entropy, spawn_key, pool_size, n_children_spawned) before spawning, so the caller's object is untouched and repeated calls reproduce the same children. Copying only entropy would drop spawn_key, which matters for a spawned child.
  2. I dropped Generator/BitGenerator from the accepted types. A generator is a stateful RNG, not a seed: its bit_generator.seed_seq is only the original lineage, so a generator advanced by a million draws would still seed the run identically to a fresh one. SPEC 7's rng semantics would instead consume it, the opposite of a reproducible seed. random_seed now takes an int, a sequence of ints, or a SeedSequence; a generator user can pass rng.bit_generator.seed_seq.

Tests cover the same-object reproducibility, n_children_spawned staying put, a spawned child with a non-empty spawn_key, and the generator rejection.

thc1006 added a commit to thc1006/RocketPy that referenced this pull request Jul 11, 2026
Addresses review feedback on RocketPy-Team#1054.

Parallel workers claimed the next index with an unlocked
keep_simulating() + increment(), so near the end of a run two workers
could both pass the count < n check and then claim sim_idx == n; the
per-index child_seeds lookup turned that into an IndexError (before, it
only wrote one extra record). Move the claim into a _claim_next_index
helper that holds the shared mutex across the check and the increment,
so each index is handed out once and the counter never overshoots. A
deterministic unit test (a barrier plus a widened check-to-increment
window) over-claims and fails if the lock is dropped.

__root_seed_sequence returned the caller's SeedSequence, and spawn()
advances its child counter, so passing the same object to simulate()
twice produced different children. Copy it from its full state instead,
which leaves the caller untouched and keeps repeated calls reproducible.
Also drop Generator/BitGenerator from the accepted types: a stateful
generator is not a seed, and reducing it to its underlying SeedSequence
ignores how far it has been consumed. random_seed now takes an int, a
sequence of ints, or a SeedSequence.

Signed-off-by: thc1006 <84045975+thc1006@users.noreply.github.com>
@thc1006 thc1006 requested a review from MateusStano July 11, 2026 01:18
thc1006 added 3 commits July 11, 2026 14:36
MonteCarlo seeded the stochastic models per worker in parallel mode (from a
fresh, unseeded SeedSequence) and once at construction in serial mode, so the
sampled inputs depended on the execution mode and the worker count, and parallel
runs were not reproducible run to run.

Add a keyword-only random_seed to simulate() (SPEC 7 style: accepts an int, a
SeedSequence, or a Generator; None keeps the previous fresh-entropy behavior).
Spawn one child seed per simulation index from that root and reseed the
stochastic models from child_seeds[i] before simulation i. SeedSequence.spawn is
prefix-stable, so index i maps to the same seed regardless of which worker runs
it, making the inputs identical across serial, parallel(2) and parallel(N). Each
index seed is split three ways so the environment, rocket and flight draw from
independent streams rather than sharing one.

The serial index field now counts from 0 to match the parallel path. Both changes
alter the numbers a fixed seed produces, so stored baselines regenerate.

Adds tests/unit/simulation/test_monte_carlo_determinism.py: serial
reproducibility, worker invariance (serial == parallel(2) == parallel(4)), and
the None-seed path.

Signed-off-by: thc1006 <84045975+thc1006@users.noreply.github.com>
The reproducible-seeding change added __root_seed_sequence and
__seed_simulation plus the per-index serial and parallel seeding, but the
only tests that reached them ran a full Monte Carlo and were marked slow,
so the coverage job (which does not pass --runslow) never executed them.

Add fast unit tests that drive the two helpers directly: every supported
random_seed type normalizes to the same root stream, None draws fresh
entropy, existing SeedSequence/Generator/BitGenerator objects are reused
rather than copied, and each child seed splits three ways so environment,
rocket and flight get independent streams.

Move the end-to-end simulate reproducibility tests into tests/integration,
next to the existing Monte Carlo simulate test. The serial reproducibility
run now lives in the non-slow suite; only the fork-based worker-invariance
test stays slow, and it imports multiprocess lazily like the library does.

Signed-off-by: thc1006 <84045975+thc1006@users.noreply.github.com>
Addresses review feedback on RocketPy-Team#1054.

Parallel workers claimed the next index with an unlocked
keep_simulating() + increment(), so near the end of a run two workers
could both pass the count < n check and then claim sim_idx == n; the
per-index child_seeds lookup turned that into an IndexError (before, it
only wrote one extra record). Move the claim into a _claim_next_index
helper that holds the shared mutex across the check and the increment,
so each index is handed out once and the counter never overshoots. A
deterministic unit test (a barrier plus a widened check-to-increment
window) over-claims and fails if the lock is dropped.

__root_seed_sequence returned the caller's SeedSequence, and spawn()
advances its child counter, so passing the same object to simulate()
twice produced different children. Copy it from its full state instead,
which leaves the caller untouched and keeps repeated calls reproducible.
Also drop Generator/BitGenerator from the accepted types: a stateful
generator is not a seed, and reducing it to its underlying SeedSequence
ignores how far it has been consumed. random_seed now takes an int, a
sequence of ints, or a SeedSequence.

Signed-off-by: thc1006 <84045975+thc1006@users.noreply.github.com>
@thc1006 thc1006 force-pushed the enh/reproducible-montecarlo-seeding branch from 9c020b6 to 3e22729 Compare July 11, 2026 06:38
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

4 participants