Negative distances at medium transitions with activated magnetic field
Depending on the energies and cuts you choose, this is a very common issue that has been observed by several people propagating electromagnetic showers:
In some cases, the propagation fails with an error of the type:
terminate called after throwing an instance of 'boost::wrapexcept<boost::math::evaluation_error>' what(): Error in function boost::math::tools::bisect<double>: No change of sign in boost::math::tools::bisect, either there is no root to find, or there are multiple roots in the interval (f(min) = -5.8245694739599418).
The reason for this error is that
corsika::proposal::ContinuousProcess::do_continuous() receives a negative distance to propagate, which PROPOSAL can not handle.
If I activate
trace logging level, the last output is
Tracking pos: (-833.118 -5.23421e-19 6.411e+06) m Tracking p: (-0.010693 -9.16556e-28 0.00439353)  GeV Tracking v: (-2.77027e+08 -2.37456e-17 1.13825e+08) m/s [corsika:debug (TrackingLeapFrogCurved.inl:110)] gyroradius 771.227 m, steplimit: 15.4235 m = 5.14975e-08 s [corsika:debug (Intersect.inl:29)] volumeNode=0x555556061ee0, numericallyInside=true [corsika:debug (Intersect.inl:47)] intersection times with currentLogicalVolumeNode: -3.444667695403443e-05 s and 3.2353023108619646e-05 s [corsika:debug (Cascade.inl:138)] transport particle by : -8.150316578976387e-06 m Medium transition after: -8.150316578976387e-06 m Decay after: inf m Interaction after: inf m Continuous limit: 1514.897747410968 m
One can see that the negative distance is produced when a particle makes a transition to a new medium. The negative distance is very small, so this looks like some numerical instability to me.
I've been running the
em_shower example on the current master with an energy cut of
5 MeV and an initial energy of
1 TeV. If the magnetic field is enabled, 4 out of 50 showers fail with this error (you can reproduce a failing shower by setting the seed to
If I set the magnetic field to zero, all showers pass. Therefore, I think this issue might have to do with TrackingLeapFrog in combination with the calculation of intersections.