Problem 2438

Summary: Wrong velocity in the first step after an optical reflection
Product: Geant4 Reporter: Pierre <pgorel>
Component: processes/opticalAssignee: Daren Sawkey <daren.sawkey>
Status: RESOLVED FIXED    
Severity: normal    
Priority: P4    
Version: 10.7   
Hardware: All   
OS: All   
Attachments: Verbose of a specific photon track showing the problem.

Description Pierre 2021-10-15 19:52:33 CEST
Created attachment 735 [details]
Verbose of a specific photon track showing the problem.

Good day,
I am simulating UV and visible photons with 10.07. I noticed that even when some photons are never leaving the inner material (total internal reflection), I would see some steps with Velocities (=stepLength/stepTime) corresponding to the outer material.

After pocking around the code, I think I found out why (and I think it is happening for all reflections):
deltaT is calculated in G4Transportation::AlongStepDoIt, using stepData.GetPreStepPoint()->GetVelocity(). At this point, the PreStepPoint has received information from the previous step PostStepPoint (Copied by G4SteppingManager).

Now, what do we have in PostStepPoint->GetVelocity() at the end of the previous step?

A reflection is handled by G4OpBoundaryProcess in a 2-step fashion.

    Step 1: Particle starts in “inner volume”. ProposeVelocity in G4OpBoundaryProcess gives the velocity of the inner volume (G4Track::CalculateVelocityForOpticalPhoton). The momentum and polarization are changed. Despite the reflection, the NextVolume is “outer volume” as G4OpBoundaryProcess cannot change the output of the transportation.
    Step 2: The Particle is in “outer volume”, pointing towards towards “inner volume”. Since the step length is 0, G4OpBoundaryProcess exits without doing much, except call CalculateVelocityForOpticalPhoton, which this time gives the velocity of the “outer volume”. Because Propose was called isVelocityChanged is set to true, we don’t recalculate (but it would not change the outcome). pPostStepPoint receives this velocity.

In the step coming after, although we are in the inner volume, PreStepPoint() has the velocity of the outer volume and the deltaT is calculated incorrectly.

To make sure the issue was not with my code, I run examples/extended/optical/OpNovice. I modified slightly SteppingVerbose to show the time and velocity used for each step, as well as the volume the photon is in. You can find the output in attachment.

I don’t quite know how to fix this issue. In G4 9.6, the problem was not apparent because G4Track::GetVelocity() was not re-calculating the Velocity. During the 0-length step, the velocity was from the inner volume, Proposed by G4OpBoundaryProcess, which was then propagated to the step after. It seems that setting G4Track.useGivenVelocity to “true” would bring back this behavior, but this needs to be done when the photons are created, which is not always possible.

Another possibility I see would be to forget about “ProposeVelocity” and call track.GetVelocity() instead of stepData.GetPreStepPoint()->GetVelocity() in G4Transportation::AlongStepDoIt. This way, the velocity will always be matching the material.
Comment 1 Daren Sawkey 2021-10-20 21:50:04 CEST
Hi,

Thank you for reporting this issue. I believe the fix is for these zero length steps to propose the velocity corresponding to the post step point (which corresponds to the inner volume for total internal reflection).

I'll submit a patch, but meanwhile if you want to see if this fixes the problem for you and provide any feedback, that will be appreciated.

https://geant4.kek.jp/lxr/source/processes/optical/src/G4OpBoundaryProcess.cc#L189

becomes:

  if(aTrack.GetStepLength() <= kCarTolerance)
  {
    theStatus = StepTooSmall;
    if(verboseLevel > 1) 
      BoundaryProcessVerbose();

    G4MaterialPropertyVector* groupvel =
      Material2->GetMaterialPropertiesTable()->GetProperty(kGROUPVEL);
    if(groupvel != nullptr)
    {    
      aParticleChange.ProposeVelocity(
        groupvel->Value(thePhotonMomentum, idx_groupvel));
    }    
    return G4VDiscreteProcess::PostStepDoIt(aTrack, aStep);
  }