Protons and neutrons in 27Al(p, x) reactions

Dear FLUKA users and expects,

I’m trying to score neutron fluence energy distribution (spectrum) born on 27Al(p, x) reaction on aluminum attenuator of a cyclotron target assembly.
I would like to know the energy distribution of the neutrons born in inelastic hadron reactions in the attenuator and than bombarding the molybdenum target sample (which contains 98Mo isotope and can undergo (n, \gamma) reaction).

There are two problems:

1. Energy of the protons

The proton beam initial energy is 18 MeV distributed according to the Gaussian law with \sigma_E \approx 0.12 MeV. It also has an annual shape (3 cm in diameter), but than the beam is collimated with the target assembly (ABSORB reg.) which has inner diameter of 0.6 cm. The beam starts to run towards the attenuator (ALFOIL3 reg.) through air and loses some kinetic energy. The first data I have to know what’s the precise proton energies just before hitting the attenuator. I use USRBDX card to score the quantity, but the result happens to be not what I’d expected. The shape of the curve looks like the Gaussian, but it is not symmetrical and has slightly distorted shape. Please, have a look at the image below:

While scoring after leaving the attenuator (boundary crossing between ALFOIL3 and AIR2) gives me pretty nice shape and symmetry:

What might cause this strange shape of the E_p curve?

2. Neutron spectrum

The most important goal of the simulation is getting the energy spectrum of the neutrons born in the inelastic nuclear reaction on 27Al and moving towards the molybdenum sample (POCKET reg.). I also use the USRBDX card to find it. I’m interested in energies of the neutrons crossing the boundary between the attenuator and air (AIR2 reg.) and also those which are moving thought air (AIR2) and Mo-sample (POCKET). Let’s focus on the first sort of them. You can see the result on the image below:

As you can see the distribution has large uncertainties and does not follow an expected shape (the outlier about 0.7 MeV). I ran 3E7 primaries for 5 cycles using 10 CPU threads (Intel i7 9th gen) and it took about 8 hours. Is the number of protons enough? I suspect no(?), than how can I increase FoM and speed the calculation up? (not including increasing the number of primaries)
Probably I have to use a more powerful machine to run it (?). May I use a variance reduction technique (and which one?)? Won’t it affect the investigation of the physical processes?

The spectrum of the neutrons hitting the Mo-sample is much more worse:

What else can I do to improve the results and get a nice proper neutrons and protons? Physics, scoring, etc …

Sincerely,
Ihor

Input files

sim_v2.flair (15.2 KB)
sim_v2.inp (12.9 KB)

Versions

FLUKA: 4-5.1
Flair: 3.4-4
OS: Ubuntu 24.04 LTS

Hello Ihor,

thank you for your patience.

  1. You actually set up a rectangular energy spread for the beam, as you can check by defining a second AIR region right in front of the source, and see what is emitted. Then, travelling through the air, the beam loses energy. What you see by the time you score it is a convolution of the rectangular window with a Landau (approximately) distribution, which is asymmetric. So no surprises on the first point.

  2. This is a perfect application of LAM-BIAS. You can change the inelastic scattering length of proton in a material, biasing the number of inelastic interactions. The current inelastic scattering length, as you can find the output file, of proton on Al is 21 cm, whereas you only have ~1mm. I recommend scaling the lambda by 0.005, so that on average each proton will undergo an inelastic interaction. Lam-biasing takes care of the weights of the particles, so you won’t negatively affect the results.

Hope this helps,

Stefano

Hello Stefano,

Thank you for the advices!

I really made a “typo” in the \Delta p function! I selected flat instead of Gauss distribution. Now it’s fixed and the protons spectrum looks well!

LAM-BIAS also helped to reduce uncertainties, you can compare these before (right) and after biassing (left) plots below.

This is another image of the new distribution (I ran 1e7 primary histories this time).

As you can see there are still few points which have the large errors and do not follow the expected shape (below 1 MeV). Maybe there is a problem with slow neutron treatment? What do you think?

sim_v2.flair (24.7 KB)
sim_v2.inp (13.1 KB)

Thanks,
Ihor

Hi Ihor,

could it just be not sufficient statistics. To speed up the simulation further, I recommend setting high EMFCUTS and delta ray production, as well as cutting the proton when it goes belo 2-3 MeV, since the production of neutrons at that point has vanished. With 5e8 primaries I do not see a big discrepancy.

Hope this helps, let me know if you manage to close the gap

Stefano

Note: I adopted the units and axes of your plot, to try to reproduce it, but you have some inconsistencies:

  • fluence is given as a differential over energy, so you should have 1/E in the units
  • usually when making the x axis logarithmic and showing a differential, it is customary to show lethargy units, i.e. multiply the function by E (just multiply the bin values by the bin centers)