Step limits

Currently there are only a few, hard-coded step limitation mechanics that apply to a particle at each step. The ordering of these can affect the performance.

Discrete interactions

Most physics processes use pre-calculated cross sections that are tabulated and interpolated: see Numeric grids. Cross sections for each process are evaluated at the beginning of the step along with range limiters.

inline StepLimit celeritas::calc_physics_step_limit(MaterialTrackView const &material, ParticleTrackView const &particle, PhysicsTrackView &physics, PhysicsStepView &pstep)

Calculate physics step limits based on cross sections and range limiters.

Processes with integral cross section rejection calculate an estimated “maximum” over the step: see PhysicsTrackView.calc_max_xs .

Todo:

For particles with decay, macro XS calculation will incorporate decay probability, dividing decay constant by speed to become 1/len to compete with interactions.

Todo:

For neutral particles that haven’t changed material since the last step, we can reuse the previously calculated cross section (but may not want to if it’s expensive).

If undergoing an interaction, the process is sampled from the stored beginning-of-step cross sections.

template<class Engine>
inline ParticleProcessId celeritas::find_ppid(MaterialView const &material, ParticleTrackView const &particle, PhysicsTrackView const &physics, PhysicsStepView &pstep, Engine &rng)

Sample the process for the discrete interaction.

  • If the particle is at rest and has an at-rest process, that is returned.

  • Sample a process by selecting from probabilities proportional to the cross sections calculated at the beginning of the step.

  • If the particle has changed energy over the step by continuous energy loss, and the process supports integral cross section rejection, then the cross section is recalculated with the new energy. The collision is sampled based by rejection sampling with the updated cross section.

    Todo:

    Support competing at-rest processes?

Returns:

Sampled process, or null if integral xs rejection cancels interaction

template<class Engine>
ActionId celeritas::select_discrete_interaction(MaterialView const &material, ParticleTrackView const &particle, PhysicsTrackView const &physics, PhysicsStepView &pstep, Engine &rng)

Choose the physics model for a track’s pending interaction.

  • Sample from the previously calculated per-process cross section/decay to determine the interacting process ID.

  • From the process ID and (post-slowing-down) particle energy, we obtain the applicable model ID.

  • For energy loss (continuous-discrete) processes, the post-step energy will be different from the pre-step energy, so the assumption that the cross section is constant along the step is no longer valid. Use the “integral

    approach” to sample the discrete interaction from the correct probability distribution (section 7.4 of the Geant4 Physics Reference release 10.6).

User-input

See celeritas::inp::TrackingLimits.

Range and energy loss

See Continuous slowing down.

Multiple scattering

See Coulomb scattering.

class UrbanMscSafetyStepLimit

Sample a step limit for the Urban MSC model using the “safety” algorithm.

This distribution is to be used for tracks that have non-negligible steps and are near the boundary. Otherwise, no displacement or step limiting is needed.

Note

This code performs the same method as in ComputeTruePathLengthLimit of G4UrbanMscModel, as documented in section 8.1.6 of [The Geant4 Collaboration, 2023] or Urbán [2006] .

class UrbanMscMinimalStepLimit

Sample a step limit for the Urban MSC model using the “minimal” algorithm.

Todo:

Here and UrbanMscSafetyStepLimit should simply calculate limit and limit_min; the caller should skip sampling if not MSC limited (or below min sampling step), and another helper (documented the hardcoded 0.1 sigma width) does the gaussian sampling.

Note

This code performs the same method as in ComputeTruePathLengthLimit of G4UrbanMscModel, as documented in section 8.1.6 of [The Geant4 Collaboration, 2023] or Urbán [2006] .

When multiple scattering is in play, the “physical” step limit must be converted to and from a “geometry” step limit before and after the propagation (movement within geometry) step limiting calculation.

class MscStepFromGeo

Convert the “geometrical” path traveled to the true path with MSC.

The “true” step length is the physical path length taken along the geometrical path length (which is straight without magnetic fields or curved with them), accounting for the extra distance taken between along-step elastic collisions.

The transformation can be written as

\[ t(z) = \langle t \rangle = -\lambda_{1} \log(1 - \frac{z}{\lambda_{1}}) \]
or \( t(z) = \frac{1}{\alpha} [ 1 - (1-\alpha w z)^{1/w}] \) if the geom path is small, where \( w = 1 + \frac{1}{\alpha \lambda_{10}}\).

Param true_path:

the proposed step before transportation.

Param gstep:

the proposed step after transportation.

Param alpha:

linearity of MSC cross section over the step

class MscStepToGeo

Convert the “true” path traveled to a geometrical approximation.

This takes the physical step limit&#8212;whether limited by range, physics interaction, or an MSC step limiter&#8212;and converts it to a “geometrical” path length, which is a smooth curve (straight if in the absence of a magnetic field).

The mean value of the geometrical path length \( z \) (the first moment) corresponding to a given true path length \( t \) is given by

\[ \langle z \rangle = \lambda_{1} [ 1 - \exp({-\frac{t}{\lambda_{1}}})] \,, \]
where \(\lambda_{1}\) is the first transport mean free path.

Since \(\lambda_{1}\) depends on the kinetic energy of the path and decreases along the step, the path length correction can be approximated as a linear function over the step:

\[ \lambda_{1} (E(t')) \approx \lambda_{1,i} (1 - \alpha t') \,, 0 \le t' \le t \]
using the tabulated MFP with the start-of-step energy \( \lambda_{1,i} \equiv \lambda_{1}(E_0)\),. In the usual case, the linearity coefficient \( \alpha\) is calculated using the calculated mean free paths from the energy at the beginning and estimated step endpoint energy:
\[ \alpha \equiv \frac{\lambda_{1,i} - \lambda_{1,f}}{t \lambda_{1,i}} \,, \]
where \( \lambda_{1,f} \equiv \lambda_{1}(E(t)) \). The following edge cases calculate alpha differently:
  • For very small steps (less than celeritas::UrbanMscParameters::min_step_transform), the cross section is assumed constant due to the small energy loss: \( \alpha = 0 \) . The step is not transformed.

  • For small steps relative to the range, the cross section is approximated as constant, but the step is transformed analytically.

  • For range-limited steps ( \( t = r = \mathrm{range} \)), \( \alpha = 1/r \) by definition since the mean free path at the end point of the step energy) is zero.

  • For other low-energy steps (kinetic energy less than mass energy), \( \alpha = 1/r \) rather than using the actual endpoint cross section. This is equivalent to approximating the endpoint MFP somewhere between the initial MFP (as with the small-step case) and zero (the range-limited case).

Since the MSC cross section should increase with decreasing energy, \( \lambda_{1,f} \) will be smaller than \( \lambda_{1,i} \) and \( \alpha \) will be positive. However, in the Geant4 Urban MSC model, different methods are used to calculate the cross section above and below 10 MeV. In the higher energy region the cross sections are identical for electrons and positrons, resulting in a discontinuity in the positron cross section at 10 MeV. This means on fine energy grids it’s possible for the cross section to be increasing with energy just above the 10 MeV threshold and therefore for \( \alpha \) is negative.

The resulting geometrical step length can be greater than 1 MFP: it’s the MSC step limiter’s job to add additional restrictions.

Note

This performs the same method as in ComputeGeomPathLength of G4UrbanMscModel of the Geant4 10.7 release.

Geometry

See Propagation.