Problem 2166 - wrong excitation for recoil nucleus in G4ParticleHPInelasticCompFS
Summary: wrong excitation for recoil nucleus in G4ParticleHPInelasticCompFS
Status: RESOLVED FIXED
Alias: None
Product: Geant4
Classification: Unclassified
Component: processes/hadronic/models/neutron_hp (show other problems)
Version: 10.5
Hardware: All All
: P4 trivial
Assignee: dennis.herbert.wright
URL:
Depends on:
Blocks:
 
Reported: 2019-06-08 20:19 CEST by Artem Zontikov
Modified: 2019-06-25 18:05 CEST (History)
1 user (show)

See Also:


Attachments

Note You need to log in before you can comment on or make changes to this problem.
Description Artem Zontikov 2019-06-08 20:19:25 CEST
Hello. 

There is an issue in particle_hp model which results in unlimited growth of G4IonTable. 

If we are dealing with an isotope which 
- does not have determined excitation levels for some interaction; 
- does not have angular of energy-angular distributions for outgoing particle; 
- has photon data represented as multiplicities and contimuum or arbitary energy distributions; 
then line #608 in G4ParticleHPInelasticCompFS.cc 
G4DynamicParticle* targ = new G4DynamicParticle( G4IonTable::GetIonTable()->GetIon ( (G4int)theBaseZ , (G4int)theBaseA , totalPhotonEnergy )  , G4ThreeVector(0) ); 
inserts a new entry into G4IonTable via consecutive calls of GetIon(), FindIon() , FindIonInMaster(), CreatIon() and InsertWorker(). 

When G4IonTable tries to find an ion it checks, among other things, the excitation energy. Sampling the energy and multiplicity of final state photons from G4NDL creates a unique "totalPhotonEnergy" which later is used as an excitation. Thus G4IonTable is expanding and worker threads are waiting longer to acquire the mutex for it. 

This behaviour could be reproduced shooting 14MeV neutrons into Mg target. The (n,p) and (n,a) reactions generate photons which are sampled from distributions. Even on a 4-core machine the speed halves in several minutes from starting the execution. 

Proposed change for line 608: 
G4DynamicParticle* targ = new G4DynamicParticle( G4IonTable::GetIonTable()->GetIon ( (G4int)theBaseZ , (G4int)theBaseA , theFinalStatePhotons[it] ? theFinalStatePhotons[it]->GetLevelEnergy() : 0 )  , G4ThreeVector(0) );
Comment 1 dennis.herbert.wright 2019-06-24 22:23:53 CEST
Hi Artem,

   Most of this seems OK, except for the case when theFinalStatePhotons[it] = 0.
Then you assign zero excitation energy.  In that case it would seem better to use totalPhotonEnergy, would it not?

Dennis
Comment 2 Artem Zontikov 2019-06-25 12:52:22 CEST
I guess that theFinalStatePhotons[it] == 0 results in thePhotons == 0 which leaves no other option but to have a zero totalPhotonEnergy. 

Also, proposed fix makes these lines redundant: 
587       G4double totalPhotonEnergy = 0.0;
588       if ( thePhotons != 0 )
589       {
590          unsigned int nPhotons = thePhotons->size();
591          unsigned int ii0;
592          for ( ii0=0; ii0<nPhotons; ii0++)
593          {
594             //thePhotons has energies at LAB system 
595             totalPhotonEnergy += thePhotons->operator[](ii0)->GetTotalEnergy();
596          }
597       }
Comment 3 dennis.herbert.wright 2019-06-25 18:05:31 CEST
Your fix has been implemented as proposed.