Stopping Power Calculations

Critical Information:

As of RUMP versions linked on 5/18/99 and later, the default stopping powers for H and He on Si and C have been changed. The new stopping powers are taken from the measurements by G. Konac, S. Kalbitzer, Ch. Klatt, D. Niemann, and R. Stoll, NIM B 136-138 (1998), 159-165. The deviations from the default Ziegler tables are significant (7-10%).

Background:

The original DOS version of RUMP used the 6th order polynomial fit to the stopping power for hydrogen and helium in elements as given in Backscattering Spectrometry, W.K. Chu, J.W. Mayer and M.-A. Nicolet (Eds.), Academic Press, NY, 1978. For backward compatibility purposes, these have been retained as the data files d_35.stp, he4_35.stp, and h_35.stp. These files were limited in use to below 3.5 MeV (as per their name). Several users have seen and used a second file high_e.dat which contained similar data for use at higher energies. The use of these files is now obsolete and should be avoided unless there is a particular need to match previous spectra simulations.

When RUMP converted from FORTRAN to C, one of the changes was to incorporate automatic generation of stopping power tables for all elements in all targets - choosing the energy range to match the needs of the simulation. The energy range of each table is determined at run time based on the beam energy being simulated. Multiple tables for different energy ranges, as well as multiple ion (for ERD), were added. Generation of the tables was based on an analytical model from J.F. Ziegler, J.P. Biersack, U. Littmark, Empirical Stopping Powers for Ions in Solids, IBM Research Report RC9250, 1982 and The Stopping and Range of Ions in Solids, Pergamon Press, New York, 1985 by the same authors..

The generation of stopping power tables happens outside of the user's view, except for several information messages. These include the creation of a new table (a table corresponds to a particle incident ion with entries for all possible target Zs):
Created new stopping power table: 2 4.000000 0.360000 3.450000
signifies a new table for 4He, covering the energy range 0.36 to 3.45 MeV. As target nuclei are required, the table is filled in.
Fitting particle 4He for Z = 14 ... zsf ... maximum deviation: 0.76%
This message specifies that the entry for Si (Z=14) was computed into one of the 4He stopping power tables (multiple tables may exist if the beam energy is changed significantly in different simulations).

Internally, RUMP continues to use a polynomial representation for the stopping power. This is necessary to efficiently implement Bragg's Law for complex target structures. The analytical models, however, are neither computationally efficient nor amenable for direct implementation of Bragg's Law. Consequently, RUMP "fits" a polynomial to the analytical models on the fly. Since the energy window changes between tables, the parameterization is always optimized for the energy range of interest. The maximum deviation (0.76% above) quantitatively shows how well the polynomial representation matches the analytical representation (neither of which is necessarily correct). Although this is an expensive computation, it need only be performed once for each ion/target pair unless the beam energy is significantly changed (either up or down).

Improving Fit Precision:

Chu, Mayer and Nicolet used a 6th order polynomial fit to the stopping power, based on linear energy. It turns out to be very difficult in this representation to handle the curvature of the stopping power function across the Bragg peak. However, representing the stopping power as polynomial in the square root of the energy maintains both the simple additivity of the polynomial expression which improving the fit precision dramatically. There are several better forms, but none can retain the computational simplicity for both evaluation and Bragg's Law. The choice of representing stopping powers in linear or root energy may be selected with a configuration command.

config stop_type linear
config stop_type sqrt

Switching from linear to square root mode will improve the precision in several ways. First, the lower energy limit in the stopping power tables is dropped by a factor of two (reflecting the fact that it is easier to match across the Bragg peak). Second, the fit is simply better. For Si, the example above would create a table extending from 0.18 MeV to 3.45 MeV and the maximum deviation (over this extended range) decreases to 0.24%.

Recommendation 1: Add the line "config stop_type sqrt return" to your RUMP.INI file to automatically enable this mode by default.

However, while the precision of fit to the analytical model is improved, the fit is still limited in accuracy to the accuracy of the underlying analytical expressions (i.e.. Ziegler model).

Incorporating newer stopping power measurements:

It is generally known (among quantitative RBS groups) that the default Ziegler stopping power for He on Si in RUMP must be scaled by 89-100% to match careful experimental values [A. Climent-Font, U. Watjen, H. Bax, NIM B 71 (1992), 81 and A. Climent-Font, M.T. Fernandez-Jimenez, U Watjen, J. Perriere, NIM A 353 (1994) 575]. More recent papers have attempted to provide quantitative measurements of several specific systems, most notably Kalbitzer's group for various ions in Si.

Stopping power measurements of 1H, 4He and 14 N in Si in the energy range of 0.02 - 1 MeV/amu, D. Niemann, G. Konac, S. Kalbitzer, NIM B 118 (1996), 11-18
Energy loss and straggling of H and He ions of keV energies in Si and C, G. Konac, S. Kalbitzer, Ch. Klatt, D. Niemann, and R. Stoll, NIM B 136-138 (1998), 159-165].

The problems with the RUMP default values were clearly pointed out in a very recent paper,

Improved stoichiometry measurement using 4He elastic backscattering: experiment and simulation, W.N. Lennard, G.R. Massoumi, T.W. Simpson and I.V. Mitchell, NIM B 152 (1999), 370-376.

Since RUMP's fitting algorithm really does not care the analytical source of the stopping power data, it was a simple matter to add the formalism from Konac's paper into an additional stopping power routine. The form is also sufficiently direct that it will permit additional ion/target pairs to be added at the user convenience. If an entry exists in the data file for the new algorithm, it will replace the Ziegler stopping subroutine. The values from the Niemann et al. and Konac et al. papers have been included by default.

The electronic component to the stopping power is expressed as a complex inverse polynomial:

S_e(E) = E^s ln (e+BE) / P(E)
P(E) = a₀ + a₁ E^1/4 + a₂ E^1/2 + a₃ E^3/4 + a₄ E + a₅ E^(1+s)

Note that this form encompasses data from both papers. For Konac et al., a₃ and a₄ are zero, while in Niemann et al.'s formalism, s = 1/2. These coefficients are maintained in a data file .../rump/data/konac.dat, the default version of which (on May 18, 1999) is shown below. In this format, stopping powers for high Z elements can also be easily incorporate - simply fit the data to the format specified.

Recommendation 3: In publications, specify the source of the stopping power data with the reference to RUMP if appropriate.

Hierarchy of stopping power calculations:

The analytical value fo the stopping power for any ion/target pair follows the following search and modification procedure:

If a user mystop.dll is loaded, attempt to determine stopping power from the routine. If not,
Look for the ion/target pair in the konac.dat tables, and calculate using those parameters. If not,
Default to using the Ziegler subroutine for all ions in all targets.

The exception stopping power table (../rump/data/except.stp) file is then scanned for this ion/target pair. If a modification expression is defined, the values from the analytical result are modified by the equation provided (The analytical stopping power is given by the variable zgl independent of how is it obtained above.) The modifed stopping powers are then fit to a polynomial form using either linear or square root of energy as the independent variable.

Reverting to old (pre 5/99) stopping power models:

If it is necessary to revert to the old stopping power models (highly discouraged), there are two options.

Delete the file .../rump/data/konac.dat. This will print an error message to warn of its absence (reminding that you are purposely using less accurate data) and not load any parameters.
Comment out the individual lines of relevance in the ../rump/data/konac.dat file.

Example konac.dat file:

# ==================================================================================================
# Stopping power fits on specific incident/target ion pairs
# This file contains data for the Konac form of the stopping power description
# for all ions. The entries in this field will be used before resorting to
# the default Ziegler stopping power calculations (scaled proton data). To
# ==================================================================================================
# Stopping power fit coefficients
#
# Primary Reference:
#   Konac 98: G. Konac, S. Kalbitzer, Ch. Klatt, D. Niemann, R. Stoll,
#              NIM B, 136-138 (1998) 159-165.
#
# Secondary data references:
#   Konac 98a: Konac 98, assuming velocity scaling from 1H and 4He
#
#===================================================================================================
# Inc Target              Coefficients as per Konac 98 reference         (a3)
# Z M    Z       s       a0        a1        a2        a3        a4        a5      beta
#===================================================================================================
1 1    6     0.39    3.40E-2 -3.93E-2   2.09E-2       0         0     7.89E-1    36.5 # Konac 98
1 2    6     0.39    3.40E-2 -3.93E-2   2.09E-2       0         0     7.89E-1    36.5 # Konac 98a
2 4    6     0.49    1.36E-3 2.47E-2 -1.81E-2       0         0     1.88E-1    36.5 # Konac 98
2 3    6     0.49    1.36E-3 2.47E-2 -1.81E-2       0         0     1.88E-1    36.5 # Konac 98a
1 1   14     0.37    4.16E-2 -1.47E-1   1.80E-1       0         0     2.79E-1    15.7 # Konac 98
1 2   14     0.37    3.66E-2 -1.24E-1   1.52E-1       0         0     2.93E-1    15.7 # Konac 98
2 3   14     0.52    8.24E-3 -2.02E-2   2.28E-2       0         0     8.18E-2    15.7 # Konac 98
2 4   14     0.52    7.93E-3 -1.93E-2   2.21E-2       0         0     8.25E-2    15.7 # Konac 98

Contact Mike Thompson if you want to discuss these issues further.