1 | \section[Ion Scattering]{Ion Scattering} |
---|
2 | |
---|
3 | The necessity of accurately computing the characteristics of interatomic |
---|
4 | scattering arises in many disciplines in which energetic ions pass |
---|
5 | through materials. Traditionally, solutions to this problem not involving |
---|
6 | hadronic interactions have been dominated by the multiple scattering, |
---|
7 | which |
---|
8 | is reasonably successful, but not very flexible. In particular, it |
---|
9 | is relatively difficult to introduce into such a system a particular |
---|
10 | screening function which has been measured for a specific atomic pair, |
---|
11 | rather than the universal functions which are applied. |
---|
12 | In many problems |
---|
13 | of current interest, such as the behavior of semiconductor device |
---|
14 | physics in a space environment, nuclear reactions, particle showers, |
---|
15 | and other effects are critically important in modeling the full details |
---|
16 | of ion transport. The process $G4ScreenedNuclearRecoil$ provides |
---|
17 | simulation of ion elastic scattering \cite{ion_scat_model}. |
---|
18 | |
---|
19 | \subsection{Method} |
---|
20 | |
---|
21 | The method used in this computation is a variant of a subset of the |
---|
22 | method described in Ref.\cite{MendenhallWellerXSection}. |
---|
23 | A very short recap of the basic material is included here. The scattering |
---|
24 | of two atoms from each other is assumed to be a completely classical |
---|
25 | process, subject to an interatomic potential described by a potential |
---|
26 | function\begin{equation} |
---|
27 | V(r)=\frac{Z_{1}Z_{2}e^{2}}{r}\phi\left(\frac{r}{a}\right)\label{VR_eqn}\end{equation} |
---|
28 | where $Z_{1}$ and $Z_{2}$ are the nuclear proton numbers, $e^{2}$ |
---|
29 | is the electromagnetic coupling constant ($q_{e}^{2}/4\pi\epsilon_{0}$ |
---|
30 | in SI units), $r$ is the inter-nuclear separation, $\phi$ is the |
---|
31 | screening function describing the effect of electronic screening of |
---|
32 | the bare nuclear charges, and $a$ is a characteristic length scale |
---|
33 | for this screening. In most cases, $\phi$ is a universal function |
---|
34 | used for all ion pairs, and the value of $a$ is an appropriately |
---|
35 | adjusted length to give reasonably accurate scattering behavior. In |
---|
36 | the method described here, there is no particular need for a universal |
---|
37 | function $\phi$, since the method is capable of directly solving |
---|
38 | the problem for most physically plausible screening functions. It |
---|
39 | is still useful to define a typical screening length $a$ in the calculation |
---|
40 | described below, to keep the equations in a form directly comparable |
---|
41 | with our previous work even though, in the end, the actual value is |
---|
42 | irrelevant as long as the final function $\phi(r)$ is correct. From |
---|
43 | this potential $V(r)$ one can then compute the classical scattering |
---|
44 | angle from the reduced center-of-mass energy $\varepsilon\equiv E_{c}a/Z_{1}Z_{2}e^{2}$ |
---|
45 | (where $E_{c}$ is the kinetic energy in the center-of-mass frame) |
---|
46 | and reduced impact parameter $\beta\equiv b/a$\begin{equation} |
---|
47 | \theta_{c}=\pi-2\beta\int_{x_{{\scriptscriptstyle 0}}}^{\infty}f(z)\, dz/z^{2}\label{theta_eqn}\end{equation} |
---|
48 | where\begin{equation} |
---|
49 | f(z)=\left(1-\frac{\phi(z)}{z\,\varepsilon}-\frac{\beta^{2}}{z^{2}}\right)^{-1/2}\label{fz_eqn}\end{equation} |
---|
50 | and $x_{{\scriptscriptstyle 0}}$ is the reduced classical turning |
---|
51 | radius for the given $\varepsilon$ and $\beta$. |
---|
52 | |
---|
53 | The problem, then, is reduced to the efficient computation of this |
---|
54 | scattering integral. In our previous work, a great deal of analytical |
---|
55 | effort was included to proceed from the scattering integral to a full |
---|
56 | differential cross section calculation, but for application in a Monte-Carlo |
---|
57 | code, the scattering integral $\theta_{c}(Z_{1},\, Z_{2},\, E_{c},\, b)$ |
---|
58 | and an estimated total cross section $\sigma_{{\scriptscriptstyle 0}}(Z_{1},\, Z_{2},\, E_{c})$ |
---|
59 | are all that is needed. Thus, we can skip algorithmically forward |
---|
60 | in the original paper to equations 15-18 and the surrounding discussion |
---|
61 | to compute the reduced distance of closest approach $x_{{\scriptscriptstyle 0}}$. |
---|
62 | This computation follows that in the previous work exactly, and will |
---|
63 | not be reintroduced here. |
---|
64 | |
---|
65 | For the sake of ultimate accuracy in this algorithm, and due to the |
---|
66 | relatively low computational cost of so doing, we compute the actual |
---|
67 | scattering integral (as described in equations 19-21 of \cite{MendenhallWellerXSection}) |
---|
68 | using a Lobatto quadrature of order 6, instead of the 4th order method |
---|
69 | previously described. This results in the integration accuracy exceeding |
---|
70 | that of any available interatomic potentials in the range of energies |
---|
71 | above those at which molecular structure effects dominate, and should |
---|
72 | allow for future improvements in that area. The integral $\alpha$ |
---|
73 | then becomes (following the notation of the previous paper) |
---|
74 | \begin{equation} |
---|
75 | \alpha\approx\frac{1+\lambda_{{\scriptscriptstyle 0}}}{30}+\sum_{i=1}^{4}w'_{i}\, f\left(\frac{x_{{\scriptscriptstyle 0}}}{q_{i}}\right)\label{alpha_eq}\end{equation} |
---|
76 | where |
---|
77 | \begin{equation} |
---|
78 | \lambda_{{\scriptscriptstyle 0}}=\left(\frac{1}{2}+\frac{\beta^{2}}{2\, x_{{\scriptscriptstyle 0}}^{2}}-\frac{\phi'(x_{{\scriptscriptstyle 0}})}{2\,\varepsilon}\right)^{-1/2}\label{lambda_eqn}\end{equation} |
---|
79 | \\ |
---|
80 | $w'_{i}\in${[}0.03472124, 0.1476903, 0.23485003, 0.1860249] \\ |
---|
81 | $q_{i}\in${[}0.9830235, 0.8465224, 0.5323531, 0.18347974]\\ |
---|
82 | Then \\ |
---|
83 | \begin{equation} |
---|
84 | \theta_{c}=\pi-\frac{\pi\beta\alpha}{x_{{\scriptscriptstyle 0}}}\label{thetac_alpha_eq}\end{equation} |
---|
85 | |
---|
86 | |
---|
87 | The other quantity required to implement a scattering process |
---|
88 | is the total scattering cross section $\sigma_{{\scriptscriptstyle 0}}$ |
---|
89 | for a given incident ion and a material through which the ion is propagating. |
---|
90 | This value requires special consideration for a process such as screened |
---|
91 | scattering. In the limiting case that the screening function is unity, |
---|
92 | which corresponds to Rutherford scattering, the total cross section |
---|
93 | is infinite. For various screening functions, the total cross section |
---|
94 | may or may not be finite. However, one must ask what the intent of |
---|
95 | defining a total cross section is, and determine from that how to |
---|
96 | define it. |
---|
97 | |
---|
98 | In Geant4, the total cross section is used to determine a mean-free-path |
---|
99 | $l_{\mu}$ which is used in turn to generate random transport distances |
---|
100 | between discrete scattering events for a particle. In reality, where |
---|
101 | an ion is propagating through, for example, a solid material, scattering |
---|
102 | is not a discrete process but is continuous. However, it is a useful, |
---|
103 | and highly accurate, simplification to reduce such scattering to a |
---|
104 | series of discrete events, by defining some minimum energy transfer |
---|
105 | of interest, and setting the mean free path to be the path over which |
---|
106 | statistically one such minimal transfer has occurred. This approach |
---|
107 | is identical to the approach developed for the original TRIM code |
---|
108 | \cite{TRIM1980}. As long as the minimal interesting energy transfer |
---|
109 | is set small enough that the cumulative effect of all transfers smaller |
---|
110 | than that is negligible, the approximation is valid. As long as the |
---|
111 | impact parameter selection is adjusted to be consistent with the selected |
---|
112 | value of $l_{\mu}$, the physical result isn't particularly sensitive |
---|
113 | to the value chosen. |
---|
114 | |
---|
115 | Noting, then, that the actual physical result isn't very sensitive |
---|
116 | to the selection of $l_{\mu},$ one can be relatively free about defining |
---|
117 | the cross section $\sigma_{{\scriptscriptstyle 0}}$ from which $l_{\mu}$ |
---|
118 | is computed. The choice used for this implementation is fairly simple. |
---|
119 | Define a physical cutoff energy $E_{min}$ which is the smallest energy |
---|
120 | transfer to be included in the calculation. Then, for a given incident |
---|
121 | particle with atomic number $Z_{1}$, mass $m_{1}$, and lab energy |
---|
122 | $E_{inc}$, and a target atom with atomic number $Z_{2}$ and mass |
---|
123 | $m_{2}$, compute the scattering angle $\theta_{c}$ which will transfer |
---|
124 | this much energy to the target from the solution of\begin{equation} |
---|
125 | E_{min}=E_{inc}\,\frac{4\, m_{1}\, m_{2}}{(m_{1}+m_{2})^{2}}\,\sin^{2}\frac{\theta_{c}}{2}\label{recoil_kinematics}\end{equation}. |
---|
126 | Then, noting that $\alpha$ from eq. \ref{alpha_eq} is a number very close to unity, one can solve for an approximate |
---|
127 | impact parameter $b$ with a single root-finding operation to find the classical turning point. |
---|
128 | Then, define the total cross section to be $\sigma_{{\scriptscriptstyle 0}}=\pi b^{2}$, |
---|
129 | the area of the disk inside of which the passage of an ion will cause |
---|
130 | at least the minimum interesting energy transfer. Because this process |
---|
131 | is relatively expensive, and the result is needed extremely frequently, |
---|
132 | the values of $\sigma_{{\scriptscriptstyle 0}}(E_{inc})$ are precomputed |
---|
133 | for each pairing of incident ion and target atom, and the results |
---|
134 | cached in a cubic-spline interpolation table. However, since the actual result isn't very critical, the |
---|
135 | cached results can be stored in a very coarsely sampled table without |
---|
136 | degrading the calculation at all, as long as the values of the $l_{\mu}$ |
---|
137 | used in the impact parameter selection are rigorously consistent with |
---|
138 | this table. |
---|
139 | |
---|
140 | The final necessary piece of the scattering integral calculation is |
---|
141 | the statistical selection of the impact parameter $b$ to be used |
---|
142 | in each scattering event. This selection is done following the original |
---|
143 | algorithm from TRIM, where the cumulative probability distribution |
---|
144 | for impact parameters is \begin{equation} |
---|
145 | P(b)=1-\exp\left(\frac{-\pi\, b^{2}}{\sigma_{{\scriptscriptstyle 0}}}\right)\label{cum_prob}\end{equation} |
---|
146 | where $N\,\sigma_{{\scriptscriptstyle 0}}\equiv1/l_{\mu}$ where $N$ |
---|
147 | is the total number density of scattering centers in the target material |
---|
148 | and $l_{\mu}$ is the mean free path computed in the conventional |
---|
149 | way. To produce this distribution from a uniform random variate $r$ |
---|
150 | on (0,1], the necessary function is\begin{equation} |
---|
151 | b=\sqrt{\frac{-\log r}{\pi\, N\, l_{\mu}}}\label{sampling_func}\end{equation} |
---|
152 | This choice of sampling function does have the one peculiarity that |
---|
153 | it can produce values of the impact parameter which are larger than |
---|
154 | the impact parameter which results in the cutoff energy transfer, |
---|
155 | as discussed above in the section on the total cross section, with |
---|
156 | probability $1/e$. When this occurs, the scattering event is not |
---|
157 | processed further, since the energy transfer is below threshold. For |
---|
158 | this reason, impact parameter selection is carried out very early |
---|
159 | in the algorithm, so the effort spent on uninteresting events is minimized. |
---|
160 | |
---|
161 | The above choice of impact sampling is modified when the mean-free-path |
---|
162 | is very short. If $\sigma_{{\scriptscriptstyle 0}}>\pi\left(\frac{l}{2}\right)^{2}$ |
---|
163 | where $l$ is the approximate lattice constant of the material, as |
---|
164 | defined by $l=N^{-1/3}$, the sampling is replaced by uniform sampling |
---|
165 | on a disk of radius $l/2$, so that\begin{equation} |
---|
166 | b=\frac{l}{2}\sqrt{r}\label{flat_sampling}\end{equation} |
---|
167 | This takes into account that impact parameters larger than half the |
---|
168 | lattice spacing do not occur, since then one is closer to the adjacent |
---|
169 | atom. This also derives from TRIM. |
---|
170 | |
---|
171 | One extra feature is included in our model, to accelerate the production |
---|
172 | of relatively rare events such as high-angle scattering. This feature |
---|
173 | is a cross-section scaling algorithm, which allows the user access |
---|
174 | to an unphysical control of the algorithm which arbitrarily scales |
---|
175 | the cross-sections for a selected fraction of interactions. This is |
---|
176 | implemented as a two-parameter adjustment to the central algorithm. |
---|
177 | The first parameter is a selection frequency $f_{h}$ which sets what |
---|
178 | fraction of the interactions will be modified. The second parameter |
---|
179 | is the scaling factor for the cross-section. This is implemented by, |
---|
180 | for a fraction $f_{h}$ of interactions, scaling the impact parameter |
---|
181 | by $b'=b/\sqrt{scale}$. This feature, if used with care so that it |
---|
182 | does not provide excess multiple-scattering, can provide between 10 |
---|
183 | and 100-fold improvements to event rates. If used without checking |
---|
184 | the validity by comparing to un-adjusted scattering computations, |
---|
185 | it can also provide utter nonsense. |
---|
186 | |
---|
187 | |
---|
188 | \subsection{Implementation Details} |
---|
189 | |
---|
190 | \label{Lobatto}The coefficients for the summation to approximate |
---|
191 | the integral for $\alpha$ in eq.(\ref{alpha_eq}) are derived from |
---|
192 | the values in Abramowitz \& Stegun \cite{AbramowitzStegunLobatto}, |
---|
193 | altered to make the change-of-variable used for this integral. There |
---|
194 | are two basic steps to the transformation. First, since the provided |
---|
195 | abscissas $x_{i}$ and weights $w_{i}$ are for integration on {[}-1,1], |
---|
196 | with only one half of the values provided, and in this work the integration |
---|
197 | is being carried out on {[}0,1], the abscissas are transformed as:\begin{equation} |
---|
198 | y_{i}\in\left\{ \frac{1\mp x_{i}}{2}\right\} \label{as_yxform}\end{equation} |
---|
199 | Then, the primary change-of-variable is applied resulting in:\begin{eqnarray} |
---|
200 | q_{i} & = & \cos\frac{\pi\, y_{i}}{2}\label{q_xform}\\ |
---|
201 | w'_{i} & = & \frac{w_{i}}{2}\sin\frac{\pi\, y_{i}}{2}\label{wprime_xform}\end{eqnarray} |
---|
202 | except for the first coefficient $w'_{1}$where the $\sin()$ part |
---|
203 | of the weight is taken into the limit of $\lambda_{{\scriptscriptstyle 0}}$ |
---|
204 | as described in eq.(\ref{lambda_eqn}). This value is just $w'_{1}=w_{1}/2$. |
---|
205 | |
---|
206 | \subsection{Status of this document} |
---|
207 | 06.12.07 created by V. Ivanchenko from paper of M.H.~Mendenhall and R.A.~Weller\\ |
---|
208 | 06.12.07 further edited by M. Mendenhall to bring contents of paper up-to-date with current implementation. |
---|
209 | |
---|
210 | \begin{latexonly} |
---|
211 | |
---|
212 | \begin{thebibliography}{99} |
---|
213 | |
---|
214 | \bibitem{ion_scat_model} |
---|
215 | M.H.~Mendenhall, R.A.~Weller, An algorithm for computing screened Coulomb |
---|
216 | scattering in Geant4, {\em Nucl. Instr. Meth. B 227 (2005) 420.} |
---|
217 | |
---|
218 | \bibitem{MendenhallWellerXSection} |
---|
219 | M.H.~Mendenhall, R.A.~Weller, Algorithms for the rapid computation of |
---|
220 | classical cross sections for screened coulomb collisions, |
---|
221 | {\em Nucl. Instr. Meth. in Physics Res. B58 (1991) 11.} |
---|
222 | |
---|
223 | \bibitem{TRIM1980} |
---|
224 | J.P.~Biersack, L.G.~Haggmark, A Monte Carlo computer program for the |
---|
225 | transport of energetic ions in amorphous targets, |
---|
226 | {\em Nucl. Instr. Meth. in Physics Res. 174 (1980) 257.} |
---|
227 | |
---|
228 | \bibitem{AbramowitzStegunLobatto} |
---|
229 | M.~Abramowitz, I.~Stegun (Eds.), Handbook of Mathematical Functions, |
---|
230 | {\em Dover, New York, 1965, pp. 888, 920.} |
---|
231 | |
---|
232 | \end{thebibliography} |
---|
233 | |
---|
234 | \end{latexonly} |
---|
235 | |
---|
236 | \begin{htmlonly} |
---|
237 | |
---|
238 | \subsection{Bibliography} |
---|
239 | |
---|
240 | \begin{enumerate} |
---|
241 | |
---|
242 | \item |
---|
243 | M.H.~Mendenhall, R.A.~Weller, An algorithm for computing screened Coulomb |
---|
244 | scattering in Geant4, {\em Nucl. Instr. Meth. B 227 (2005) 420.} |
---|
245 | |
---|
246 | \item |
---|
247 | M.H.~Mendenhall, R.A.~Weller, Algorithms for the rapid computation of |
---|
248 | classical cross sections for screened coulomb collisions, |
---|
249 | {\em Nucl. Instr. Meth. in Physics Res. B58 (1991) 11.} |
---|
250 | |
---|
251 | \item |
---|
252 | J.P.~Biersack, L.G.~Haggmark, A Monte Carlo computer program for the |
---|
253 | transport of energetic ions in amorphous targets, |
---|
254 | {\em Nucl. Instr. Meth. in Physics Res. 174 (1980) 257.} |
---|
255 | |
---|
256 | \item |
---|
257 | M.~Abramowitz, I.~Stegun (Eds.), Handbook of Mathematical Functions, |
---|
258 | {\em Dover, New York, 1965, pp. 888, 920.} |
---|
259 | |
---|
260 | \end{enumerate} |
---|
261 | |
---|
262 | \end{htmlonly} |
---|