Created attachment 535 [details] Patch implementing the described improvements. See also https://github.com/Geant4/geant4/pull/5 [I have also submitted this as https://github.com/Geant4/geant4/pull/5 because I don't know which method is preferred.] Added calculation of the Fermi density effect using atomic data as described in R.M. Sternheimer et al. "Density Effect For The Ionization Loss of Charged Particles in Various Substances" Atom. Data Nucl. Data Tabl. 30 (1984) 261-271. This calculation gives values for the density effect 'delta' which are much more accurate than that provided by the Sternheimer 3-part approximation, introduced in his 1952 paper, Phys.Rev. 88 851-859. This is particularly true when using the general method for producing parameter values given in his 1971 paper, Phys.Rev. B3 3681-3692, but doing the full calculation is still a significant improvement over using the tabulated values given in the 1984 paper for elements and a selected set of mixtures and compounds. The purpose of the 3-part approximation was to avoid heavy computation from the perspective of sliderules (1952) or minicomputers (1971). It is no longer necessary. I ran tests using the NOvA simulation framework and could measure no difference in performance between using the approximation and the full treatment. A major advantage of this method is that the stopping power of a substance does not abruptly change when the user modifies it so that it is no longer a tabulated material. In our case, we noticed this effect when we changed iron to steel (98% iron, 2% other). Previously, iron got the 1984 tabulated values while steel got the 1971 general form, which gives up to a 1.3% difference in mean muon range, with the largest difference at about 200MeV. I have tested the code with muons in iron and verified that the mean range changes as expected between the approximate and exact forms. Implementation notes: The calculation involves finding roots and can fail. In this case, the code falls back to the approximation and prints a message saying it has done so. I have not found any cases that prompt a failure, but have only tested with materials found in and around the NOvA detector. Attempted to make clear in places that refer to the parameterization, such as G4IonisParamMat::SetDensityEffectParameters(), that these calls only have any effect if the full calculation fails. I found a separate implementation of the Sternheimer parameterization in source/processes/electromagnetic/highenergy/src/G4mplIonisationModel.cc I don't know why this parallel implementation exists, so I simply changed a comment to make it more clear that it was the approximate form. G4IonisParamMat::DensityCorrection() is no longer inlined. Removed debugging print statements in G4IonisParamMat::ComputeDensityEffect().
Dear Matthew, Sorry, but Geant4 is using this ( Atom. Data Nucl. Data Tabl. 30 (1984) 261-271) parametrization since 2009 - see G4DensityEffectData class. VI
Hi Vladimir, Please take a closer look. My patch does not implement the parameterization from the 1984 paper. It implements the full calculation that the paramterization approximates. That is, it is not a way to find the values {x0, x1, k, C, delta0}, but rather it calculates the exact value of delta at any energy without using those parameters at all. Thanks, Matthew
Created attachment 537 [details] Improved patch This improved patch fixes two bugs in the previous version and gives a more informative error message if the calculation fails. * Zero out the sternf array before using it. Otherwise we're adding to uninitialized values. * Don't divide by Z when counting the number of electrons. This is already taken into account so was leading to wrong values (sometimes subtle, sometimes obvious -- my earlier test cases were in the subtle category). * Give a better message when calculation fails. One failure mode that can be reproduced by hand is if the density of a substance is arbitrarily increased. This raises the plasma energy until Sternheimer 1984 Eq (8) has no solution. Pure hydrogen at several g/cc is such a case.
Created attachment 556 [details] 3rd version of this patch, with conditional warnings Wrap the warnings that the new exact calculation of the Fermi density effect can print in checks for the value of G4NistManager::Instance()->GetVerbose() so that it is possible for the user to suppress them.
The first implementation of this proposal is included in the 10.5ref08. We can close this ticket.
