Message: Re: Electrons misteriously killed Not Logged In (login)
 Next-in-Thread Next-in-Thread
 Next-in-Forum Next-in-Forum

None Re: Electrons misteriously killed  

Forum: Fields: Magnetic and Otherwise
Re: Question Electrons misteriously killed (Lafont Fabien)
Re: Question Re: Electrons misteriously killed (Lafont Fabien)
Date: 14 Apr, 2014
From: John Apostolakis <John Apostolakis>

Hi Fabien,

I think that you have found a feature of propagation in field - which can kill particles either verbosely or silently - but one that you can control using existing methods built into the G4Transportation process. *** (see footnote)

When the propagation in field undertakes a lot of integration steps and geometry inquiries without encountering a boundary - or reaching the location of a physics process - it will give up, and signal its failure to complete the requested tracking step within the (large) budget of trial steps to the Transportation process.

In response the Transportation process 
   - can either issue a warning, 
   - if the energy of the particle is lower than a threshold (default is 250 MeV I believe) it will not issue it and will kill the track silently.

There are a few member variables in G4Transportation which govern this behaviour: 

 // Thresholds for looping particles: 
  // 
     G4double fThreshold_Warning_Energy;     //  Warn above this energy
     G4double fThreshold_Important_Energy;   //  Hesitate above this
     G4int    fThresholdTrials;              //    for this no of trials
       // Above 'important' energy a 'looping' particle in field will 
       //   *NOT* be abandoned, except after fThresholdTrials attempts.

These parameters can be changed - either by revising their default values in the Geant4 source code (G4Transportation.cc) or better by using the public methods of G4Transportation created for this purpose: 

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

The default values are set in the constructor of G4Transportation:
    fThreshold_Warning_Energy( 100 * MeV ),  
    fThreshold_Important_Energy( 250 * MeV ), 
    fThresholdTrials( 10 ), 

The reasons for these default are to avoid a large number of warnings in high energy simulations - but clearly you will need to activate the warnings, and increase the number of trials in order to make sure that tracks are not abandoned except when you cannot any longer afford to track them.

I expect that you may also be interested to use some different steppers, in case your computation time is problematic.  Please look at the Geant4 documentation and consider using G4HelixExplicitEuler  or potentially G4HelixSimpleRunge if this is the case.

Best regards,
John

Note ***: In case you are using parallel navigation, due to the presence of additional geometries, it may be a bit different.  Please signal if this is the case.
===================================================
John Apostolakis,  PH Department, CERN
SFT (SoFTware for Experiments) Group


On Apr 14, 2014, at 4:52 PM, Lafont Fabien wrote:

> 
> *** Discussion title: Fields: Magnetic and Otherwise
> 
> Sorry for the log file format ...
> 
>>> AlongStepDoIt (process by process):    Process Name = Transportation
> ++G4Step Information
> Address of G4Track    : 0x34f80b40
> Step Length (mm)      : 13.75406869033371
> Energy Deposit (MeV)  : 0
> -----------------------------------------------------------------------
> StepPoint Information               PreStep            PostStep
> -----------------------------------------------------------------------
> ...
> Total Energy (MeV)  :   0.9772813196679687  0.9772813196679687
> Kinetic Energy (MeV):   0.4662824096679687  0.4662824096679687
> ...
> Process defined Step:            Undefined      Transportation
> -----------------------------------------------------------------------
> !Note! Safety of PostStep is only valid after all DoIt invocations.
> ++G4ParticleChange Information
> 
> -------------------------------------------------------------
> Visit this GEANT4 at hypernews.slac.stanford.edu message (to reply or unsubscribe) at: 
> http://hypernews.slac.stanford.edu/HyperNews/geant4/get/emfields/213/1.html 

Inline Depth:
 1 1
 All All
Outline Depth:
 1 1
 2 2
 All All
Add message: (add)

1 None: Re: Electrons misteriously killed   (Lafont Fabien - 16 Apr, 2014)
(_ None: Re: Electrons misteriously killed   (John Apostolakis - 16 Apr, 2014)
(_ Question: Re: Electrons misteriously killed   (Gumplinger Peter - 16 Apr, 2014)
(_ None: Re: Electrons misteriously killed   (John Apostolakis - 17 Apr, 2014)
 Add Message Add Message
to: "Re: Electrons misteriously killed "

 Subscribe Subscribe

This site runs SLAC HyperNews version 1.11-slac-98, derived from the original HyperNews


[ Geant 4 Home | Geant 4 HyperNews | Search | Request New Forum | Feedback ]