Helium-3 neutron capture cross section doesn't seem right

Physics Processes, Models and Cross Sections
1

Dear Geant4 experts,
I am tryng to simulate a He3-based neutron counter. However, in a test application, Helium-3 neutron cross section appears to be two orders of magnitude lower than expected: ~11 barn versus ~5300 barn.

I have constructed He3 G4Material like so:

//He3 element
G4double atomicMass=3.016*g/mole;
G4Isotope* he3 = new G4Isotope("he3", 2, 3, atomicMass);
G4Element* He3 = new G4Element("He3", "He3", 1);
He3->AddIsotope(he3,100.*perCent);

//He3 gas
G4double    pressure = 5*bar,
            temperature = 293*kelvin,
            molar_constant = CLHEP::Avogadro*CLHEP::k_Boltzmann,
            density = (atomicMass*pressure)/(temperature*molar_constant);
G4Material* He3Gas = new G4Material("He3_gas", density,	1, kStateGas, temperature, pressure);
He3Gas->AddElement(He3, 100.*perCent);

Then, I used this material to build a simple detector, a tube along z axis:
изображение
The primary particles are set to be neutrons, and are launched with thermal (0.025 eV) energy along z axis.

I have used a modular physics list, where I registered the following physics (which, to my understanding, contain the most accurate low energy neutron processes):

myPhysicsList::myPhysicsList():G4VModularPhysicsList(){
int ver = 0;
this->RegisterPhysics( new G4EmStandardPhysics(ver) );
G4String evaluation("ENDF/BVII.1");
this->RegisterPhysics( new G4HadronElasticPhysicsLEND(ver,evaluation) );
this->RegisterPhysics( new G4HadronPhysicsShieldingLEND(ver) );
this->RegisterPhysics( new G4ThermalNeutrons(ver) );
}

The following log entry appears in physics summary at runtime:

                       Hadronic Processes for neutron

  Process: hadElastic
        Model:             hElasticCHIPS: 19.5 MeV ---> 100 TeV
        Model:               LENDElastic: 4 eV  ---> 20 MeV
        Model: NeutronHPThermalScattering: 0 eV  ---> 4 eV 
     Cr_sctns: NeutronHPThermalScatteringData: 0 eV  ---> 4 eV 
     Cr_sctns:   LENDElasitcCrossSection: 0 eV  ---> 20 MeV
     Cr_sctns:        G4NeutronElasticXS: 0 eV  ---> 100 TeV

  Process: neutronInelastic
        Model:                      FTFP: 3 GeV ---> 100 TeV
        Model:            BertiniCascade: 19.9 MeV ---> 6 GeV
        Model:             LENDInela    stic: 0 eV  ---> 20 MeV
         Cr_sctns: LENDInelasitcCrossSection: 0 eV  ---> 20 MeV
         Cr_sctns: JENDLHEInelasticCrossSection: 0 eV  ---> 100 TeV
         Cr_sctns:  BarashenkovGlauberGribov: 0 eV  ---> 100 TeV

  Process: nCapture
        Model:               LENDCapture: 0 eV  ---> 20 MeV
        Model:               nRadCapture: 19.9 MeV ---> 100 TeV
     Cr_sctns:   LENDCaptureCrossSection: 0 eV  ---> 20 MeV
     Cr_sctns:        G4NeutronCaptureXS: 0 eV  ---> 100 TeV

  Process: nFission
        Model:               LENDFission: 0 eV  ---> 20 MeV
        Model:                G4LFission: 19.9 MeV ---> 100 TeV
     Cr_sctns:   LENDFissionCrossSection: 0 eV  ---> 20 MeV
     Cr_sctns:          GheishaFissionXS: 0 eV  ---> 100 TeV

With this, I get ~1200 hits per 2E6 neutrons, which, at 10cm He3 volume at 5bar pressure, equals 11.5 barns.
Is there a cross section not defined? If so, how can I define it?

Any help will be highly appreciated.

2

Attached a macro for example Hadr03 and the result.
akim.mac.txt (184 Bytes)
akim.out.txt (2.9 KB)

1 Like
3

Thank you for your reply, @maire!

It was very helpful of you to supply a concrete instruction for me to validate physics on my end. One should be careful, however, as the macro contains incorrect density for helium: /testhadr/det/setIsotopeMat He3 2 3 177.43 mg/cm3 will set the density to about 3 orders more than it needs to be (0.1346*mg/cm3). This can’t influence the per-atom cross section that Hadr03 calculates, and I guess that is why it came unnoticed.

I was able to solve the problem in the end. There is, however, a peculiar aspect to the detection. I am writing this for anybody who may encounter a similar feature.

In my setup described above, the He3 cylinder is set to be a SensitiveDetector. I measured the cross section “manually”: if N primary neutrons cause the detector to register M hits, then, given that the capture cross section is dominating at thermal energies, we can say that, approximately, M/N=1-exp(-snL), where s is capture cross section, n is atomic density (2.68*10^19 1/cm3 here), L is in-material path (10 cm here).
My SensitiveDetector was, effectively, writing down all non-zero occurrences of
energy = step->GetTotalEnergyDeposit()
in ProcessHits method. This way, I was getting about 9 hits per 10000 shots (which equals about 3.36 barn). Then, I switched to writing all more-than-zero
energy = step->GetPreStepPoint()->GetKineticEnergy() - step->GetPostStepPoint()->GetKineticEnergy()
Now, each 10000 primary neutrons yields about 7600 hits – and that gives about 5400 barn, as it should be!

Can anyone shed light on why GetTotalEnergyDeposit() does not include capture?
Other than that, problem solved.

4

That may be due to production of secondary particles. In my very limited experience with Geant4, step->GetPreStepPoint()->GetKineticEnergy() - step->GetPostStepPoint()->GetKineticEnergy() = step->GetTotalEnergyDeposit() + (all secondary->GetKineticEnergy()).

This is because Geant4 uses production thresholds via range.
If a neutron generates a secondary which cannot travel the range cut (usually 1mm), then the secondary energy is output as GetTotalEnergyDeposit(). Otherwise, the secondary will track through the volume (now being considered as a primary) and deposit energy in the same way.

Secondary production is discrete, while GetTotalEnergyDeposit() is continuous.

If you remove the primary check, and register all total energy deposit, it should work:
instead of

if(primary)
{
energy = step->GetTotalEnergyDeposit()
}

simply

energy = step->GetTotalEnergyDeposit()

ought to work.
Since all produced secondaries can then track out the energy through continuous processes - which are accounted in GetTotalEnergyDeposit.

If you are already using the second option, it may be that your produced secondaries are leaving the He3 cylinder (carrying energy with them) and decreasing the amount of energy logged.

Hope this helps.
Akhil

5

Thank you for the detailed explanation, @akhilsadam. My understanding of this problem is now complete. All issues are resolved.