Message: Re: Strange behaviour for AtRest processes Not Logged In (login)
 Next-in-Thread Next-in-Thread
 Next-in-Forum Next-in-Forum

None Re: Strange behaviour for AtRest processes 

Forum: Event and Track Management
Re: Warning Strange behaviour for AtRest processes (Daniel Barna)
Date: 28 Nov, 2012
From: Marc Verderi <Marc Verderi>

Dear Daniel,

I came back on the issue you exposed, and I believe you're right. Please 
submit a bug report on this issue to keep track of the problem and of 
subsequent fix.

Best regards,
Marc

On 11/02/2012 02:07 PM, Daniel Barna wrote:
> *** Discussion title: Event and Track Management
>
> Hello, According to the documentation, when a particle is at rest, the
> rest-process proposing the shortest interaction length will be called,
> the other rest processes associated with this particle not (unless they
> are forced).
>
> However, this does not work. If I associate these two processes below
> with a particle, it is always only MyProcessA that should be called,
> since this is the one with the shorter interaction length. Both are
> called, as can easily be tested. If I invert the GPIL values (MyProcessA
> returning 20, and MyProcessB returning 10), then only MyProcessB is
> called. This depends on the ordering of the processes, as described
> below this code...
>
> [code]
>
> class MyProcessA : public G4VRestProcess
> {
> public:
>      G4bool IsApplicable(const G4ParticleDefinition& ) { return true; }
>      MyProcessA() : G4VRestProcess("MyProcessA",fHadronic) {}
>      G4double AtRestGetPhysicalInteractionLength(const G4Track& , G4ForceCondition* condition)
>          {
>              *condition = NotForced;
>              cerr<<"MyProcessA::AtRestGetPhysicalInteractionLength returning 10"<<endl;
>              return 10;
>          }
>      G4VParticleChange *AtRestDoIt(const G4Track &track, const G4Step &)
>          {
>              aParticleChange.Initialize(track);
>              aParticleChange.SetNumberOfSecondaries(0);
>              aParticleChange.ProposeTrackStatus(fStopAndKill);
>              cerr<<"MyProcessA::AtRestDoIt"<<endl;
>              return &aParticleChange;
>          }
>      G4double GetMeanLifeTime(const G4Track&, G4ForceCondition*) { return 0; }
> };
>
> class MyProcessB : public G4VRestProcess
> {
> public:
>      G4bool IsApplicable(const G4ParticleDefinition& ) { return true; }
>      MyProcessB() : G4VRestProcess("MyProcessB",fHadronic) {}
>      G4double AtRestGetPhysicalInteractionLength(const G4Track& , G4ForceCondition* condition)
>          {
>              *condition = NotForced;
>              cerr<<"MyProcessB::AtRestGetPhysicalInteractionLength returning 20"<<endl;
>              return 20;
>          }
>      G4VParticleChange *AtRestDoIt(const G4Track &track , const G4Step &)
>          {
>              aParticleChange.Initialize(track);
>              aParticleChange.SetNumberOfSecondaries(0);
>              aParticleChange.ProposeTrackStatus(fStopAndKill);
>              cerr<<"MyProcessB::AtRestDoIt"<<endl;
>              return &aParticleChange;
>          }
>      G4double GetMeanLifeTime(const G4Track&, G4ForceCondition*) { return 0; }
> };
> [/code]
>
> The problem lies in G4SteppingManager::InvokeAtRestDoItProcs(), see my
> comments in this code:
>
> [code]
> void G4SteppingManager::InvokeAtRestDoItProcs()
> {
>     // Select the rest process which has the shortest time before
>     // it is invoked. In rest processes, GPIL()
>     // returns the time before a process occurs.
>     G4double lifeTime, shortestLifeTime;
>
>     fAtRestDoItProcTriggered = 0;
>     shortestLifeTime = DBL_MAX;
>
>     unsigned int NofInactiveProc=0;
>     for( size_t ri=0 ; ri < MAXofAtRestLoops ; ri++ ){
>       fCurrentProcess = (*fAtRestGetPhysIntVector)[ri];
>       if (fCurrentProcess== 0) {
>         (*fSelectedAtRestDoItVector)[ri] = InActivated;
>         NofInactiveProc++;
>         continue;
>       }   // NULL means the process is inactivated by a user on fly.
>
>       lifeTime =
>         fCurrentProcess->AtRestGPIL( *fTrack, &fCondition );
>
>       if(fCondition==Forced && fCurrentProcess){
>         (*fSelectedAtRestDoItVector)[ri] = Forced;
>       }
>       else{
>         (*fSelectedAtRestDoItVector)[ri] = InActivated;
>         if(lifeTime < shortestLifeTime ){
>            shortestLifeTime = lifeTime;
>            fAtRestDoItProcTriggered =  G4int(ri);
>
>            // ----------  MY COMMENTS ---------------------------
>            // this line below should be outside of the loop!!!
>            // If the processes happen to return decreasing lifetimes,
>            // all processes will be set to NotForced in sequence, and
>            // therefore be called in the later part of this function!
>            (*fSelectedAtRestDoItVector)[fAtRestDoItProcTriggered] = NotForced;
>         }
>       }
>     }
>     ......;
> }
> [/code]
>
> Did I miss something or is it indeed a (quite serious) bug? I am using
> geant4.9.5.p01. I see that patch-02 is out now, but this code is the
> same there.
>
> -------------------------------------------------------------
> Visit this GEANT4 at hypernews.slac.stanford.edu message (to reply or unsubscribe) at:
> http://hypernews.slac.stanford.edu/HyperNews/geant4/get/eventtrackmanage/1091.html

 Add Message Add Message
to: "Re: Strange behaviour for AtRest processes"

 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 ]