Mismatch of the CSDA energy loss depositions using the printed dE/dx values and FLUKA simulations

Physics, Transport and Magnetic Fields
1

Dear FLUKA experts,

I tried to simulate carbon ions penetrated in a water phantom.
My aim is to check the consistency between the printed dE/dx values from the DELTARAY card and the depth-energy relation curves from FLUKA simulation.

I used the following settings in the .inp file:

STEPSIZE is fixed to 0.1cm (same as the voxel size in the z-direction)
IONFLUCT is switched off (no energy fluctuation for hadrons & electrons in a step)
EMFCUT for electrons transport is 10 GeV (electrons deposit energy locally)
MCSTHRES is switched off (no MCS, only along the z-direction)
THRESHOL is set to a very large value (no inelastic and elastic nuclear interactions)
PAIRBREM, inhibit both

In this case, to my understand, only CSDA and ionization (with delta electrons) happen, and carbons only move along the z-axis.

Here I performed a simulation of 100000 0.1 GeV carbon ions in a water phantom and plotted the depth-energy curves of carbon and electron, with the resolution of 0.1cm in the z-direction (same as the stepsize).

Second, I used the printed dE/dx values (from the DELTARAY card) and made interpolations to construct the restricted stopping power function L(E).

Third, I set up my own calculations, with initial energy 1200 MeV, the L(E) above, stepsize=1mm. The energy loss of CSDA in each step E_loss was calculated by Eq. (D2) of Ref. Kawrakow, Medical Physics 27 (2000) 485, which is proved to be precise enough.

image

Note that the energy loss of each step in my own calculation is actually

E_total = E_loss + E_electron,

where the latter is extracted from the electrons’ depth-energy curve of step 1.

However, results showed that my own calculations of the energy loss in each step (using the printed dE/dx values) is in general smaller than the FLUKA simulations. I’m not sure whether my settings in the .inp file are correct. My aim is just to check the CSDA and ionization processes without any other physics processes. I have attached the .inp file and my own calculations (using Mathematica).

My own calculations, columns from left to right are
Depth(cm), current E (own), E_loss (own), current E (FLUKA), E_loss (FLUKA), E_electron (FLUKA)

image
image602Ă—644 21.1 KB

image

Carbon_test.inp (1.8 KB)

3

Hello Kaiwen,

Sorry for the shamefully late reply. Almost a year delay. I am now certainly confirmed in the hall of shame of this forum.

My aim is to check the consistency between the printed dE/dx values from the DELTARAY card and the depth-energy relation curves from FLUKA simulation.

Understood. Let’s proceed.

I used the following settings in the .inp file:

STEPSIZE is fixed to 0.1cm (same as the voxel size in the z-direction)
IONFLUCT is switched off (no energy fluctuation for hadrons & electrons in a step)
EMFCUT for electrons transport is 10 GeV (electrons deposit energy locally)
MCSTHRES is switched off (no MCS, only along the z-direction)
THRESHOL is set to a very large value (no inelastic and elastic nuclear interactions)
PAIRBREM, inhibit both

I have two observations concerning your course of action:

  • In the input file, EMFCUT is commented out and hence inactive. I henceforth assume this was a mild hiccup and that you had it instead activated in your exercise.

  • You are not passing a DELTARAY card governing the production of delta rays by the primary ion, which is then ruled by default settings - implicitly requested by your DEFAULTS card with SDUM=PRECISIOn, which reigns supreme, and implies that your 12C are still producing delta rays with energies of 100 keV and above. While your continuous energy losses are fluctuation-free as per the IONFLUCT card you pass, the energy of the produced delta ray is not. That is, not all delta rays have the same energy - it is sampled from a distribution! Maybe the overall effect is minor, but strictly speaking you are not doing a continuous slowing down approximation. A minor note in passing: since you have a discrete interaction mechanism active (delta-ray production) with potentially short steps, it may well be that your STEPSIZE card was not as effective as you think (but i did not cross-check this further as this is not central to your post).

In any event, for your purposes what you wanted to check is whether when you put the simulation in conditions sufficiently close to the continuous slowing down approximation (CSDA).

I have the impression that your binning in z / your step sizes of 0.1 cm are too coarse to accurately check the agreement. The dE/dx curve varies very happily in a log scale (a lot).

I repeated your exercise (with a DELTARAY card blocking explicit production of delta rays by 12C, and now with 400 bins in z from 0 to 4 cm - and accordingly short STEPSIZEs) and obtained the following curve (adapting the energy so as to match precisely one of the FLUKA dE/dx tabulation points):
pl

The mild difference you see is because of the coarse integration i used (trapezoidal rule)… It is then sufficient for me to relax the stepsize to 0.1 cm to obtain a deviation:

pl_yourstepsize

So, in a nutshell, it’s a precision problem. Put thinner stepsizes and binnings, and you should be alright.

Here goes a tarball with the .flair and a quick and dirty fortran snippet to produce the curves. Normally it suffices to go into csda/ and type “make”. It’s not at all intended for production purposes, it’s just a quick snippet:

20230628_C_CSDA.tgz (3.9 MB)

I did (will) not do a cursory cross check of the expression you’re using to perform the integrals. Could it be that there is a further misinterpretation of its intended use? In any event, a trapezoidal rule suffices (at the cost of a thinner grid, sure).

Cesc

3 Likes