[joshmann-inversionsemi]

Makes the additions to inversion-fbpic:

• Add mathematical expression interpreting, MathExpression (basic Python operators via ast, as well as numpy and scipy methods)
• Add a ExpressionDensityProfile which takes a MathExpression as the density profile
• Add a _ParticleBeam class for PWFA (or other) simulations, covering fbpc.lpa_utils.bunch; includes a ExpressionParticleBeam class which takes in a MathExpression
• Add FBPIC particle-beam diagnostics and beam-aware moving-window beta resolution.
• Permit laser-only and beam-only simulations; only background gas/plasma-only simulations are excluded
• Add a PWFA downramp-injection demo with post-processing.

A couple notes:
For the purposes of permitting users to use general methods to describe density and beam profiles, a limited MathExpression was preferred over an alternative native Python method + pickling scheme. Pickling would prevent the method from being readable or modifiable in its serialized form. It also brings potential security issues alla arbitrary code execution. I think we could and should extend support for this scheme, either with additional security measures or an insistent warning to the user, for the case that more complex logic is desired for the generation of these profiles. In the meantime, the MathExpression scheme should cover most practical cases.
A gap in the logic for determining a sufficient element configuration for a running simulation, where lacking a _DensityProfile is permitted, was left. In PR #37 , a simulation_length parameter was introduced. Once PR #37 is merged, this parameter is sufficient to replace the otherwise required _DensityProfile class. Presently valid simulations will still be valid whether this PR or PR #37 is merged first.
Simulations require at least one driver (a beam or a laser), but _DensityProfiles are optional if simulation_length is set.

[joshmann-inversionsemi]

@claude , please review this PR.

[claude]

Claude finished @joshmann-inversionsemi's task in 3m 2s —— View job

---

PR Review: General profiles & particle beams

• Gather context (CLAUDE.md, diff)
• Review math_expression.py
• Review particle_beam.py / density_profiles.py
• Review simulation.py changes
• Review tests & demo

Overall this is a well-structured, carefully-documented PR — the attrs validators, the .pyi stubs, the hash-payload update, and the duplicate-species-name consolidation are all solid. Two issues stand out, one of which I consider blocking.

---

🔴 1. MathExpression sandbox is bypassable → arbitrary file write & RCE

inversion_fbpic/lib/math_expression.py — The validator allows calling any public attribute of numpy / scipy / math (_is_public_module_attribute only checks that no attribute in the chain starts with _ and that the root is an allowed module). NumPy exposes filesystem and deserialization primitives that defeat the stated security goal:

• np.ndarray.tofile(np.array([...], np.uint8), <path-bytes>) → write arbitrary bytes to an arbitrary path.
• np.load(<path>, allow_pickle=True) → unpickle → arbitrary code execution.

Both pass _validate_expression because np.load / np.ndarray.tofile are "public module attributes" and the allow_pickle=True keyword is an allowed ast.keyword. Chaining them (write a pickle with tofile, then np.load(..., allow_pickle=True)) yields self-contained RCE from a config file — exactly the "arbitrary code execution" the PR description cites as the reason to avoid pickling. (I confirmed this by analysis of the validator; I could not execute a PoC here because the test sandbox blocked running Python — see note at the bottom.)

The root cause is that the allowlist is at module granularity. A robust fix is to allowlist specific callables/attributes (e.g. an explicit set like {"sin","cos","exp","where","sqrt","special.jv", ...}) rather than "any non-underscore attribute of these modules", or at minimum to maintain a denylist for the known-dangerous surface (load, loads, save, savez*, fromfile, tofile, tobytes, memmap, fromregex, genfromtxt, DataSource, ctypeslib, ndarray.view, …). A denylist is fragile; an explicit allowlist is the defensible design here. Fix this →

---

🟠 2. Laser-only / beam-only simulations pass validation but crash at setup

simulation.py now permits configs with no density profile (good), but Simulation.z_extent (line ~421) is computed solely from self.densities and returns (None, None) when empty. That value flows into:

• setup_simulation: L_interact = self.z_extent[1] - self.z_extent[0] + right_buffer → TypeError: unsupported operand type(s) for -: 'NoneType'.
• calculate_mean_group_beta (via resolve_beta_window for a laser-only sim): np.linspace(self.z_extent[0], ...) with None bounds.

So a beam-only or laser-only simulation — the headline capability of this PR — will fail once it reaches setup_simulation, unless a density also happens to be present. The assembly check and the runtime extent logic are now out of sync. z_extent should fall back to the laser (get_z_extent) and beam (p_zmin/p_zmax, or Gaussian extents) bounds when there are no densities. The new tests only exercise assembly and resolve_beta_window directly, so they don't catch this. Fix this →

---

🟡 Minor notes

• ExpressionParticleBeam / FileParticleBeam.__attrs_post_init__ don't call super().__attrs_post_init__(). ExpressionDensityProfile does. If _ParticleBeam/SerializableConfig ever gains post-init logic this silently skips it — worth adding the super() call for consistency.
• _density_function(self.dens_func) in ExpressionParticleBeam.__attrs_post_init__ is called only to validate (result discarded) and then rebuilt in add_to_simulation. Fine, but a short comment noting it's a validation call would help, since it reads as dead code.
• resolve_beta_window vs calculate_mean_group_beta: both now return the max beta. For calculate_mean_group_beta() with no args this is a behavior change (previously "first laser with a wavelength"). The test covers it, but confirm downstream callers expect the fastest rather than the first.
• Demo run_simulation.py: beam_data = list(flattop_beam[:7]) then overwrites each entry with concatenated flattop+downramp — clear, but a one-line comment on why index 7 (flattop_beam[7]) is passed through unconcatenated would aid readers.

---

⚠️ Note: I was unable to run the test suite or a security PoC — the execution sandbox in this job blocked python/pytest invocations (permission prompts I can't approve non-interactively). Findings above are from static analysis of the diff; the sandbox-escape in #1 follows directly from the validator's allow-any-public-attribute logic. If you'd like, I can push fixes for #1 and #2 — just say the word.
• IG-029-general-profiles-and-particle-beams

[joshmann-inversionsemi]

@claude , I have addressed each of your concerns. Please provide an updated review of this PR.

[claude]

Claude finished @joshmann-inversionsemi's task in 1m 15s —— View job

---

Updated PR review: general profiles and particle beams

• Gather context (latest commits bfb37ed … 6997537)
• Re-check Ig 00 init #1 MathExpression sandbox
• Re-check Ig 000 init #2 laser-only and beam-only extents
• Re-check minor notes
• Look for new issues

You fixed both blocking issues from my last review, so I no longer see anything blocking. The remaining items below are small cleanups.

---

✅ #1 MathExpression sandbox: resolved

The validator now allowlists exact callable paths, and eval only receives SimpleNamespace proxies built from that allowlist, never the real modules. I tried to find a way around it and couldn't:

Attack	Result
Access a non-allowlisted attribute of np, numpy, sp, scipy, mt or math (np.load, np.ndarray.tofile, np.ctypeslib…, scipy.io…)	Rejected: the path isn't in _ALLOWED_SYMBOLS, and the proxy doesn't have it anyway
Reach attributes through call results or subscripts (np.array([z]).tofile, z[0].x)	Rejected: _attribute_path returns None
Use dunder attributes (z.__class__, np.sin.__class__)	Rejected: the path isn't allowlisted
Call anything that isn't an exact allowlisted path (z(...), a subscripted callable)	Rejected by the _ALLOWED_CALLABLES membership check
Lambda, comprehension, walrus, */** unpacking, strings, IfExp, BoolOp	All rejected by the closed node allowlist
Change the shared _ALLOWED_NAMESPACES	Not possible: eval mode has no assignment, setattr isn't reachable, and __builtins__ is empty

The ufunc and math sets contain only pure numerical functions, and the parametrised rejection test covers the vectors I reported before. Good fix.

✅ #2 laser-only and beam-only crash: resolved by deferring the feature

Raising whenever there is no density profile (075d5cc) means z_extent can no longer return (None, None) into setup_simulation or calculate_mean_group_beta. That's reasonable until #37 lands. Three follow-ups:

• Unreachable code. The if not (self.lasers or self.beams): ... block in inversion_fbpic/lib/simulation.py:411-419 now sits after an unconditional raise, so it never runs. Either delete it and keep the TODO comment, or reorder it so it's clearly the code to restore after Ig 022 exploration campaign 1 #37.
• Missing space in the error message. The two adjacent string literals at simulation.py:404-405 and :408-409 join into "...simulation extent.Provide a density profile."
• Stale PR description. It still lists "Permit laser-only and beam-only simulations" as a bullet. Update it so reviewers and the changelog don't assume the feature ships with this PR.

Fix this →

✅ Minor notes from last time: resolved

• super().__attrs_post_init__() is now called in both beam subclasses. The base class has a no-op hook, and a test asserts the call order.
• The validation-only _density_function(...) call now has a comment.
• The "index 7 is the particle charge" comment was added to the demo.
• I take the new calculate_mean_group_beta() behaviour (max across lasers) as intentional, since test_mean_group_beta_uses_all_laser_wavelengths pins it.

---

🟡 New, low-severity (optional)

1. out= lets an expression overwrite its input arrays. All ufuncs are allowlisted and any identifier keyword is accepted. So an expression like np.multiply(z, 0, out=z) passes validation and changes the z/r arrays FBPIC passes into the density function. That is a correctness hazard rather than a security one. Rejecting the out keyword in the ast.keyword branch of _validate_expression (math_expression.py:317) closes it cheaply.
2. Unbounded expressions can hang a run. Something like 9**9**9**9 (Python big-int power) or np.arange(1e12) validates and then hangs or runs out of memory when the profile is evaluated. A config author can already exhaust resources through grid sizes, so this isn't a new trust boundary. It's worth a sentence in the MathExpression docstring in case configs ever come from less-trusted sources.
3. Missing constants. np.inf and np.nan aren't exposed. They come up in density cutoffs, for example np.where(z < z0, np.nan, …). Consider adding them next to np.e and np.pi at math_expression.py:174.

---

⚠️ Python and pytest are still blocked in this job, so this review is from reading the code and I haven't run the tests. If you'd like CI-independent verification, adding Bash(python -m pytest:*) to the action's --allowedTools would let me run the suite next time.
• IG-029-general-profiles-and-particle-beams

[joshmann-inversionsemi]

Addressed each of Claude's concerns. Ready for human review!

[joshmann-inversionsemi]

Merged the changes from #40 ; @claude , please review this merge.

[claude]

Claude finished @joshmann-inversionsemi's task in 1m 21s —— View job

---

Review of the #40 merge (7ed6906)

• Identify the merge commit and what Ig 028 fix restart iteration bounds #40 brought in
• Check the conflict resolution in simulation.py
• Check the merged test_simulation.py and README_LIB.md
• Check how the restart logic interacts with the new particle beams
• Post findings

I found nothing blocking. The merge keeps both sides intact, and the restart changes from #40 slot in correctly after the new beam-injection block.

---

What the merge changed on this branch

Comparing 610f60b (this branch before the merge) with 7ed6906, only #40's changes came across:

• simulation.py
  • Electron tracking moved from inside the density loop to after the restart attempt, guarded by not restart_loaded.
  • restart_from_checkpoint now runs before lasers are added and before the moving window is set. If no checkpoints exist, it falls back to a fresh start.
  • Lasers are only added on a fresh start.
  • run_simulation now steps num_steps - iteration + 1, clamped to 0 with a warning.
• test_simulation.py: four new unit tests, iteration=0 on the mocked backend, and the restart-equivalence integration test.
• README_LIB.md: the new section on the boosted-frame restart limitation.

None of the PR's own changes were dropped or altered:

• ExpressionDensityProfile
• beam diagnostics
• the density-required check
• resolve_beta_window
• the MathExpression changes

✅ Ordering in setup_simulation is correct

The resolved order (simulation.py:996-1150) is:

1. Densities and their diagnostics
2. Beams and their diagnostics
3. restart_from_checkpoint
4. Tracking
5. Lasers
6. Moving window
7. Checkpoints

This is the right place for the beam block:

• Species order matches. By default, FBPIC's restart_from_checkpoint maps checkpoint species to sim.ptcl by index ("species 0", "species 1", …). Beams are appended after the density species on both fresh and restarted runs, so the order is the same and beam particles are restored from the checkpoint.
• Beam diagnostics are registered before the restart. That's the same pattern Ig 028 fix restart iteration bounds #40 uses for density diagnostics.
• elecs is still defined when the tracking loop moved out of the density loop, since it's declared at :996 before the loop.

🟡 Low-severity notes (optional)

1. Beam space-charge is computed and then discarded on restart. The add_particle_bunch* calls in particle_beam.py keep FBPIC's default direct_space_charge=True. On a restart, each beam is generated and its initial fields are solved and deposited, and then restart_from_checkpoint overwrites both. The result is still correct. It's wasted work for large FileParticleBeam/openPMD beams, though, and the extra particles briefly use memory. A cheap fix is to skip generation on restart and add only an empty species. That is more invasive, so a short comment near :1065 explaining that the restart overwrites these may be enough.
2. No restart coverage for beams. Ig 028 fix restart iteration bounds #40's unit tests and the test_restarted_diagnostic_dumps_match_uninterrupted_run integration test only use [hyp, density, laser]. Adding a beam case, or a mock test that asserts beam.add_to_simulation runs before restart_from_checkpoint, would pin the species-index ordering. A future refactor that moves beams after the restart would otherwise silently leave them at their fresh initial state. Fix this →
3. track_electrons ignores beam species. This was true before the merge too. Beam particles are never tracked, even when tracking is on. For PWFA, where you often want to tag the driver or witness, it may be worth documenting that the flag covers background-plasma electrons only, or adding a per-beam track option later.

Still open from my last review (not merge-related)

These are still present on 7ed6906:

• Unreachable block: the if not (self.lasers or self.beams) block at simulation.py:411-419 comes after the unconditional raise at :407, so it never runs.
• Missing space: the string literals at simulation.py:404-405 and :408-409 still join into "...simulation extent.Provide a density profile."
• Stale PR description: the description still lists "Permit laser-only and beam-only simulations".

Fix the first two →

---

⚠️ I reviewed this by reading the code only. I didn't run pytest or the integration tests: Python/pytest were blocked in earlier runs of this job, and I didn't attempt them this time. The CI results for 7ed6906 are the authoritative check.
• IG-029-general-profiles-and-particle-beams

[joshmann-inversionsemi]

Now that #37 will be merged soon, I will want to take more time to merge those changes here; hold off on reviewing, all.

[moellep]

The PR looks good and will help simplify lume-fbpic a lot. When this PR and #37 are on main, I can rework lume-fbpic.

[cedoss-inversion]

#37 is now merged!

[joshmann-inversionsemi]

• Merged from Ig 022 exploration campaign 1 #37
• added np.clip as a MathExpression function (among others)
• fixed the long expression issue and added test coverage for it
• added a SimulationHyperparameters.use_cuda option (defaulting to True to be consistent with previous simulations)
• added a new check for the recorded git hash consistency (in case a future contributor forgets to run pre-commit install)

@claude , please review this updated PR.

[claude]

Claude finished @joshmann-inversionsemi's task in 1m 49s —— View job

---

Review of the updated PR (#37 merge, np.clip, long YAML expressions, use_cuda, git-hash test)

• Gather context (commits ccfad42 … 7af4892, diff vs origin/main)
• Review Ig 022 exploration campaign 1 #37 merge and the density-free simulation path in simulation.py
• Review MathExpression additions (np.clip etc.)
• Review long-expression YAML fix
• Review use_cuda hyperparameter
• Review git-hash consistency test
• Re-check items still open from last review
• Post findings

Most of the update looks good. One real bug remains: a laser-only simulation now passes validation but should crash in setup_simulation. Everything else is minor.

---

🟠 1. A laser-only simulation with simulation_length set should crash in resolve_beta_window

3d1e85b allows [hyp(simulation_length=...), laser] with no density. On that path, setup_simulation calls resolve_beta_window() (simulation.py:873), which calls calculate_mean_group_beta(laser.wavelength). That runs:

zs = np.linspace(self.z_extent[0], self.z_extent[1], num_samples, endpoint=False)  # simulation.py:711

With no densities, z_extent returns (None, None), so this should raise TypeError. Fixing the extent alone wouldn't help either: sum([]) at :714 gives the int 0, and len(betas) at :721 would then fail.

A laser's group velocity in vacuum is c, so calculate_mean_group_beta should return 1.0 early when self.densities is empty.

Notes on the other paths:

• Beam-only works because it only uses gamma0.
• The L_interact path is safe: simulation_length is required in this case, so z_extent isn't used at :866.
• Why the tests don't catch it: test_explicit_simulation_length_adds_right_buffer monkeypatches resolve_beta_window, and no test builds a laser-only simulation with simulation_length. A test that calls resolve_beta_window() on [hyp(simulation_length=4e-6), laser] would cover it.

Fix this →

🟡 2. Small cleanups in the new validation (simulation.py:408-433)

• Redundant inner check. The if not (self.lasers or self.beams) block at :418-426 repeats the check at :427. Its message, "Provide a density profile or at least one driver", is also misleading: a driver is always required, and a density alone is no longer enough. Deleting the inner block would leave the single, accurate check at :427.
• Stale test names. test_laser_requires_density_profile and test_particle_beam_requires_density_profile (test_simulation.py:432, :453) now actually test "no density and no simulation_length". Renaming them and adding the positive case from Ig 00 init #1 would make the new rule clear.

The missing-space and unreachable-code items from my last review are resolved, and the PR description is updated.

✅ MathExpression additions

np.clip and the other additions are pure numerical functions, and the symbols are still bound through the exact-path allowlist and SimpleNamespace proxies. I didn't find any new escape. Two small things:

• np.empty and np.empty_like return uninitialised memory, so an expression that uses them as a base gives non-deterministic profiles. Unless you need them, I'd drop them; zeros and full cover the same use.
• out= is still accepted (from my last review, optional). math_expression.py:339 accepts any identifier keyword, so np.multiply(z, 0, out=z) can still overwrite the z/r arrays FBPIC passes in. The new view-returning functions (ravel, reshape, transpose, …) can also be out= targets. Rejecting node.arg == "out" there is a one-line fix.

✅ Long-expression YAML fix

width=float("inf") in dump_yaml_with_comments (commented_yaml.py:466) fixes the root cause, which was PyYAML wrapping a plain scalar so the injected comment broke it. It also covers any other long string field, not just MathExpression. The YAML and JSON round-trip tests cover it.

✅ use_cuda

It is wired correctly into FBPICSimulation(use_cuda=...), defaults to True, and the stub and docstring are synced.

Optional: use_cuda is an execution setting, but it now goes into config_hash() through hyparams.to_dict(). The same physics run on CPU and on GPU will therefore get different hashes, so skip_if_hashed won't treat them as the same run. The recorded git hash already changes per commit, so this matters little. If you want to keep execution settings out of the hash, you could exclude use_cuda, along with use_mpi.

🟡 Git-hash consistency test (test_record_git_hash.py:20)

The test works in CI: pip install -e records the hash of the checked-out merge commit, and fetch-depth: 0 is set. Locally it has some rough edges:

• Missing file. If git_hash.txt doesn't exist (fresh clone, not yet installed), read_text raises FileNotFoundError instead of showing your helpful assertion message.
• No .git directory. In an sdist or a copied source tree, git rev-parse raises CalledProcessError. A pytest.skip would suit that case better.
• The hint only covers one cause. The message says to run pre-commit install. But README_LIB.md:466 notes that git reset isn't covered by the hooks, so a developer who already has the hooks installed can still hit this. Mentioning python tools/record_git_hash.py in the message would cover both cases.

Fix this →

Also checked

• Ig 022 exploration campaign 1 #37 merge: simulation_length and right_buffer combine correctly at :864-867.
• Beam diagnostics: still registered before restart_from_checkpoint.
• Beam restart test: the new test (4a91a1d) covers the restart species-ordering concern from my last review.
• Pre-commit hook regex: .pre-commit-config.yaml now includes math_expression and particle_beam.

---

⚠️ This job's sandbox blocked Python and pytest again, so I couldn't run the tests or reproduce #1. I traced it by reading the code, from setup_simulation → resolve_beta_window → calculate_mean_group_beta → np.linspace(None, None, …). A quick test would confirm it.
• IG-029-general-profiles-and-particle-beams

[joshmann-inversionsemi]

Addressed each of your concerns, ready for final review from you @claude .

[claude]

Claude finished @joshmann-inversionsemi's task in 49s —— View job

---

Final review of PR #41

• Gather context (commits 5d7316d … 3541632)
• Re-check the laser-only calculate_mean_group_beta fix
• Re-check the validation cleanup and test renames
• Re-check MathExpression (empty, out=)
• Re-check the git-hash test
• Look for new issues
• Post findings

You've addressed everything from my last review, and I found nothing blocking. From my side this is ready to merge. Two optional notes are at the end.

---

✅ Laser-only crash: fixed (5d7316d)

• Early return. calculate_mean_group_beta now returns 1.0 when self.densities is empty (simulation.py:685). This runs before z_extent is used, so np.linspace(None, None) and sum([]) can't be reached anymore.
• Other z_extent uses. I checked the only other caller, at :859. It only runs when simulation_length is None, and validation already requires densities in that case. So z_extent can't return (None, None) into any code path.
• Test coverage. test_laser_only_simulation_uses_vacuum_group_beta tests both calculate_mean_group_beta and resolve_beta_window on [hyp(simulation_length=4e-6), laser]. This is the test I asked for.

✅ Validation cleanup: done (95664d1, 3c4a141)

• Redundant check. The duplicate inner check is gone. The two remaining checks are independent and their messages are accurate:
  • "a density or simulation_length"
  • "at least one driver"
• Tests. The renamed tests describe what they actually test now. There are new positive tests for a laser-only config and for a beam-only config.

✅ MathExpression: done (ee3a04e, 3da662f)

• empty and empty_like are removed from the allowlist.
• out= is rejected. The check runs in the ast.keyword branch, so it covers every call whether the callable is a ufunc or not. test_function_rejects_out_keyword_to_preserve_input_arrays covers it. "outer" is a function name, not a keyword, so the check doesn't affect it.

✅ Git-hash test: done (3541632)

• No git or not a checkout. It now skips when the git binary is missing (FileNotFoundError) or when rev-parse fails.
• Missing file. If git_hash.txt doesn't exist, the test fails with an assertion message instead of raising FileNotFoundError.
• Hints. Both failure messages now mention python tools/record_git_hash.py as well as pre-commit install.
• pytest is imported at the top of the file.

---

🟡 Optional, non-blocking

1. Stale beta_window docstring. The docstring for SimulationHyperparameters.beta_window (simulation.py:151-153, and the copies in simulation.pyi:209 and :295) still describes the old rule: "the first LaserPulse object with a wavelength… If no wavelength… 1.0". The actual rule is in resolve_beta_window:

   • it takes the largest beta across all lasers and beams (sqrt(1 - gamma0**-2) for a beam);
   • it uses vacuum (1.0) for lasers when there's no density.

   Syncing the docstring with resolve_beta_window would keep the generated docs accurate.

2. Beam particles can't be selected in particle diagnostics. _selected_recordings in diagnostics.py:212-237 only builds its targets from densities. A beam's particle_name is written to disk, but load_particles can't select it: you get "Unknown recorded … target". For PWFA, people will want to analyse the driver and witness, so this is a natural follow-up. It doesn't need to block this PR.

---

⚠️ I reviewed this by reading the code only. Earlier runs of this job blocked Python and pytest, so I didn't run the tests. The CI results for 3541632 are the authoritative check.
• IG-029-general-profiles-and-particle-beams

[joshmann-inversionsemi]

Updated the beta_window docstring and added support for applying (post-simulation) diagnostic analyses to particle beams.

Ready for another review, @moellep !

[moellep]

Looks good - thanks for the changes!