Secondaries reverse time with wrapped RadioactiveDecay

Particles, Track, Event, Run and Biasing
1

I am using a G4Wrapper to bias the secondaries (direction and number). In the overridden AtRestDoIt I call the pRegProcess-> AtRestDoIt and receive the below output for some decays. Many decays work perfectly fine.

G4WT1 >  G4VParticleChange::CheckSecondary :
G4WT1 >  The global time of secondary goes back compared to the parent !!
G4WT1 >  ParentTime(ns)=3.5172e+16 SecondaryTime(ns)= 4.01313e+14 Difference(ns)=3.47707e+16
G4WT1 >  Ekin(MeV)=9.63804e-06Zr90[2319.000] created by model_RDM_IT
G4WT1 >       -----------------------------------------------
G4WT1 >         G4VParticleChange Information
G4WT1 >         TrackID             : 135
G4WT1 >         ParentID            : 51
G4WT1 >         Particle            : e-
G4WT1 >         Kinetic energy (MeV): 0.0021923056
G4WT1 >         Position (mm)       : (1.0334841e-311,371.31385,136.82683)
G4WT1 >         Direction           : (0.347514,-0.48819155,-0.80056419)
G4WT1 >         PhysicsVolume       : PV-LV_innerbox

In these instances I never make it past my call the the pRegProcess->AtRestDoIt so that code is
likely irrelevant to this issue.

Doing a stack trace I end up the DecayIt with a Mean lifetime of: 188.99310000000000 and a globalTime for the track as 401312856727909.75.Which doesn’t seem to match the error.

This error doesn’t appear if I turn off the G4Wrapper around RadioaciveDecay Physics.

Below is some of my implementation:

G4VParticleChange* RadioactiveDecayWrapper::AtRestDoIt(const G4Track& track, const G4Step& step){

       //get the particle change
       G4VParticleChange* radDecayChange = pRegProcess->AtRestDoIt(track, step);
       G4int noTargets = targetMap.size();
       G4int noSplits = splittingMap.size();

       [...]

        return radDecayChange;
}

And in my physisList

voidPhysicsList::WrapRadioactiveDecay()
{
        G4ProcessManager* processManager = G4GenericIon::GenericIon()- 
        >GetProcessManager();
        G4VProcess* ionRDecay = processManager->GetProcess("RadioactiveDecay");
        if (!ionRDecay) 
        {
            return;
        }

        decayWrap->RegisterProcess(ionRDecay);
        processManager->RemoveProcess(ionRDecay);
        processManager->AddProcess(decayWrap,3, -1, -1);
        G4cout << "Wrapping RadioactiveDecay for generic ions"  << G4endl;
}

_Geant4 Version:_v11.1.2
Operating System: win11
_Compiler/Version:_MSVC v14.39.33519
_CMake Version:_3.28.0

2

Update:

Using break points, it looks like G4Wrapper or my RadioacitveDecayWrapper is not thread safe.

I added some G4out lines to output so I can get the info right before G4RadioactiveDecay takes over:

G4WT4 > globalTime:  0.805504
G4WT4 > Dy particle: Pt181
G4WT4 > Track ID:    31
G4WT4 > Parent Id:   1
G4WT4 >
G4WT4 > globalTime:  1.16208e+10
G4WT4 > Dy particle: Ir181[1680.760]
G4WT4 > Track ID:    32
G4WT4 > Parent Id:   31
G4WT4 >
G4WT4 > globalTime:  1.16208e+10
G4WT4 > Dy particle: Ir181
G4WT4 > Track ID:    361
G4WT4 > Parent Id:   32
G4WT4 >
G4WT4 > globalTime:  3.93938e+11
G4WT4 > Dy particle: Os181
G4WT4 > Track ID:    1787
G4WT4 > Parent Id:   361
G4WT4 >
G4WT2 > globalTime:  1.69613
G4WT2 > Dy particle: Au187
G4WT2 > Track ID:    40
G4WT2 > Parent Id:   1
G4WT2 >
G4WT4 > globalTime:  1.98084e+12
G4WT4 > Dy particle: Re181[432.430]
G4WT4 > Track ID:    2044
G4WT4 > Parent Id:   1787
G4WT4 >
G4WT0 > globalTime:  0.158701
G4WT0 > Dy particle: Pt177
G4WT0 > Track ID:    32
G4WT0 > Parent Id:   1
G4WT0 >
G4WT2 >  G4VParticleChange::CheckSecondary :
G4WT2 >  The global time of secondary goes back compared to the parent !!
G4WT2 >  ParentTime(ns)=0.158701 SecondaryTime(ns)= 4.54839e+11 Difference(ns)=-4.54839e+11
G4WT2 >  Ekin(MeV)=0.101763Ir183 created by model_RDM_Alpha
G4WT2 >       -----------------------------------------------
G4WT2 >         G4VParticleChange Information
G4WT2 >         TrackID             : 32
G4WT2 >         ParentID            : 1
G4WT2 >         Particle            : Pt177
G4WT2 >         Kinetic energy (MeV): 0
G4WT2 >         Position (mm)       : (-10.297004,-3.9912292,-34.6215)
G4WT2 >         Direction           : (0.33584882,-0.45288854,-0.82589197)
G4WT2 >         PhysicsVolume       : PV-LV_innerbox
G4WT2 >         Material            : Average_SV_int
G4WT2 >       -----------------------------------------------
G4WT2 >         # of secondaries    :                    0
G4WT2 >       -----------------------------------------------
G4WT2 >         Energy Deposit (MeV):                    0
G4WT2 >    NIEL Energy Deposit (MeV):                    0
G4WT2 >         Track Status        :         StopButAlive
G4WT2 >         TruePathLength (mm) :                    0
G4WT2 >         Stepping Control    :                    0
G4WT2 >        First step in volume
G4WT2 >       -----------------------------------------------
G4WT2 >     G4ParticleChangeForDecay proposes:
G4WT2 >     Proposed local Time (ns):          0.011576413
G4WT2 >     Initial local Time (ns) :          0.011576413
G4WT2 >     Initial global Time (ns):           0.15870094
G4WT2 >     Current global Time (ns):           0.15870094
G4WT2 >
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : TRACK001
      issued by : G4VParticleChange::CheckSecondary()
Secondary with illegal time and/or energy and/or momentum
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

G4WT2 > globalTime:  0.159008
G4WT2 > Dy particle: Ir183
G4WT2 > Track ID:    41
G4WT2 > Parent Id:   40
G4WT4 >
G4WT4 > globalTime:  1.98084e+12
G4WT4 > Dy particle: Re181[356.720]
G4WT4 > Track ID:    2362
G4WT4 > Parent Id:   2044
G4WT4 >
G4WT0 >  G4VParticleChange::CheckSecondary :
G4WT0 >  The global time of secondary goes back compared to the parent !!
G4WT0 >  ParentTime(ns)=1.98084e+12 SecondaryTime(ns)= 8.83298e+09 Difference(ns)=1.97201e+12
G4WT0 >  Ekin(MeV)=4.46861e-05Ir177 created by model_RDM_BetaPlus
G4WT0 >       -----------------------------------------------
G4WT0 >         G4VParticleChange Information
G4WT0 >         TrackID             : 2362
G4WT0 >         ParentID            : 2044
G4WT0 >         Particle            : Re181[356.720]
G4WT0 >         Kinetic energy (MeV): 0
G4WT0 >         Position (mm)       : (142.37867,116.7671,-28.453447)
G4WT0 >         Direction           : (0.52273975,-0.80202952,-0.28894948)
G4WT0 >         PhysicsVolume       : PV-LV_innerbox
G4WT0 >         Material            : Average_SV_int
G4WT0 >       -----------------------------------------------
G4WT0 >         # of secondaries    :                    0
G4WT0 >       -----------------------------------------------
G4WT0 >         Energy Deposit (MeV):                    0
G4WT0 >    NIEL Energy Deposit (MeV):                    0
G4WT0 >         Track Status        :         StopButAlive
G4WT0 >         TruePathLength (mm) :                    0
G4WT0 >         Stepping Control    :                    0
G4WT0 >        First step in volume
G4WT0 >       -----------------------------------------------
G4WT0 >     G4ParticleChangeForDecay proposes:
G4WT0 >     Proposed local Time (ns):          0.013354366
G4WT0 >     Initial local Time (ns) :          0.013354366
G4WT0 >     Initial global Time (ns):        1.9808403e+12
G4WT0 >     Current global Time (ns):        1.9808403e+12
G4WT0 >
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : TRACK001
      issued by : G4VParticleChange::CheckSecondary()
Secondary with illegal time and/or energy and/or momentum
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

Looking at the output from the wrapper the globalTime is 0.158701 and the dynamic particle is Pt-177. When G4RadioactiveDecay gets ahold of the track it is all messed up. It’s getting the correct track (ID: 32) but it’s on the wrong thread. The timing makes sense but not why it’s throwing an error. With break points we can look at the data in the track 32:
Initial Global time: 0.158701 (parentTime)
final Global time: 8.8329e+09 (secondaryTime)
Parent Track ID: 1

If we continue on, we can see Re181[356.720] is also on the wrong thread, and that it has the secondaryTime from the traced Pt-177 (8.8329e+09). I’ve also included the entire decay chain for the Re181[356.720] starting at the Pt-181 output from the wrapper.

The final global time/secondary time is being attributed to the incorrect track on the wrong thread. In some runs I even receive a segfault, as it tries to add a secondary to a particleChange object on another thread.

It was my understating per this page that each thread handles an entire event.

3

Final update and solved. The issue arose in the physicsList where I construct the RadioactiveDecayWrapper:

I created the member variable RadioactiveDecayWrapper* decayWrap in my header, then built in my constructor. This made one decayWrap for all the threads and it wasn’t thread safe.

The fix was to create a local variable in the WrapRadioactiveDecay() function (see below).

void PhysicsList::WrapRadioactiveDecay()
{
        G4ProcessManager* processManager = G4GenericIon::GenericIon()- 
        >GetProcessManager();
        G4VProcess* ionRDecay = processManager->GetProcess("RadioactiveDecay");
        if (!ionRDecay) 
        {
            return;
        }

        RadioactiveDecayWrapper* decayWrap = new RadioactiveDecayWrapper(<parameters>);
        decayWrap->RegisterProcess(ionRDecay);
        processManager->RemoveProcess(ionRDecay);
        processManager->AddProcess(decayWrap,3, -1, -1);
        G4cout << "Wrapping RadioactiveDecay for generic ions"  << G4endl;
}

There were some code refactoring required such as passing in the needed parameters into the constructor of the RadioactiveDecayWrapper but it works now.

4

This topic was automatically closed 7 days after the last reply. New replies are no longer allowed.