Physics list for detector in nuclear reactor

I’m working on simulated a detector within a nuclear reactor that generates electrical current through nuclear reactions (mostly (n,Beta), (n,y)(y,e-) and (y,e-) reactions). My current physics list is based heavily off of the activation example. My primary particles are a spectrum of neutrons and gamma rays (generated using GPS), in both cases the maximum energy I’m looking at is ~2 MeV. The problem I’m having is that I’m currently getting what appears to be way more electrons from gamma rays than I should be. Any recommendations would be greatly appreciated! My physics list is as follows:

// Hadron Elastic scattering
RegisterPhysics( new HadronElasticPhysicsHP(verb) );

// Hadron Inelastic Physics
RegisterPhysics( new G4HadronPhysicsQGSP_BIC_HP(verb));

// stopping Particles
RegisterPhysics( new G4StoppingPhysics(verb));

// Gamma-Nuclear Physics
RegisterPhysics( new GammaNuclearPhysics(“gamma”));

// EM physics
RegisterPhysics( new G4EmStandardPhysics_option4());

// Decay
RegisterPhysics( new G4DecayPhysics());

// Radioactive decay with gamma rays turned on
RegisterPhysics( new G4RadioactiveDecayPhysics());
G4RadioactiveDecay* radioactiveDecay = new G4RadioactiveDecay();
radioactiveDecay->SetHLThreshold(-1.*s);
radioactiveDecay->SetICM( true );
radioactiveDecay->SetARM( true );

Where hadron elastic is:
G4HadronElasticPhysics::ConstructProcess();
GetNeutronModel()->SetMinEnergy(0.00001eV);
G4HadronicProcess
process = GetNeutronProcess();
G4ParticleHPElastic* model1 = new G4ParticleHPElastic();
process->RegisterMe(model1);
process->AddDataSet(new G4ParticleHPElasticData());

and GammaNuclear is:
G4ProcessManager* pManager = G4Gamma::Gamma()->GetProcessManager();
G4PhotoNuclearProcess* process = new G4PhotoNuclearProcess();
G4CascadeInterface* bertini = new G4CascadeInterface();
bertini->SetMaxEnergy(10*GeV);
process->RegisterMe(bertini);
pManager->AddDiscreteProcess(process);

Thank you!

1 Like

You instantiate RadioactiveDecay process twice : first within G4RadioctiveDecayPhysics constructor and a second time directly. This can make trouble.
In fact I do not understand why you instantiate G4RadioactiveDecay explicitly; you do not need.
It is enough to RegisterPhysics( new G4RadioactiveDEcayPhysics());

Thank you very much. The problem with just using RegisterPhysics( new G4RadioactiveDecayPhysics()); is there there are no gamma rays emitted by radioactive decay by default, which does make a difference in my detector. Do you have any advice on properly enabling gamma ray emission?

Thanks again!

It would be very surprising that radioactive decay does not generate gamma.
Which Geant4 version do you use ?

You are correct! 10.6 now has gamma emission by default. However, I’m still seeing far too many electrons produced by reactions with gamma rays in my detector. Do you have any other suggestions for changes to the physics list?

Once again, thanks so much. You help is greatly appreciated!

We have not enough informations to discuss fruitfuly your problem.
Please try to explain what you do and why you believe to see too much e-

An additionnal question : why do you need to use G4EmStandardPhysics_option4 ? why not the default ElectromagneticPhysics of Activation example ?

Thanks so much for working with me on this! I use sensitive detector and hit collections to track electron generation and movement across my volumes. The results of a run for a vanadium detector in a thermal reactor (both neutrons and prompt gamma spectra for primary particles) show that current generated in the detector (defined as electrons entering collector from the emitter which is the neutron sensitive portion of the detector) has the following results:

Net number of electrons entering collector: 258886
Electrons entering collector that originated in the emitter: 1180359
Electrons entering collector from (n,y) reaction in emitter: 63666
Electrons entering collector from (n,B) reaction in emitter: 5195
Electrons entering collector from reactor gamma reaction in emitter: 1111498

Literature has current in a detector of this type in a thermal reactor being 99% delayed ( (n, beta) reaction), with 1% of the signal being prompt (defined as (n,y)(y,e-) reaction series or (y,e-) reaction if the primary is a gamma ray) which you can see the signal is nothing approaching that. Thus the worry.

About the option4, that was a result of me playing around with different physics packages. Do you recommend the standard option?

Once again, thank you so much for your help. It is greatly appreciated.

I have not big experience with nuclear reactor physics …
What I have understood : you have a volume called emitter made of a mixture of Yttrium and Boron ( ratio ?)
You wish to compare beta decay of nuclear products induced by thermal neutrons and/or gamma (energies of gamma ?)
For debugging, I would keep the electromagneticPhysics constructor of Activation example but comment (eg. do not register) gamma processes (photoelectric, Compton, conversion) to prevent e- production from electromagnetic interactions. Keep also big cut (see PhysicsList::SetCuts()) to prevent e- production from ionisation.
Make a run with neutrons alone and another with gamma alone. And check if your results are coherent …
Additionally you can count the number of neutrinos at creation.

Thanks again for your help! My detector is essentially a semiconductor that produces current through nuclear reactions. I have a volume called the emitter made up of vanadium (natural abundance, so 0.25% 50V and 99.75% 51V). I’m wanting to track all electrons that are generated in this volume that make it out through a second volume (insulator) and to a third volume (collector) so unfortunately all electrons are important (whether they be from neutron capture or gamma processes because the prompt gamma flux in the reactor should be producing a sizable current). What I was seeing that was worrisome was that I was getting too much current from gamma rays compared to whats in literature for this detector. I’ll go ahead and make those changes to see what happens. However, do you still recommend removing those processes? I suspect that they should be there given that will be a source of current in the real detector.

For time being, the question is not to do a realistic simulation but to try to understand what your program does.
I suggest to keep only e- coming directly from a beta decay. They should not exceed the number of anti-neutrinos created.

On the other hand, I have run Activation example. Attached macro and output.

vanadium.mac.txt (231 Bytes)
vanadium.out.txt (5.1 KB)

I am not sure to understand how you count electrons.
Here, an exercice which may be of interest for you.
With a lighty modified Hadr06 (suppress cut 10 km, use G4EM option4), shoot thermal neutrons in a block of natural vanadium, plot e- spectrum at creation. macro in attachment.

  • green : beta decay alone (58 000 entries)
  • blue : the previous one plus e- from (Compton, photoelectric, conversion) processes. 303 484 entries
  • red : the previous one plus e- from atomic deexcitation (Auger cascade). 1 041 994 entries
    There is a lot of Auger elections, but at low energy. Can they reach your collector ?
    neutron.mac.txt (613 Bytes)