2D PIC simulations of collisional transport of relativistic electrons in dense plasma

We report results of a series of simulations about electron transport in plasma close to solid density performed with a collisional 2D3V PIC code and compare the results to published ones obtained using hybrid codes. We show that the dispersion of energetic particles remains similar to the one observed in the collisionless case and discuss and compare our results in the light of hybrid codes. The transport of fast electrons in dense cold plasmas remains a challenging problem in many situations of interest. The main tool for the study of kinetic effects in the laser-plasma interaction, i.e. PIC code, is at pain in this regime. First it must be extended in the collisional regime, which has proven to be doable but add a significant burden in term of computing time. Second the stability of the model imposes – a priori – severe constraints on the mesh size and the time step usable. In order to alleviate these constraints several kinds of “reduced models” aiming at the replacement of “brute force” by physical insight have been proposed and used. The basic physical assumption of all these models is that for a given temperature of the background plasma, above some value of the density, collisional damping prevails and suppresses collective modes of the plasma [1]. The use of this assumption has given rise to various flavours of codes. So-called hybrid models ignore the laser/plasma interaction and merely describe the transport of energetic particles in a background plasma assuming quasi-neutrality. The energetic particles are described by PIC particles, while the background plasma is handled by fluid equations. The electric field is obtained by Ohm law applied to the thermal plasma and the Faraday law then gives the magnetic field. One of the difficulties with this kind of model is the initialization of the high energy electrons, whose distribution function is by essence arbitrary. Several studies showed that the results are dependent on the injection used. Attempts have been made to alleviate these constraints, either by coupling this kind of model with a PIC code [2] or by building directly a PIC code presenting an adequate behaviour for dense plasmas. Two of these codes have been recently published. In one of them [3] the collective behaviour is progressively ignored as the density becomes larger than some critical value ndamp , while in the other one [4] Maxwell equations are replaced by MHD equations-as in an hybrid codeabove some critical value ndamp. The study of the published literature shows that the validation of such models is limited to peculiar cases, that a lot of caution has to be exercised in their use and in particular that the hypothesis about collisional damping must remain valid during all the simulation including after heating of the background plasma by energetic particles. Finally the same study of literature shows that similar problems handled with different hybrid codes give different results [2, 5]. In the present paper we report results based on a series of reference PIC simulations, in the sense that all standard conditions for stability have been fulfilled. The use of high order spline interpolation combined with fourth order smoothing has allowed us to increase the mesh up to values large as compared to the Debye length (20 d) but this point in itself has been validated by comparing ae-mail: anne.heron@cpht.polytechnique.fr This is an Open Access article distributed under the terms of the Creative Commons Attribution License 2.0, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited. Article available at http://www.epj-conferences.org or http://dx.doi.org/10.1051/epjconf/20135917015 EPJ Web of Conferences Figure 1. Low frequency magnetic field (unit: 1 corresponds to 104 Tesla). results with smaller value of the mesh size. In all the simulations that have been performed, no self heating is measurable during the length of the simulation. We have used as many particles as possible in order to reduce numerical collisions and improve the statistic. The code uses an original collisional model due to K. Nanbu [6, 7] which is an extension of the model of T. Takizuka and H. Abe [8] usually used. This model was extended to the relativistic case by L. Gremillet [9]. As in [8] it is a binary collision model that mimics a Fokker Planck operator [10]. It uses the same technique for pairing particles as in [8], but computes the cumulative scattering angle of a succession of N small-angle collisions undergone by a particle in plasma (N → ∞). It offers the advantage over the standard model not to be restricted to small deviation imposing small time step for collisions. It naturally gives isotropization if s = νDt > 1 without artificial cut-off. For a given value of the so called isotropization parameter s, it provides the probability function of the angle of deviation. The geometry of the simulation is a plasma slab with a sharp density step (no density gradient) of n/nc = 100 irradiated by a focalized laser presenting a flat top temporal profile and a focal spot of FWHM = 10 m. The results that we present have been obtained for a laser flux of 1019 W/cm2. The size of the system is 72 × 125 m with a vacuum region of approximately 25 m in front of the system and 5 m at the rear of it. To eliminate the role of particle boundary effects the simulations have been stopped as soon as energetic particles reach the boundary. The initial temperature of the plasma is 115 eV. The target is supposed to be fully ionized Al and we choose to set the value of the Coulomb logarithm to 7.7 which corresponds to νei/ p = 1.9. The mesh size is k0dx = 0.03(dx = 20 D) and the time step is 0 dt = 0.015(k0 and 0 are the wave number and the pulsation of the laser). Note that the mesh size is such that p dx/c = 0.3 which means that the skin depth is fully resolved. There are 70 particles per cell for a total of 3.51010 particles. The physical run time is approximately 350 fs. A run necessitates 2048 Blue Gene Nodes or 4 TB of memory and 2106 hours of CPU time. A common diagnostic of the interaction is the low frequency magnetic field. It is expected that the propagation of the beam of energetic particles inside the plasma is associated with a confining toroïdal (in 3D) magnetic field creating a pinching effect of the particles. This has been demonstrated several times with the help of hybrid codes. Figure 1a shows a blow up of magnetic structure at the surface of interaction 25 fs after the impact of the laser. A highly filamented magnetic field is clearly visible, the scale length (and the longitudinal extension) is exactly p/c. We obtained exactly the same results in the collisionless case [11]. Note the intensity of this field that already reaches 5000 T (1 in the units of these pictures is 104 T or 100 MG). Figure 1b shows the same picture but on a larger scale 100 fs after the laser impact. The color scale has been saturated to a value much smaller than the preceding one in order to show the details inside the plasma. One can clearly see at the edge of the beam the reminiscence of what should be a toroïdal field, but it is already filamented. Finally Figure 1c shows

The transport of fast electrons in dense cold plasmas remains a challenging problem in many situations of interest.The main tool for the study of kinetic effects in the laser-plasma interaction, i.e.PIC code, is at pain in this regime.First it must be extended in the collisional regime, which has proven to be doable but add a significant burden in term of computing time.Second the stability of the model imposes -a priori -severe constraints on the mesh size and the time step usable.In order to alleviate these constraints several kinds of "reduced models" aiming at the replacement of "brute force" by physical insight have been proposed and used.The basic physical assumption of all these models is that for a given temperature of the background plasma, above some value of the density, collisional damping prevails and suppresses collective modes of the plasma [1].The use of this assumption has given rise to various flavours of codes.So-called hybrid models ignore the laser/plasma interaction and merely describe the transport of energetic particles in a background plasma assuming quasi-neutrality.The energetic particles are described by PIC particles, while the background plasma is handled by fluid equations.The electric field is obtained by Ohm law applied to the thermal plasma and the Faraday law then gives the magnetic field.One of the difficulties with this kind of model is the initialization of the high energy electrons, whose distribution function is by essence arbitrary.Several studies showed that the results are dependent on the injection used.Attempts have been made to alleviate these constraints, either by coupling this kind of model with a PIC code [2] or by building directly a PIC code presenting an adequate behaviour for dense plasmas.Two of these codes have been recently published.In one of them [3] the collective behaviour is progressively ignored as the density becomes larger than some critical value n damp , while in the other one [4] Maxwell equations are replaced by MHD equations-as in an hybrid code-above some critical value n damp .The study of the published literature shows that the validation of such models is limited to peculiar cases, that a lot of caution has to be exercised in their use and in particular that the hypothesis about collisional damping must remain valid during all the simulation including after heating of the background plasma by energetic particles.Finally the same study of literature shows that similar problems handled with different hybrid codes give different results [2,5].In the present paper we report results based on a series of reference PIC simulations, in the sense that all standard conditions for stability have been fulfilled.The use of high order spline interpolation combined with fourth order smoothing has allowed us to increase the mesh up to values large as compared to the Debye length (20 d ) but this point in itself has been validated by comparing a e-mail: anne.heron@cpht.polytechnique.frThis is an Open Access article distributed under the terms of the Creative Commons Attribution License 2.0, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.results with smaller value of the mesh size.In all the simulations that have been performed, no self heating is measurable during the length of the simulation.We have used as many particles as possible in order to reduce numerical collisions and improve the statistic.The code uses an original collisional model due to K. Nanbu [6,7] which is an extension of the model of T. Takizuka and H. Abe [8] usually used.This model was extended to the relativistic case by L. Gremillet [9].As in [8] it is a binary collision model that mimics a Fokker Planck operator [10].It uses the same technique for pairing particles as in [8], but computes the cumulative scattering angle of a succession of N small-angle collisions undergone by a particle in plasma (N → ∞).It offers the advantage over the standard model not to be restricted to small deviation imposing small time step for collisions.It naturally gives isotropization if s = νDt > 1 without artificial cut-off.For a given value of the so called isotropization parameter s, it provides the probability function of the angle of deviation.
The geometry of the simulation is a plasma slab with a sharp density step (no density gradient) of n/n c = 100 irradiated by a focalized laser presenting a flat top temporal profile and a focal spot of FWHM = 10 m.The results that we present have been obtained for a laser flux of 10 19 W/cm 2 .The size of the system is 72 × 125 m with a vacuum region of approximately 25 m in front of the system and 5 m at the rear of it.To eliminate the role of particle boundary effects the simulations have been stopped as soon as energetic particles reach the boundary.The initial temperature of the plasma is 115 eV.The target is supposed to be fully ionized Al and we choose to set the value of the Coulomb logarithm to 7.7 which corresponds to ν ei / p = 1.9.The mesh size is k 0 dx = 0.03(dx = 20 D ) and the time step is 0 dt = 0.015(k 0 and 0 are the wave number and the pulsation of the laser).Note that the mesh size is such that p dx/c = 0.3 which means that the skin depth is fully resolved.There are 70 particles per cell for a total of 3.510 10 particles.The physical run time is approximately 350 fs.A run necessitates 2048 Blue Gene Nodes or 4 TB of memory and 210 6 hours of CPU time.
A common diagnostic of the interaction is the low frequency magnetic field.It is expected that the propagation of the beam of energetic particles inside the plasma is associated with a confining toroïdal (in 3D) magnetic field creating a pinching effect of the particles.This has been demonstrated several times with the help of hybrid codes.Figure 1a shows a blow up of magnetic structure at the surface of interaction 25 fs after the impact of the laser.A highly filamented magnetic field is clearly visible, the scale length (and the longitudinal extension) is exactly p /c.We obtained exactly the same results in the collisionless case [11].Note the intensity of this field that already reaches 5000 T (1 in the units of these pictures is 10 4 T or 100 MG). Figure 1b shows the same picture but on a larger scale 100 fs after the laser impact.The color scale has been saturated to a value much smaller than the preceding one in order to show the details inside the plasma.One can clearly see at the edge of the beam the reminiscence of what should be a toroïdal field, but it is already filamented.Finally Figure 1c shows IFSA 2011 the magnetic field on an even larger scale 250 fs after the impact.Note that the scale is 2.5 times larger than in Figure1b.The toroïdal structure has completely disappeared leaving place to intense filaments that diverge radially.Close to the surface of interaction some magnetic activities persist with very small wavelength.The large filaments of Figure 1c are absent in collisionless simulation.Further analysis is still necessary to understand their origin.
One of the basic assumptions of hybrid codes is the fact that the current of energetic particles is balanced by a return current due to thermal ones.Therefore we have studied the distribution functions inside the slab.Because we are interested specifically in the thermal particles we have computed the projection f(v x ), f (v y ), f(v z ) rather than their equivalent in p = v.These distributions are local in the sense that they are computed only in the center of the focal spot on a width of few microns.Four slices are presented moving from the right of the plasma at respectively 40 m, 25 m, 15 m, 6 m of the plane of interaction.These distribution functions are plotted as f(v i ) versus v 2 i (i = x, y, z).In this representation a Maxwellian distribution function is a straight line and one gets readily the temperature from the slope.Three cuts are represented in Figure 2a (40, 25, 15 m) 250 fs after the impact of the laser.The narrow triangle in dotted line corresponds to 40 m, the three components all fully superimposed.It corresponds to the limit of penetration of energetic particles and basically there is no plasma response, the slope gives the initial temperature.The next triangle in dashed line corresponds to 25 m, one starts to see the plasma response.Heating is visible but while the two transverse components remain superimposed and therefore indistinguishable, the longitudinal component f(v x ) (thicker line) presents an asymmetry with the slope on the negative side being smaller than on the positive one.This shows the existence of a return current.This trend is even more visible on the last set which is a cut at 15 m from the surface of interaction with a clear asymmetry of the distribution function of the longitudinal component f(v x ).
Figure 2b shows the same distribution functions on a larger scale for the cut performed at 6 m.Again the distribution function f(v y ) (thin line) and f(v z ) (dashed line) are indistinguishable except in the tail (which is less collisional) and are essentially Maxwellian with a temperature of approximately 2.5 keV.The component f(v x ) (thick line) however is different presenting not only a difference between the positive and negative slopes but also a curvature.This shows that in presence of a large electric field collisions are not sufficient to maintain a Maxwellian distribution function and that rather than a return current carried by a drifting Maxwellian one has to think in term of deformation of distribution function.This may be of importance in term of kinetic instabilities.We have also computed the repartition of current density as function of v x in Figure 2c.It corresponds to the integral of v x f(v x ) in the case of Figure 2b.The dotted line corresponds to the value at half maximum.It is clearly visible that half of the return current is carried by particles with energy above 2.5 keV (v x /c ∼ 0.1).Conversely half of the forward current is carried by particles with energy below 10 keV (v x /c ∼ 0.2).The black triangle corresponds to energy of 100 keV which shows that at most 20% of relativistic particles contribute to the forward current.A significant retention of generated energetic particles occurs in the vicinity of the surface of interaction.
Finally Figure 3 shows at the same time, the structure of the forward current density i.e. the current density associated with a positive value of v x (3a) as a function of space, the structure of the total current density (3b), the mean kinetic energy of all the particles with an energy below 10 keV -an arbitrary estimate of the thermal population-(3c) and the structure of the low frequency electric field i.e. filtered with respect to p (3d).The total current density shows that there is locally a slight excess of forward current over the return one.This excess is of the order of 1% which is small in relative value but one has to keep in mind that in these figures the unit of current density is n 0 ec which is 4.810 4 A/ m 2 and that such a current is therefore capable of generating large magnetic fields.We checked that rotB reflects accurately the structure of the total current.
Finally the spatial repartition of the kinetic energy is clearly correlated with the structure of the electric field suggesting heating by Joule effect.The nature of this electric field is definitively of resistive origin; absolutely no charge density is detectable.Note the scale is such that a value of 1 corresponds to field of the order of 3 × 10 12 V/m, which means that these fields are very large and capable of stopping particles of several keV over few microns.

Figure 2 .
Figure 2. Distribution function and repartition of current density.

Figure 3 .
Figure 3. Current density for v x > 0 (a), total current density (b), kinetic energy (c) and longitudinal component of low frequency electric field (d).