Hadronic physics fails in my custom physics list after migrating from Geant4 10.6.1 to 10.6.3

This seems like déjà vu. I recently upgraded from Geant4 10.6.1 to 10.6.3. To check that everything was installed correctly, I recompiled and ran a program that was working fine under 10.6.1. The program models scintillation in different scintillator/PMT combinations irradiated by various hadronic particle and gamma ray beams of low energy (<20 MeV). It outputs a file of the number of scintillation and PMT photons for each event and it has been working for several years. It broke when upgrading from Geant4 10.5.1 to 10.6.1 but was fixed to run under the latter. (See forum item “Hadronic physics fails to work in my user physics list after migrating from Geant4 10.5.1 to 10.6.1” in Physics Lists for a history of that.)

The program runs OK when the incident particles are gamma rays, electrons, alphas or ions other than a hydrogen ion. When I use neutrons or protons, the program fails by throwing one or more exceptions and a core dump. For example, for 1000 MeV primary neutrons incident on an NE213 scintillator, the error message is:

$ ./OpNovice -m Neutron.mac -t 8 > Neutron_1000keV_pos0_0_774_iso_1000_OpNovice1.out
G4WT5 >
-------- EEEE ------- G4Exception-START -------- EEEE -------
*** G4Exception : had001
issued by : G4CrossSectionDataStore::GetIsoCrossSection
No isotope cross section found for proton off Element Hydrogen in NE213 Z= 1 A= 1 E(MeV)= 0.367108

*** Fatal Exception *** core dump ***
G4WT5 > G4Track (0x7fc8083a8d30) - track ID = 3, parent ID = 1
G4WT5 > Particle type : proton - creator process : hadElastic, creator model : Undefined
G4WT5 > Kinetic energy : 367.108 keV - Momentum direction : (-0.685377,0.363223,-0.631132)
G4WT5 > Step length : 0 fm - total energy deposit : 0 eV
G4WT5 > Pre-step point : (-0.721052,4.28328,7.65669) - Physical volume : PhysicalScintillator (NE213)
G4WT5 > - defined by : not available
G4WT5 > Post-step point : (-0.721052,4.28328,7.65669) - Physical volume : PhysicalScintillator (NE213)
G4WT5 > - defined by : not available
G4WT5 > *** Note: Step information might not be properly updated.
G4WT5 >
-------- EEEE -------- G4Exception-END --------- EEEE -------

G4WT5 >
G4WT5 > *** G4Exception: Aborting execution ***
Aborted (core dumped)

The output also says:
NeutronHP: /Inelastic file for Z = 1, A = 1 is not found and CrossSection set to 0. (This message occurs regardless of the primary particle.)

So it seems to be unable to find cross sections for protons (presumably the recoils from the elastic neutron scattering) on hydrogen. However, I checked that the appropriate cross section libraries are installed, and they are:
$env | grep G4 | sort | sed ‘s/=/ /’ | column -t
G4ABLADATA /home/jmcfee/geant4/geant4.10.6.p03-install/share/Geant4-10.6.3/data/G4ABLA3.1
G4ENSDFSTATEDATA /home/jmcfee/geant4/geant4.10.6.p03-install/share/Geant4-10.6.3/data/G4ENSDFSTATE2.2
G4INCLDATA /home/jmcfee/geant4/geant4.10.6.p03-install/share/Geant4-10.6.3/data/G4INCL1.0
G4LEDATA /home/jmcfee/geant4/geant4.10.6.p03-install/share/Geant4-10.6.3/data/G4EMLOW7.9.1
G4LENDDATA /home/jmcfee/geant4/geant4.10.6.p03-install/share/Geant4-10.6.3/data/LEND_GND1.3_ENDF.BVII.1
G4LEVELGAMMADATA /home/jmcfee/geant4/geant4.10.6.p03-install/share/Geant4-10.6.3/data/PhotonEvaporation5.5
G4NEUTRONHPDATA /home/jmcfee/geant4/geant4.10.6.p03-install/share/Geant4-10.6.3/data/G4NDL4.6
G4NEUTRONHP_DO_NOT_ADJUST_FINAL_STATE 1
G4NEUTRONHP_NEGLECT_DOPPLER 0
G4NEUTRONHP_SKIP_MISSING_ISOTOPES 1
G4NEUTRONHP_USE_ONLY_PHOTONEVAPORATION 0
G4PARTICLEHPDATA /home/jmcfee/geant4/geant4.10.6.p03-install/share/Geant4-10.6.3/data/G4TENDL1.3.2
G4PARTICLEXSDATA /home/jmcfee/geant4/geant4.10.6.p03-install/share/Geant4-10.6.3/data/G4PARTICLEXS2.1
G4PIIDATA /home/jmcfee/geant4/geant4.10.6.p03-install/share/Geant4-10.6.3/data/G4PII1.3
G4PROTONHPDATA /home/jmcfee/geant4/geant4.10.6.p03-install/share/Geant4-10.6.3/data/G4TENDL1.3.2/Proton
G4RADIOACTIVEDATA /home/jmcfee/geant4/geant4.10.6.p03-install/share/Geant4-10.6.3/data/RadioactiveDecay5.4
G4REALSURFACEDATA /home/jmcfee/geant4/geant4.10.6.p03-install/share/Geant4-10.6.3/data/RealSurface2.1.1
G4SAIDXSDATA /home/jmcfee/geant4/geant4.10.6.p03-install/share/Geant4-10.6.3/data/G4SAIDDATA2.0

(The data files are the latest available and are indeed at the locations indicated.)

I must stress that this example worked perfectly a few weeks ago under 10.6.1. I have made no changes to my code or macro files. The 10.6.2 and 10.6.3 Release Notes do not show any changes since 10.6.1 that would explain the problem.

Am I missing a library? As always any help is gratefully appreciated. Thanks.

My physics list attached.OpNovicePhysicsList.cc (25.1 KB)

I upgraded to Geant4 10.7 and the problem remains unchanged. I can circumvent the problem by writing a physics list that just addresses the optical physics and using the modular physics list method to load the EM and hadronic physics production lists together with my optical physics list . (I have completed this and am currently testing it. The problem in this post has disappeared but I am not completely sure that my optical physics list is being used.) However, the present problem is a real one, since it’s disconcerting when code breaks everytime there is a minor distribution upgrade.

Hello, John,

I am sorry, but this Physics List is problematic. Our recommendation is to use reference Physics Lists (see $G4INSTALL/examples/extended/hadronic/Hadr00 or Hadr01). All reference Physics Lists are working fine in these releases. Any Physics Lists organized in a way like you did will required very detail study of Geant4 modifications between different versions.

If you prefer do not use reference Physics Lists, then our recommendation to use G4VModularPhysicsList interface and combine components from Geant4 distribution. In particular, in your method OpNovicePhysicsList::ConstructHad() obsolete GHEISHA cross section classes are used, also I am not sure about model combination. OpNovicePhysicsList::ConstructEM() is also obsolete. I would instead use all components from Geant4 distribution and make your own optical physics if you disagree with Geant4 optical builder.

VI

Hello Vladimir,

Thank you for replying. You are correct about the problems with building my own complete physics list from scratch and I have abandoned that approach. It was too much of a headache to keep debugging it and rewriting it.

But, I needed my own optical physics process to handle the particle-dependent fast/slow scintillation component yield ratios, so I built my own optical physics process module. I just finished it and succeeded in adding it to an existing reference physics list using a Modular Physics list approach. It seemed to work, with a few glitches.

Before I could iron those out, Geant4 10.7.0 was released with its enhanced time constant methods and properties. By rewriting my MaterialsPropertyTable with the new methods available, I no longer need my own custom Physics List and I can use the normal reference ones, such as are in OpNovice.cc.

A big thank you to the optical physics development team.