Problem 2063

Summary: Transportation is prematurely killing charged particles expected to travel a long distance in the field
Product: Geant4 Reporter: Krzysztof Genser <resnegfk>
Component: processes/transportationAssignee: John Apostolakis <John.Apostolakis>
Status: RESOLVED FIXED    
Severity: major CC: Gabriele.Cosmo
Priority: P4    
Version: 10.4   
Hardware: All   
OS: All   

Description Krzysztof Genser 2018-06-04 05:45:51 CEST
It looks like the Transportation process is killing high energy charged particles
with a long expected physical step in the field when the default number of
integration steps reach its limit.  

The Transportation then sets the track status to StopAndKill and the
step status as Geom Limit. The status is somewhat misleading and, under normal
verbosity and non debug conditions, there are no warnings indicating
that transportation is killing a particle due to the number of integration steps limit.

This behavior can cause e.g. particles traveling in a field in a vacuum to be
killed before they decay or reach a volume boundary.

The problem is often masked by short physical steps, small volumes or
by a step limiter.

Below is an example when 200 MeV/c momentum (or 104MeV KE) pi- 
in a 3.5Tesla uniform magnetic field is killed after traveling ~283mmm

The shown printout excerpts were obtained with a high verbosity levels and with the toolkit 
compiled with many of the G4DEBUG_... flags.

Stepper precision parameters were:
g4epsilonMin        5e-05
g4epsilonMax        0.001
g4DeltaOneStep      0.05 mm
g4DeltaIntersection 0.05 mm
g4DeltaChord        0.05 mm
The box volume half length was ~50m and the material was G4_Galactic

Krzysztof

*********************************************************************************************************
* G4Track Information:   Particle = pi-,   Track ID = 1,   Parent ID = 0
*********************************************************************************************************

Step#    X(mm)    Y(mm)    Z(mm) KinE(MeV)  dE(MeV) StepLeng TrackLeng  NextVolume ProcName
    0        0        0        0       104        0        0         0 BoxInTheWorld initStep
...
 >>DefinePhysicalStepLength (List of proposed StepLengths): 
    ++ProposedStep(PostStep ) =    100000 : ProcName = StepLimiter (No ForceCondition)
    ++ProposedStep(PostStep ) = 1.1887e+27 : ProcName = pi-Inelastic (No ForceCondition)
    ++ProposedStep(PostStep ) = 1.03063e+27 : ProcName = hadElastic (No ForceCondition)
    ++ProposedStep(PostStep ) =   1720.59 : ProcName = Decay (No ForceCondition)
    ++ProposedStep(PostStep ) = 2.70417e+30 : ProcName = CoulombScat (No ForceCondition)
    ++ProposedStep(PostStep ) = 1.79769e+308 : ProcName = hPairProd (No ForceCondition)
    ++ProposedStep(PostStep ) = 1.79769e+308 : ProcName = hBrems (No ForceCondition)
    ++ProposedStep(PostStep ) = 1.01185e+24 : ProcName = hIoni (No ForceCondition)
    ++ProposedStep(PostStep ) = 1.79769e+308 : ProcName = Transportation (Forced)
    ++ProposedStep(AlongStep) = 4.17169e+26 : ProcName = hIoni (CandidateForSelection)
    ++ProposedStep(AlongStep) =   1720.59 : ProcName = msc (NotCandidateForSelection)
...
  Difficult track - taking many sub steps.
       Current Position  and  Direction 
Step#       s        X(mm)      Y(mm)      Z(mm)    N_x     N_y     N_z   Delta|N|   StepLen  StartSafety   PhsStep         NextVolume 
Step taken was 0.282814 out of PhysicalStep = 1720.59
Final safety is: 49999.9
Chord length = 0.257601
...
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomNav1002
      issued by : G4PropagatorInField::ComputeStep()
  Killing looping particle  after 1000 field substeps  totaling 282.776 mm  in *volume* BoxInTheWorld
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

    ++ProposedStep(AlongStep) =   282.776 : ProcName = Transportation (CandidateForSelection)
 G4Transportation is killing track that is looping or stuck 
   This track has 104.315 MeV energy.
   Number of trials = 0   No of calls to AlongStepDoIt = 1
...
    ++G4Track Information 
      -----------------------------------------------
        G4Track Information  
      -----------------------------------------------
        Step number         :                    1
        Position - x (mm)   :                    0
        Position - y (mm)   :              -0.0107
        Position - z (mm)   :               0.0632
        Global Time (ns)    :                 1.15
        Local Time (ns)     :                 1.15
        Momentum Direct - x :                    0
        Momentum Direct - y :               -0.332
        Momentum Direct - z :                0.943
        Kinetic Energy (MeV):                  104
        Polarization - x    :                    0
        Polarization - y    :                    0
        Polarization - z    :                    0
        Track Length        :                  283
        Track ID #          :                    1
        Parent Track ID #   :                    0
        Next Volume         :        BoxInTheWorld 
        Track Status        :          StopAndKill
        Vertex - x (mm)     :                    0
        Vertex - y (mm)     :                    0
        Vertex - z (mm)     :                    0
        Vertex - Px (MomDir):                    0
        Vertex - Py (MomDir):                    0
        Vertex - Pz (MomDir):                    1
        Vertex - KineE (MeV):                  104
        Creator Process     :      Event Generator
      ----------------------------------------------

    1        0  -0.0107   0.0632       104 9.35e-24      283       283 BoxInTheWorld Transportation
Comment 1 John Apostolakis 2018-06-07 16:40:22 CEST
Are you aware that there are already parameter in G4Transportation which govern the values of 
  - the energy above which a user is warned when a particle is abandoned ( ThresholdWarningEnergy ) 
  - the energy above which particles are not abandoned at all ( ThresholdImportantEnergy )
  - the number of trials before a particle is abandoned.

The default values for these parameters are likely not appropriate for your application - but they can be tuned as desired.

Here is an excerpt from G4Transportation.hh with the relevant methods:

     inline G4double GetThresholdWarningEnergy() const; 
     inline G4double GetThresholdImportantEnergy() const; 
     inline G4int GetThresholdTrials() const; 

     inline void SetThresholdWarningEnergy( G4double newEnWarn ); 
     inline void SetThresholdImportantEnergy( G4double newEnImp ); 
     inline void SetThresholdTrials(G4int newMaxTrials ); 

     // Get/Set parameters for killing loopers: 
     //   Above 'important' energy a 'looping' particle in field will 
     //   *NOT* be abandoned, except after fThresholdTrials attempts.
     // Below Warning energy, no verbosity for looping particles is issued

The current ‘default’ values are very high:
-     fThreshold_Warning_Energy( 100 * MeV ),  
-     fThreshold_Important_Energy( 250 * MeV ), 
I will try to reduce them significantly:
+     fThreshold_Warning_Energy( 0.0 ),
+     fThreshold_Important_Energy( 1.0 * MeV ), 

I will also improve the current ‘drab’ warning, adding information about the type of particle and position, and also providing some advice about what to do to avoid the problem.

Please let me know if there are any additional things that would remain after using these (existing) method - and after improving the printout.
Comment 2 John Apostolakis 2018-06-07 16:45:20 CEST
I note that some sort of protection of this type is necessary in case a particle acquires a momentum with a very small component in the direction of a strong magnetic field (in a vacuum).

Otherwise the particle could 'drift' forever in computational time and go a very distances.

Even using a helical stepper *may* not help in those cases if the track is near a boundary - which could force smaller steps in order to ensure that a track does not intersect a volume boundary.
Comment 3 John Apostolakis 2018-06-09 00:33:57 CEST
Correction:  particles above the 'important' energy can be killed, but only after being given repeated changes (the GetThresholdTrials parameter governs this.)
Comment 4 John Apostolakis 2018-11-02 16:53:36 CET
Preparing a revision, to ensure that only stable particles are killed by default - to address stated requirements.

In particular this should ensure that muons, but also pions and kaons, are propagated until they range out or leave the detector by default.
Comment 5 John Apostolakis 2019-03-06 12:45:17 CET
The change to stops stable particles (e-, proton) being killed by the Tranportation processes was part of Geant4 release 10.5