Monte Carlo calculation of energy loss of hydrogen and helium ions transmitted under channelling conditions in silicon single crystal

Nuclear Instruments and Methods in Physics Research B 268 (2010) 1361–1366

Contents lists available at ScienceDirect

Nuclear Instruments and Methods in Physics Research B journal homepage: www.elsevier.com/locate/nimb

Monte Carlo calculation of energy loss of hydrogen and helium ions transmitted under channelling conditions in silicon single crystal O. El Bounagui *, H. Erramli Faculty of Sciences Semlalia, University Cadi Ayyad, Marrakech, Morocco

a r t i c l e

i n f o

Article history: Received 10 June 2009 Received in revised form 7 January 2010 Available online 18 January 2010 Keywords: Channelling Electronic energy loss Simulations code Silicon lattice

a b s t r a c t In this work, we report on calculations of the electronic channelling energy loss of hydrogen and helium ions along Sih1 0 0i and Sih1 1 0i axial directions for the low energy range by using the Monte Carlo simulation code. Simulated and experimental data are compared for protons and He ions in the h1 0 0i and h1 1 0i axis of silicon. A reasonable agreement was found. Computer simulation was also employed to study the angular dependence of energy loss for 0.5, 0.8, 1, and 2 MeV channelled 4He ions transmitted through a silicon crystal of 3 lm thickness along the h1 0 0i axis. Ó 2010 Elsevier B.V. All rights reserved.

1. Introduction Accurate values of stopping powers are important not only for depth profiling of elements in materials [1] but also for the determination of the deposited energy and defect production [2] by the bombarding ion. Experiments that involve ions channelling along the major directions of a crystalline target require an even more detailed description of the energy loss mechanisms as a function of the impact parameter. A full treatment including all basic energy loss processes in a large energy range is still missing although the first steps have been made for simple systems [3–7]. Hence, empirical and semi-empirical procedures have mostly been adopted by the ion-beam community in order to estimate the energy loss for practical applications. In the early 1960s, it had been observed [8–11] that the energy loss of positively charged particles incident on crystalline materials at an angle close to a low index axis or plane consists of only a fraction of the energy loss in a random direction of incidence. This was explained by Lindhard [12] with the low charge density in the channel, leading to a lower energy loss. According to this theory, a value of 0.5 is expected for E, the ratio of the energy loss for channelled versus randomly incident particles. Values for this ratio reported for many crystals in a variety of experiments, range from 0.3 to 0.7. Detailed calculations of channelling effects usually require a Monte Carlo simulation of trajectories of an ion in a crystal. These * Corresponding author. Tel.: +212 60402741. E-mail addresses: [email protected], [email protected] (O. El Bounagui). 0168-583X/$ - see front matter Ó 2010 Elsevier B.V. All rights reserved. doi:10.1016/j.nimb.2010.01.010

calculations can be time consuming and difficult to interpret. On the other hand, continuum models give simple views on channelling that can lead to analytical results. In this picture the natural variable of importance is the transverse energy, and it is natural to interpret the behaviour of channelled ions as depending on the transverse energy. Several codes were developed, however each experiment has its own special features that demand extra effort to produce a modified version of the program. The main purpose of the present paper is to develop a faster, smaller, and more accurate new Monte Carlo code. This program permit the calculation of the channelling energy loss and taking into account the energy straggling due to inelastic collisions of ions moving under axial channelling conditions and thermal motion of crystal atoms which are not taken of in Ref. [13]. In this study, we calculate the electronic energy loss for a system of practical interest, namely 1H and 4He ions channelled along h1 0 0i and h1 1 0i directions of the silicon crystal. In the following, we shall give a brief description of the channelling method and the basic theory used to calculate the electronic energy loss under channelling conditions. In Section 4, we present our results for the calculation of the stopping power of different channelled light ions in a thin silicon crystal and compare the calculations to existing experimental data. The variations of the ratio of channelled to random stopping power for energy range 0.5–2.0 MeV for 4He ions along the silicon h1 0 0i axis as a function of the incidence angle was also studied. Good agreement has been found between our results and those obtained experimentally and theoretically by other works and with a high efficiency for our program.

1362

O. El Bounagui, H. Erramli / Nuclear Instruments and Methods in Physics Research B 268 (2010) 1361–1366

Our program is written in Fortran 90 under Windows, it runs at 180 ion collisions per second. This translates into 1 min of CPU time for about 104 ions passing through the target.

model that the nuclear charge of the atoms in a row (or plane) is uniformly averaged along the row (or plane). It is well known [12] that the critical channelling angle in the case of axial channelling is given by the relation

2. Theory

" " ##1=2 2 w Ca wc ðqÞ ¼ p1ffiffiffi Ln þ1 r min 2

2.1. Channelling

where:

Whenever energetic positive ions penetrate into a crystal lattice within a small angle wi < wc of a major crystallographic axis the ions are steered away from axis of atoms (Fig. 1) and channelling occurs. It arises from the gradually increasing Coulomb repulsive force between an incident particle and successive lattice nuclei. The main consequence of such steered motion is that the channelled ion never approaches closer to a lattice atom than a minimum distance of approach (rmin), where rmin is approximately equal to the mean thermal vibration amplitude (0.1 Å) of lattice atoms. For such a trajectory all collision processes such as large-angle scattering, nuclear reactions and inner-shell X-ray excitation are for the most part eliminated. Channelling exists for all positive ions, from protons to uranium ions, for energies in the keV until the 100 GeV regions. Channelling gives rise to a deep ‘‘tail” in implantation profiles of heavy ions, due to a reduced nuclear stopping power [8] and also to a reduced energy loss of light ions by a reduction in the electronic stopping power [14,15]. Channelling can be used for lattice site location studies, for the determination of thermal vibration amplitudes, for depth profiling of damage (for instance due to implantation) and annealing studies of damage, for the study of unilateral strain in epitaxially grown layers, and also to study fundamental aspects of the propagation of fast ions in matter. This summing-up contains the main applications, but is certainly not exhaustive.

w1 ¼

rin

ψc

Z1

ð3Þ

where U(r), is the interaction potential of the particle of the beam. Relation (3) gives the distance of minimum approach regardless of the specific potential that is used (e.g., Molière’s potential, Lindhard’s potential and so on). It is obvious that the results will differ according to the specific potential used but they do not differ markedly. The thermal smearing of the atom positions sets a lower limit to the minimum distance for which a row can provide the necessary correlated sequence of scatterings required for the channelling condition. The most useful first approximation to the critical angle is obtained by substituting rmin = q in Eq. (1), where q2 is twothirds of the mean square thermal vibration amplitude. 2.2.2. Electronic stopping powers In the energy range of interest, ions loss their energy essentially by interacting with electrons. This energy loss occurs through two types of collisions: long range collision with a low energy transfer (10 eV) and close collisions with high energy transfer (100 eV to 1 keV). The close collisions are directly proportional to the electrons density. For random irradiation, electrons density appears to be constant. On the other hand, well channelled ions have

d

x2 + r 2

ð2Þ

Uðr min Þ ¼ Ew2c

2.2.1. Critical angle Channelling is characterized by a critical angle, which is the maximum angle between ion and channel for a glancing collision to occur. Analytic models for ion channelling date back to the 1960’ [10]. Since then, several models have been published in the literature [12,16]. In these models formulas for calculating critical angles and minimum yields for axial and planar channels have been presented. Previously, these analytic models were widely used in the analysis of experimental channelling results. The trajectory of a channelled ion is such that the ion makes a glancing angle impact with the axes (axial channelling) or planes (planar channelling) of the crystal and is steered by small-angle scattering collisions at distances greater than 0.01 nm from the atomic cores. Since the steering of the channelled particle involves collisions with many atoms, one may consider in a continuum

r

 1=2 2Z 1 Z 2 e2 Ed

C2 is usually taken as equal to 3, a is the Thomas–Fermi screening distance, rmin is the minimum distance of approach, Z1 and Z2 are the charges of the beam particles and of the atoms of the nuclei of the target, respectively, e is the electron charge, E is the kinetic energy of the beam particles and d is the spacing of the atoms in the row. The critical channelling angle, as it has been established by Lindhard, means that if a particle is incident upon a row of atoms at an angle less than that specified by relation (1) the continuum model holds, which means that the potential due to the atomic row is considered as uniformly smeared out along the atomic row of the crystal. In the experimental works, a particle that emerges from a crystal is considered as channelled if the exit angle is less than that specified by relation (1). Thus the theoretical value of the critical distance of approach can be found if we use

2.2. Basic theory

x

ð1Þ

Z2

rmin r ψ=0

E ⊥ = U ( rmin )

ψc

E ⊥ = E.ψ c2

E ⊥ = E.ψ c2 Fig. 1. Deflection of ion by a string of atoms in the continuum approximation.

O. El Bounagui, H. Erramli / Nuclear Instruments and Methods in Physics Research B 268 (2010) 1361–1366

1363

trajectories in low electronic density regions and particles that are steered along channels do not approach the lattice atoms in the axial rows and planes closely enough to undergo a wide angle elastic scattering processes. Consequently stopping power decreases and the particles can go deeper in the matter. Analytical use of channelling techniques needs accurate data for stopping power.

gration step. The scattering angle is assumed to be normally distributed with the mean square deviation proportional to the local energy loss. In this work, a Gaussian distribution was used for sampling the deflection angle [29]. The azimuthal angle is assumed to be uniformly distributed between 0 and 2p. The stopping due to valence electrons is calculated as follows:

3. Calculation methods

Sv al ðrÞ ¼

The use of Monte Carlo (MC) simulations as a tool to analyse and interpret channelling measurement results has now become a common practice. In the simulations of ion channelling, the ion-atom scattering can be treated classically; a screened Coulomb interaction potential is used (based on the Thomas–Fermi model). The numerical approximation of the Thomas–Fermi screening function by Molière [17] is commonly used together with the formula for the screening length presented by Firsov [18] or by Biersack and Ziegler [19]. Ion channelling MC simulation programs have been developed, for example, by Barrett [20], Smulders and Boerma [21] and Dygo and Turos [22]. In these programs, two different approaches are used to calculate the deflections of ions travelling in the crystal; the binary collision model (e.g. Barrett [20]) and Lindhard’s continuum model approximation [12]. In the binary collision (BC) model, ions are presumed to interact only with the closest neighbour atoms. In the continuum potential approximation by Lindhard [12], the charge distribution induced by lattice atoms is presumed to be continuously distributed in the plane perpendicular to the atom strings. The effect of adjoining atom rows on the ion trajectories is therefore easily taken into account, contrary to the BC model. On the other hand, lattice vibrations are more difficult to include in the continuum potential. 3.1. Computer simulations The computer simulation of ion channelling is done by calculating many trajectories of ions impinging on a crystal at randomly chosen positions but in well-defined direction. The calculation and simulation programs have been described in the literature [20–23]. Comparison of measured values with computer simulation data is currently still the most flexible and reliable method for the analysis of large variety of channelling measurements. These programs are devised to perform the simulations and to produce data which can be directly compared with measurements. The simulation method used in this work is based on the calculation of individual trajectories for particles penetrating a single crystal in the direction close to one of the crystal axes. The crystal is treated as a set of static atomic strings, according to the Lindhard model for axial channelling [12]. The string potential chosen is the Molière–Erginsoy potential [23] with the standard set of parameters [24]. Individual particle trajectories are calculated by solving numerically the classical equation of motion in the force field of the surrounding atomic strings and taking into account the inelastic stopping by target electrons. To calculate the energy loss of a channelled particle, we assumed that the contributions to the stopping power from valence and core electrons of the target are independent. The stopping due to valence electrons is calculated using the electron gas model with a realistic distribution of the valence electron density over the channel [24,25]. The contribution of core electrons is obtained using the impact-parameter dependent mean energy loss calculated within the framework of the semiclassical approximation [26]. The Z 31 or Barkas corrections were taken into account in calculating both contributions [27,28]. Multiple scattering of moving particles on target electrons is taken into account after each inte-

  4pZ 21 e4 8p e2 xp Z ½ ð1  a Þ þ aq ðrÞ L ð v ; x Þ 1 þ Z 0 p 1 v al e mv 2 10 v 3 ð4Þ

Here r is the distance to a string, Z1e and v are the projectile charge and velocity, Zval is the number of valence electrons per atom, q(r) is the relative local electron density averaged over the corresponding channel direction, a = 0.5 and xp is the plasma frequency. The stopping number L0(v, xp) is calculated by following the Lindhard–Winther formalism [30]. The last term in Eq. (4) is the Z 31 or Barkas correction in the approximate form suggested by Sung and Ritchie [31]. The above described procedure is applied for calculating the particle trajectory until exit from the crystal if at any moment its distance from any atomic string pffiffiffi is greater than the minimum distance of approach r > rmin = 2u? , where u? is the thermal vibration amplitude in the transverse plane. In these calculations, ions were chosen to be incident on the crystal at angles lower than the critical angle wc, and its distance from any atomic rows is greater than the minimum distance of approach rmin. The incident particle enters in the crystal lattice under an angle w ¼ nwc ; where n is a random number (0 6 n < 1), generated by using a linear congruential random number generator [32]. The values used for the parameters, such as the modulus, the multiplier, the increment and the seed that specify the generator, are identical to those recommended by Lewis et al. [33]. For 1.0 MeV He+ incident on Sih1 1 0i at room temperature, wc = 0.65°. The impact parameter, r, is evaluated after each collision. A particle is assumed to be dechannelled if its distance from an atomic string is smaller than rmin. The electronic density is calculated by using Ziegler data [34], and thermal vibration effects were included in the calculation by averaging each string potential over a distribution function for atomic transverse displacements [35]. 3.2. Calculation of transmitted energy Transmission techniques are used for different kinds of analysis, of which measuring stopping powers is one of the most well known. By measuring the energy loss of particles passing through thin self-supporting foils of which the real density is known exactly, the stopping power can directly be determined [36,37]. The described model was used for the determination of the channelled energy loss of protons and He ions in an energy region from 0.5 to 2 MeV transmitted through a silicon crystal under axial channelling conditions. The calculation was made for a uniformly thick, perfect single crystal. The transmitted energy, Ef, in channelling conditions was calculated through a 3 lm thick Si crystal with incident energy Ei. By exploiting those energies values, i.e. Ei and Ef, one can determine the energy loss for channelled ions. Indeed, we have:

SðEÞchannelled ¼

Ei  Ef x

ð5Þ

where x, is the thickness of the silicon crystal. To see the effect of the channelled stopping on the energy loss, we compared channelled and random energy losses. We have used

1364

O. El Bounagui, H. Erramli / Nuclear Instruments and Methods in Physics Research B 268 (2010) 1361–1366

the electronic energy loss formula for comparison since the main mechanism of the stopping is due to the electronic energy losses at high energies. Thus, we propose the following relation for the conversion between the channelled and random stopping powers:

SðEÞchannelled ¼ eSðEÞrandom

ð6Þ

Here e is a constant ratio and S(E)random is the specific energy loss in a random direction. 4. Results and discussion The simulation program is conceived in the Fortran 90 language under Windows. Random pulling is carried out with the assistance of the Monte Carlo method. The program simulates the trajectories of the incident particles in matter, with the assistance of the random procedure. It enables us to allot to any ion any medium energy (500 keV–5 MeV), a direction, and to follow the evolution of these parameters when the ion crosses the matter by using the Monte Carlo method. We have calculated the ratio, e, of channelled to random stopping powers for protons and 4He ions moving along the h1 1 0i and h1 1 0i directions through the silicon single crystal. The crystal should not be so thick that dechannelling substantially destroys the channelling effect. We have chosen a crystal thickness of 3 lm. The target thicknesses used in our calculations (=3 lm) were chosen so as to cause substantial energy loss and multiple scattering of the particles. We also wanted to have (for w < wc) a sufficient number of particles remaining channelled until the exit from the crystal. The code can be run on any modern PC (512 Kb RAM and 2 GHz Pentium IV processor are recommended) with MS Windows operating system. In this study a total of 106 ion histories were run for each calculation. Table 1 shows the obtained data by calculation for protons and helium ions along the h1 0 0i and h1 1 0i channel axes of silicon. It is found that the stopping power for channelled particles is smaller than for random ones. The channelled stopping power is (in average) 70% of the random stopping power. The statistical error is quoted as a fractional standard deviation, fd. If r is the standard deviation for the ratio of channelled and random energy loss and e is the mean value of this ratio, the fractional standard deviation, fd, is:

r fd ¼ ¼ e

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Pn 2 ffi 1 2 i¼1 ei  ne nðn1Þ

e

The fractional standard deviation is just a number, it does not have any unit and it is often quoted as a percent. Table 2 compares the experimental and theoretical data obtained by other works. The value of the ratio of channelled to random stopping powers was also obtained theoretically by Vos et al. [38], Grande and Schiwietz [39], and experimentally by Appleton et al. [13], Jin and Gibson [40] and Erramli et al. [41]. These results obtained respectively along h1 1 0i axis for 2.3 MeV energy protons

Table 1 Values of the ratio e of channelled to random stopping power. Ions

Channels

Ions Initial energy Ei (MeV)

Ratio e between channelled and random stopping powers

1

h1 0 0i hl 1 0i hl 0 0i hl 1 0i

1 1 2 2

0.66 ± 0.06 0.67 ± 0.02 0.63 ± 0.03 0.69 ± 0.05

H+ H+ 4 He+ 4 He+ 1

in silicon crystal of 8.5 lm thickness, for 4He+ ions channelling along the Sih1 1 0i and h1 1 1i directions, for 9 MeV transmitted protons through 48 lm thick silicon h1 1 0i crystal, h1 0 0i axis for 2 MeV 4He+ ions in a thin Si (100) crystal of 0.88 lm thickness and for transmitted ions trough 3 lm thin silicon crystal. We noticed a reduction in the order of 30% to 40% between the energy losses of the channelled to random along the major axis. Good agreements have been found between data obtained by this method and those gotten experimentally by Jin and Gibson and Erramli et al. and theoretically by Grande and Schiwietz and Vos et al. The differences observed between our calculation and those obtained by Appelton et al. can be explained by the fact that they used thick silicon targets and particles of great energy. Indeed, channelling is only observed on some micrometers and the beam of ions is gradually dechannelled. The interactions that are the origin of the dechannelling cause oscillations of the particles during their movement in the channel. The amplitude of these oscillations increases progressively in target and leads to dechannelling. Thus the larger the crossed thickness is, the more the dechannelling is important and the more e tends towards 1. The results of the ratio, e, of channelled to random stopping power given by our Monte Carlo program and those obtained by MSCPI code [13] were in good agreement with a substantial reduction in computational time and variance. Indeed, our program shows a large Figure of Merit which is 10 times higher than the MSCPI one. The Figure of Merit is defined by: T:1r2 , where T is the computer time. Since T is proportional to the number of trials n, and r2 is proportional to 1/n for large values of n. The Figure of Merit does not depend on n. In the second part of this work, we have studied the angular dependence of energy loss of channelled helium ions along Sih1 0 0i axis. This study was carried out for incident angles randomly distributed in the intervals 0 6 w 6 1:5 : The angular scan for the negative region of the tilt angle w was obtained by symmetry. The crystal thickness used in the current calculation was 3 lm. The calculations results of the angular dependent energy loss of helium channelled along Sih1 0 0i are presented in Fig. 2. In the same Figure, the experimental data of Azevedo et al. [42] at 0.8 and 2 MeV are also shown. We notice a good agreement between our results obtained by computer simulation study and those obtained by using the experimental channelling [40,42]. It can be observed that the ratio between channelled and random stopping powers at w = 0° (well channelled ions) decreases monotonically from 0.87 at 0.5 MeV to 0.63 at 2.0 MeV. This reflects the fact that for these energies above the maximum stopping power (around 0.5 MeV), the angular dependence of the stopping power is due to projectile electron loss and capture process from silicon inner-shell electrons. The reduction of the ratio, e, of channelled to random stopping power is due to the fact: when the 4He ions energies increase, the contribution of inner-shell silicon electrons to random stopping power increase. The contribution of the silicon inner-shell electrons to the total random stopping power becomes larger with increasing projectile energy (10 of 14 electrons are inner-shell electrons). This contribution is almost completely suppressed under channelling conditions and thus, the ratio (channelling to random) decreases significantly at larger energies. Another important feature in Fig. 2 is the width of the angular dependent energy loss. These widths are about half of the channelling dips obtained from the dechannelling profiles. Thus, there are channelling trajectories where 4He ions probe an energy loss rate as large as in a random incident angle. This is a direct consequence of different impact parameter dependencies on the electronic and nuclear energy transfers.

O. El Bounagui, H. Erramli / Nuclear Instruments and Methods in Physics Research B 268 (2010) 1361–1366

1365

Table 2 Comparison of values of e. Ions

Channels

Vos et al. [38]

Grande and Schiwietz [39]

Appleton et al. [14]

Jin and Gibson [40]

Erramli et al. [41]

This work

1

hl 0 0i hl 1 0i hl 0 0i hl 1 0i

– 0.65 – –

– –

– 0.80

– – 0.64 ± 0.02 –

0.67 ± 0.05 – 0.69 ± 0.05 –

0.66 ± 0.06 0.67 ± 0.02 0.63 ± 0.03 0.69 ± 0.05

+

H H+ He+ 4 He+ 1 4

0.62 ± 0.06

Fig. 2. Ratio e of channelled to random stopping power for 0.5, 0.8, 1, 1.5, and 2 MeV He ions transmitted trough a 3 lm Si single crystal along the h1 0 0i direction as a function of the incident angle. The open circles at 0.8 and 2 MeV are the experimental results of Ref. [42].

5. Conclusion In this work, the energy loss of hydrogen and helium transmitted under channelling conditions in silicon single crystal have been calculated by applying the numerical Monte Carlo method. As expected, a decrease has been observed in the channelled energy loss along the h1 0 0i and h1 1 0i silicon axes. One interesting result is that we estimated roughly 30% decrease in the stopping power due to the channelled ions. By using the obtained results, an expression for quick transformation of the random and channelled energy loss is proposed. We also calculated the energy loss of 0.5–2.0 MeV for 4He ions channelled along the silicon h1 0 0i axis as a function of the inci-

dence angle. We noticed that for energies higher than the energy (0.5 MeV) corresponding to the maximum of stopping power, the reduction of the ratio, e, is due primarily to the interactions with the valence electrons of Si. In other hand, if the incident energies increase, the contribution of inner-shell silicon electrons to random stopping power increases. This contribution is almost completely suppressed under channelling conditions and thus, the ratio (channelling to random) decreases significantly at larger energies. The calculation code developed in this work has the advantage to be faster (1 min of CPU time for about 104 ions passing through the target), smaller, more accurate and is a good tool for calculating channelled stopping power of various crystal materials and along all directions.

1366

O. El Bounagui, H. Erramli / Nuclear Instruments and Methods in Physics Research B 268 (2010) 1361–1366

This calculation shows the importance of using the correct stopping power in the experiments combining channelling and other nuclear analytical method. References [1] J.W. Mayer, E. Rimini, Ion Beam Handbook for Material Analysis, Academic Press, New York, 1977. [2] K.B. Winterbon, P. Sigmund, J.B. Sanders, Mat. Fys. Medd. K. Dan. Vidensk. Selsk. 37 (14) (1970). [3] P. Bauer, F. Kastner, A. Arnau, A. Salin, P.D. Fainstein, V.H. Ponce, P.M. Echenique, Phys. Rev. Lett. 69 (1992) 1137–1139. [4] G. Schiwietz, Phys. Rev. A 42 (1990) 296–306. [5] P.L. Grande, G. Schiwietz, Phys. Rev. A 44 (1991) 2984–2992. [6] G. Schiwietz, P.L. Grande, Nucl. Inst. Meth. Phys. Res. B 69 (1992) 10–17. [7] P.L. Grande, G. Schiwietz, Phys. Rev. A 47 (1993) 1119–1122. [8] G.R. Piercy, F. Brown, J.A. Davies, M. McCargo, Phys. Rev. Lett. 10 (1963) 399– 400. [9] R.S. Nelson, M.W. Thomson, Phil. Mag. 8 (94) (1963) 1677–1690. [10] M.T. Robinson, O.S. Oen, Phys. Rev. 132 (1963) 2385–2398. [11] H. Erramli, G. Blondiaux, Appl. Radiat. Isot. 46 (7) (1995) 413–418. [12] J. Lindhard, K. Dan, Vidensk. Selsk. Mat.-Fys. Medd. 34 (14) (1965). [13] M.A. Misdaq, R. Elassali, Nucl. Instr. Meth. Phys. Res. B 119 (1996) 325–330. [14] B.R. Appleton, C. Erginsoy, W.M. Gibson, Phys. Rev. 161 (1967) 330–349. [15] H. Erramli, Thesis, University Cadi Ayyad, Marrakech, 1995. [16] D.S. Gemmell, Rev. Mod. Phys. 46 (1974) 129–227. [17] G. Molière, Z. Naturforschung 2a (1947) 133–145. [18] O.B. Firsov, Soviet Phys. JETP 6 (1957) 534. [19] J.P. Biersack, J.F. Ziegler, Nucl. Instr. Meth. Phys. Res. 194 (1–3) (1982) 93–100. [20] J.H. Barrett, Phys. Rev. B 3 (1971) 1527–1547. [21] P.J.M. Smulders, D.O. Boerma, Nucl. Instr. Meth. Phys. Res. B 29 (1987) 471– 489. [22] A. Dygo, A. Turos, Phys. Lett. A 127 (1988) 281–284. [23] C. Erginsoy, Phys. Rev. Lett. 15 (1965) 360.

[24] L.L. Balashova, N.M. Kabachnik, V.I. Shulga, Ch. Trikalinos, J. Phys.: Condens. Matter 4 (1992) 4883. [25] K. Lenkeit, Ch. Trikalinos, L.L. Balashova, N.M. Kabachnik, V.I. Shulga, Phys. Status Solidi B 161 (1990) 513. [26] N.M. Kabachnik, V.N. Kondratiev, O.V. Chumanova, Phys. Status Solidi B 145 (1988) 103. [27] Ch. Trikalinos, L.L. Balashova, N.M. Kabachnik, V.I. Shulga, Nucl. Instr. Meth. B 69 (1992) 94. [28] L.L. Balashova, N.M. Kabachnik, V.N. Kondratiev, Phys. Status Solidi B 161 (1990) 113. [29] B. Rossi, High-Energy Particles, Prentice-Hall, New York, 1952. Chapter 2: Theory of Electromagnetic Interactions, pp. 10–89. [30] J. Lindhard, A. Winther, K. Dan, Vidensk. Selsk, Mat. Fys. Medd. 34 (4) (1964). [31] C.C. Sung, R.H. Ritchie, Phys. Rev. A 28 (1983) 674. [32] D.E. Knuth, The Art of Computer Programming (TAOCP), Addison–Wesley, 1969. [33] P.A.W. Lewis, A.S. Goodman, J.M. Miller, Pseudo-random number generator for the system/360, IBM Sys. J. 8 (2) (1969) 136–146. [34] J. Ziegler, J.P. Biersack, U. Littmark, The Stopping and Ranges of Ions in Solids, vol. 1, Pergamon Press, New York, 1985. [35] H.F. Krause, S. Datz, P.F. Dittner, J. Gomez de1 Campo, P.D. Miller, C.D. Moak, N. Neskovic, P.L. Pepmiller, Phys. Rev. B 33 (1986) 9. [36] H.H. Andersen, J.F. Bak, H. Knudsen, B.R. Nielsen, Phys. Rev. A 16 (1977) 1929– 1940. [37] A. L’Hoir, D. Schmaus, Nucl. Instr. Meth. Phys. Res. B 4 (1984) 1–12. [38] M. Vos, D.O. Boerma, P.J.M. Smulders, Nucl. Instr. Meth Phys. Res. B 30 (1988) 38–43. [39] P.L. Grande, G. Schiwietz, Nucl. Instr. Meth. Phys. Res. B 136–138 (1998) 125– 131. [40] H.S. Jin, W.M. Gibson, Nucl. Instr. Meth. Phys. Res. B 13 (1986) 76–80. [41] H. Erramli, M.A. Misdaq, G. Blondiaux, Nucl. Instr. Meth. Phys. Res. B 170 (2000) 362–368. [42] G. de, M. Azevedo, M. Behar, J.F. Dias, P.L. Grande, J.H.R. dos Santos, R. Stoll, Chr. Klatt, S. Kalbitzer, Instr. Meth. Phys. Res. B 136–138 (1998) 132–136.