Collisionless current sheet equilibria

Current sheets are important for the structure and dynamics of many plasma systems. In space and astrophysical plasmas they play a crucial role in activity processes, for example by facilitating the release of magnetic energy via processes such as magnetic reconnection. In this contribution we will focus on collisionless plasma systems. A sensible first step in any investigation of physical processes involving current sheets is to find appropriate equilibrium solutions. The theory of collisionless plasma equilibria is well established, but over the past few years there has been a renewed interest in finding equilibrium distribution functions for collisionless current sheets with particular properties, for example for cases where the current density is parallel to the magnetic field (force-free current sheets). This interest is due to a combination of scientific curiosity and potential applications to space and astrophysical plasmas. In this paper we will give an overview of some of the recent developments, discuss their potential applications and address a number of open questions.


Introduction
In many laboratory or natural plasma systems, current sheets play important roles for understanding their structure and their dynamics (see e.g. [1,2]). For example in space and astrophysical plasmas, current sheets are crucial for facilitating energy release and conversion processes such as magnetic reconnection. The focus of this paper will be on collisionless plasmas described by the Vlasov-Maxwell (VM) equations.
For the investigation of many phenomena, for example waves or instabilities, but also for general modelling purposes it is useful to start with configurations that are in equilibrium. VM equilibrium theory has a long history (see the discussion and references in, for example [2][3][4]) and is, of course, not restricted to Cartesian geometry (see e.g. [3] for some references, but also e.g. [5,6]), to spatial variations in just one direction (see e.g. [2,4,[7][8][9][10][11][12][13]) or to systems that are not charge-neutral (e.g. [14]). In particular, collisionless current sheet equilibria with different properties than those considered in this paper can be found by making use of another (adiabatic) constant of motion (see e.g. [10,[15][16][17][18][19]).
Equilibria of collisionless current sheets can be modelled theoretically as configurations which depend on only one Cartesian spatial coordinate, which we will take to be the z-coordinate in this paper. We will only discuss non-relativistic cases (for examples of relativistic treatments see e.g. [4]). We shall also assume that the plasma systems we discuss are quasi-neutral (i.e. the charge density vanishes, but there may be a non-vanishing electric field) or exactly neutral (both the charge density and the electric field vanish). We will focus on the so-called 'inverse' problem in which the microscopic part of a VM equilibrium, i.e. the distribution functions are determined by using information about the macroscopic magnetic field profile. In particular the emphasis of this paper will be on force-free magnetic fields and the 'inverse' problem associated with finding distribution functions for such fields.
The structure of this paper is as follows. In section 2 we will present the basic theoretical background for VM equilibrium theory as used in this paper. In section 3 we will present a number of examples for solutions to the 'inverse' problem for force-free fields, concentrating on recent work and including a completely new example as well. In section 4, a discussion and our conclusions will be presented.

Basic Theory
The non-relativistic VM equations are as follows (e.g. [2]): where s indicates the particle species, f s the distribution function of species s, q s and m s the charge and mass of the species, E the electric field, B the magnetic induction (which we will call the magnetic field from now on), and 0 m and 0  are the permeability of free space and the vacuum permittivity, respectively. The symbols r  and v  indicate gradient operators in configuration space and in velocity space.
In this paper, we are interested in equilibrium solutions of this set of equations (hence t 0 ¶ ¶ = ). We will also assume that our equilibrium solutions are (quasi-)neutral and that the charge density vanishes (this point will be discussed later in a bit more detail). Because we are interested in current sheets we shall assume that the magnetic field B depends only on one Cartesian coordinate, which we will take to be z, and that its z-component vanishes, i.e. that To satisfy equation (5) we can write B in the usual way as the curl of a vector potential A, B A r = ´. Under the assumptions made it is without loss of generality possible to take For later reference we give the form that Ampèreʼs law (2) takes on when expressed in terms of the vector potential: The (stationary) Vlasov equation itself can in general be solved by the method of characteristics, which are identical to the particle orbits in the equilibrium electromagnetic fields. We usually do not know the particle orbits a priori, either because the electromagnetic fields are not known until we solve the equilibrium Maxwell equations or, even if the fields are known, as we will assume later in this paper, because the particle orbits are not explicitly known in full. Therefore, although under the assumptions we have made above the particle orbit problem is in principle integrable this is in practice not of much use and instead one makes use of the constants of motion associated with the symmetries of the problem: time-invariance leading to a constant Hamiltonian Although this is what could be called the standard approach to determining collisionless equilibria it does not cover all possible cases for even one-dimensional equilibria. This problem has already been pointed out by Grad [20] and as we will discuss in a bit more detail later, the method is, for example, unable to represent multiple current sheet cases or certain classes of asymmetric current sheets. Especially in the case of asymmetric current sheets there have been recent developments to calculate VM equilibria of asymmetric collisionless current sheets that are based on taking some aspects of the characteristics into account (see e.g. [21,22]). One can define two different approaches to the VM equilibrium problem. The first one, which we will call the 'forward' approach starts by specifying the equilibrium distribution functions F H p p , , s s xs ys ( ) . One can then calculate the charge and current densities as functions of the electric and vector potentials, x y s s s s xs ys 3 Instead of solving Gauss' law explicitly the assumption of quasi-neutrality is an excellent approximation if the typical length scales of the plasma are much larger than the Debye length. Mathematically this means that one sets the charge density to zero and determines Φ as a function of A x and A y (see e.g. [2,23]; also [3] for a detailed discussion). Using A A , x y F( )in the expressions for j x and j y , Ampèreʼs law in the form of equations (10) and (11) becomes a set of two coupled second order ordinary differential equations, which can be solved by standard methods subject to appropriate boundary conditions, thus determining A z x ( ) and A z y ( ) and by differentiation B z x ( ) and B z y ( ). The second approach, which we will call the 'inverse' approach, starts from a given magnetic field profile B z x ( ) and B z y ( ) and tries to determine from this information compatible equilibrium distribution functions F H p p , , s s xs ys ( ) . Two problems are immediately apparent. As already briefly discussed, there are magnetic field profiles that cannot be generated by distribution functions of the assumed form. We mentioned a couple of examples above, but it would be important to have a method to find out a priori whether or not a magnetic field profile is compatible with the assumed type of distribution function. We will come back to this point a bit later. The second problem is uniqueness, because there are many different distribution functions of the form F H p p , , s s xs ys ( )which result in the same magnetic field profile. As an example, we mention the Harris sheet (B z z L tanh x µ ( ) ( )) which in its original form [24] has distribution functions of the form where T s is the temperature of species s) and u ys a constant which is equal to the bulk flow speed of species s in the y-direction. However, other forms of distribution functions for the Harris sheet magnetic field are known (e.g. [25]).
In this paper we will focus on the 'inverse' approach and we will choose a method which can address both of the problems discussed above. The method we will use has been suggested by Channell [26] (for earlier work using a similar approach, see [27]; see also [23,28]). The method makes use of the fact that the current density can be directly linked to one component of the pressure tensor, which in the coordinate system we use is the P zz component, defined as zz x y s s z s s xs ys Here, we have already taken into account that Channellʼs method imposes not only quasi-neutrality, but exact neutrality ( 0 F = ) and hence P zz will be a function of Ax and A y only. The current density can generally be expressed as (see e.g. [3,23]) leading to Ampèreʼs law taking on the form This set of coupled ordinary differential equations is equivalent to the equation of motion of a particle in a potential, with the role of 'time' played by the coordinate z, the position of the particle given by A z x ( ) and A z y ( ), and the potential in which the particle moves given by P zz 0 m . These properties of the 1D VM equilibrium problem have been noticed before by many authors (see e.g. [3] for a discussion and references).
As an example we show a surface plot for P A A , Because the potential P zz depends only on one of the coordinates A x , A y , the particle will be able to move in the A x direction at constant 'velocity' corresponding to a constant guide field B B y y 0 = . The value of this constant is determined by the initial conditions of the particle motion.
We discuss this analogy to particle motion in a potential in so much detail because it is very useful in addressing the first of the two problems mentioned above, namely to find out a priori what types of magnetic field profiles are inconsistent with the 'inverse' approach as specified so far. As an example we will use an asymmetric current sheet without guide field, i.e. with only B z 0 . Then, as in the Harris sheet case P zz will only be a function of A y . For a particle trajectory to represent a current sheet, i.e. for the magnetic field to reverse direction, the particle trajectory has to have a turning point where it stops and turns around. However, on a surface depending only on one single variable the return branch of the trajectory has to be the same as the inward part of the trajectory, hence ruling out any potential asymmetry. We remark that this changes if a guide field, even if constant, is introduced, because that opens up the possibility of making P zz dependent on A x as well (see e.g. [27,29,30]). Reasoning along the same lines one can, for example, rule out representing configurations such as double current sheets with this approach.
The uniqueness problem can be addressed to some extent by restricting the class of distribution function under consideration. In his method Channell [26] uses equilibrium ) and n s 0 a constant with the dimension of a particle density. Here, g s is an unknown function of the canonical momenta. Before charge neutrality is imposed, P zz has the form For the method to work one has to ensure that [3,31,32]) for all possible values of A x , A y . This imposes additional conditions which the parameters of the distribution function have to satisfy. The neutral P zz is then given by Using the canonical momenta instead of the velocity components as integration variables and using equation (27) ) a known function of A x and A y , this is a Fredholm integral equation of the first type that has to be solved for the unknown function g p p , s xs ys ( ). We remark that different functional dependences of the distribution function on the Hamiltonian H s are possible, leading to integral equations with different kernels. As far as the authors are aware this possibility has not yet been investigated in much detail (except in [33]). One interesting possibility would be to consider an energy dependence leading to Kappa-type distribution functions as in e.g. [25] for the Harris sheet magnetic field.

Examples
In his paper, Channell [26] gives a number of examples of how his method could be used to find distribution functions for given functions N A A , x y ( ). If one would like to start from a given magnetic field profile z B z B z B , ,0 x y = ( ) ( ( ) ( ) ) there is an additional step that one has to carry out, namely to determine N A A , x y ( )and hence P A A , zz x y ( ) from the magnetic field profile. The first step in this process should of course be to determine the vector potential components A z x ( ) and A z y ( ) for the given magnetic field profile. One can determine how P zz varies with z along the given trajectory A z A z , Integrating with respect to z and calling the integration constant P T , one gets Constructing a function of two variables P A A , zz x y ( ) defined over the complete A x -A y -plane from knowledge of P z zz ( ) along one single trajectory A z A z , x y ( ) ( ) is sometimes problematic due to the intrinsic non-unique nature of this task, but as we will see later this can also provide opportunities for finding more convenient solutions. We see from equation (30) that the absolute value of P zz along the trajectory is determined only up to a constant, because the value of P T is arbitrary apart from demanding that P z zz ( ) has to be positive for all z. We remark that equations (20) and (21) are constraints that also have to be satisfied, and imply that along the known trajectory not only P zz is determined, but also its gradient. Apart from these constraints the surface P A A , zz x y ( ) can in principle be deformed arbitrarily (obviously one needs to ensure P 0 zz > ). For the remainder of this section we will focus on a particular type of current sheet magnetic fields, the so-called force-free current sheets. Force-free magnetic fields have a current density that is parallel to the magnetic field (j B  ), i.e. j B 0 m a = .
Here α can be either a constant, leading to linear force-free fields, or a function of position, leading to nonlinear force-free fields. Under the assumptions made we can have at most z a a = ( ). In the context of the 'inverse' approach as discussed above the question of what type of P A A , zz x y ( ) functions are consistent with force-free magnetic fields has been answered in [3]. For force-free fields which together with the force balance equation (29) gives the condition implying that the trajectory A z x ( ), A z y ( ) representing a forcefree solution of Ampèreʼs law has to coincide with a contour of the P zz surface. This is, for example, satisfied by P zz surfaces that are equivalent to attractive central potentials depending on A A x y 2 2 + . Such potentials have circular contours and admit circular trajectories thus satisfying the condition for force-free fields (see figure 2), actually leading to one-dimensional linear force-free fields z z z B sin , cos , 0 A distribution function leading to a nonlinear forcefree field was first found by Harrison and Neukirch [31]. The magnetic field profile for that case was ), with B 0 the constant value of the magnitude of the magnetic field. This magnetic field profile has been named the force-free Harris sheet, because it has the same B z x ( ) as the Harris sheet, but a B z y ( ) which renders it force-free. The form of P zz found in [31] was with L a typical length scale (half-width) of the current sheet configuration and b representing a constant background pressure. The background pressure is needed to keep P zz positive (see figure 3). The form of the distribution functions for the pressure (33) can be found and is Here a s , b s , u xs and u ys are constant parameters of the problem (the other parameters have been defined before). As discussed in [31,32]  = ( ) is the thermal gyroradius of species s. If we assume that all parameters except L and u xs are fixed, one sees that a decrease in the current sheet width L corresponds to an increase in u xs , which will eventually lead to a violation of condition (35) and hence multiple maxima in the v x direction (for a detailed discussion see [32]).
A number of important extensions have been made to this first nonlinear force-free collisionless current sheet equilibrium. The non-uniqueness of the 'inverse' approach was shown by Wilson and Neukirch [33] who showed that the P zz given in equation (33) can be obtained with distribution functions that have a different dependence on the Hamiltonian H s than in equation (34), but the same dependence on the canonical momenta. It was also shown in [40] how the forcefree Harris sheet distribution function (34) can be generalised to the relativistic regime.
An important extension to previous work was made by Abraham-Shrauner [41], who generalised the approach to a whole family of magnetic field profiles of the form x k sn ; ( ) and x k cn ; ( ) are Jacobian elliptic functions, with k 0 1   the modulus. This magnetic field profile includes the previously known cases of the linear force-free field in the limit k 0  and of the force-free Harris sheet in with a is constant parameters. Not surprisingly, the velocity space structure of these distribution functions is of similar complexity to that of the force-free Harris sheet.
Another important generalisation to the case of the forcefree Harris sheet was made by Kolotkov et al [43]. Distribution functions of the type (34) as well as those derived in [33] not only lead to a constant P zz along the force-free solution, but also to a constant particle density and consequently to a constant temperature (if defined by the ratio of P zz and the particle density). However, a constant P zz at the macroscopic level could also be achieved by letting both the particle density and the temperature vary, but keeping their product constant. This is exactly what was achieved in [43] by using modified distribution functions of the form = -+ + - with an additional parameter γ introduced to allow different parts of the distribution function to have different energy dependence (for a similar approach in a different context see e.g. [28,44]). The macroscopic form of P A A , zz x y ( ) given in equation (33) does not change, but some of the relations between microscopic and macroscopic parameters are modified due to the additional dependence of the distribution function on γ.
One of the shortcomings of all nonlinear force-free cases discussed so far is that they all have a plasma beta that is larger than unity ( 1; pl b > for discussion see e.g. [45]) independently of the choice of parameters, as long as the constraint F 0 s > is satisfied. This is unsatisfactory, because force-free fields are usually associated with 1 pl b < . One can, however, make use of the fact that P A A , zz x y ( ) can be changed along the lines described above. For force-free cases a mathematical formulation of this property has been given in [3] showing that a function P A A , where P ff is the constant value of P zz on the force-free contour and x y ( ) is a function that is arbitrary apart from the constraint that the expression on the right hand side of equation (41) has to be positive. As noted in [45,46] a transformation of the form with P 0 a positive constant pressure, leads to (using equation (41)) a transformed value of P zz on the force-free contour of P P ff 0 = , i.e. one can in principle reduce pl b for the force-free solution to any value below one. Starting from equation (33), as a function of Ax and A y , P zz becomes [45] This form of P zz leads to a much more complicated form of equation (28) that needs to be solved, resulting in quite complex mathematical problems. In [45,46] a solution in the form of infinite series of Hermite polynomials is found. The authors also show that the series converges for all values of pl b , in particular 1 pl b < , and under which conditions positivity of the resulting distribution function is to be expected. So far, in all parameter regimes accessible to numerical investigation it has been found that the distribution functions have single maxima in velocity space and are not equal to, but similar to Maxwell distributions. This is a very interesting difference to the distribution function (34). It should be emphasised, however, that only a part of all possible parameter combinations has been explored so far.
To illustrate the transformation method, we will present a different transformation which also leads to force-free solutions with 1 pl b < , but is mathematically less demanding. The transformation we will use is simply Solving equation (28)    = and the macroscopic expression for P zz matches the expression calculated from the v z 2 moment of the distribution function (48).
One has to ensure that F 0 s > for all p xs , p ys and that leads to certain restrictions on the parameter values if one also wants to ensure that P 2 ff < . For example, one has to have b 3 2 < and u v xi i 2 th, 2 has to be approximately smaller than 1.2. Contrary to the work by Allanson and co-workers [45,46], there is no issue with convergence of an infinite series. A preliminary study of the distribution functions for this case for 1 pl b < has only found distribution functions which have a single maximum in velocity space and seem to be quite close to Maxwellian distribution functions (see figure 4 for an example). Further and more detailed investigations will have to be carried out to see whether this is the case for all parameter combinations with 1 pl b < or whether there are cases with 1 pl b < and distributions functions with multiple maxima in velocity space. We would like to emphasise that while the transformation consists simply of squaring P zz , it is not the case that the distribution function is simply squared as well. By trial and error, we have found some examples of multi-peaked distribution functions for transformed cases with 1 pl b > , but the plasma beta had to be greater than approximately 2.5. Due to the larger number of terms in the distribution function for the transformed case it is not at all clear whether an analytical criterion such as equation (35) for the original, untransformed case, can also be found for the transformed case. Also, no immediate conclusions for the distribution function of the transformed case can be drawn from the fact that the distribution function for the untransformed case can have multiple maxima. However, for the case shown in figure 4, the parameter b s in equations (34) and (35) corresponds to a s 6 in equation (48). If we choose the same parameters as used in figure 4 for the transformed case in the original distribution function (34) (setting b a s s 6 = ), we would find that the original distribution function would also only have a single maximum in v x . As stated above, the possible parameter choices for the transformed case are to some extent influenced by the fact that the distribution functions have to be positive and one could speculate that this restriction might play a role in the possible velocity space structures of the distribution functions. It is definitely interesting that so far all equilibrium distribution functions that have been investigated for the cases of the force-free Harris sheet magnetic fields with a 1 pl b < have a relatively simple structure in velocity space, but it is too early to conclude that this is always the case.

Discussion and conclusions
In this paper we have presented a general overview of relatively recent work on the 'inverse' problem for collisionless current sheet equilibria, focussing on one-dimensional forcefree equilibria in Cartesian coordinates. A natural question is whether this can be extended to other geometries, e.g. rotationally symmetric equilibria representing flux tubes and some work in this direction has been undertaken [5,6].
An important application for collisionless current sheet models are plasma and magnetic field structures in planetary magnetospheres and a large amount of work on collisionless equilibria in general has been done in this area (see e.g. [2,8,10,29] for overviews and [5,10,[47][48][49][50][51] for recent work).
Collisionless plasma equilibria are also used as a starting point for numerical investigations of dynamical plasma processes such as magnetic reconnection (e.g. [52]). In many cases the Harris sheet with a constant guide field is used as the underlying equilibrium but force-free magnetic field configurations have recently been used more often as starting point for investigations, for example using particle-in-cell (PIC) simulations. PIC simulations starting with linear force-free fields and the corresponding exact equilibrium distribution functions have been carried out by a number of authors (see e.g. [39,[53][54][55][56]). Simulations have also been done for the force-free Harris sheet, but usually with initial distribution functions that are Maxwell-Boltzmann distributions shifted by the macroscopic bulk flow of the particle species (e.g. [57][58][59][60][61][62][63]. Other authors used the bi-Maxwellian selfconsistent equilibrium distribution function for linear forcefree fields (e.g [34,35,39]) to initialise their simulations (e.g. [64,65]). The equilibrium distribution function for the force-free Harris sheet found in [31] has been used as initial equilibrium for PIC simulations in [66] (see also [56]). It is interesting that the results of most of these simulations with regards to the properties of collisionless reconnection agree despite the differences in the simulation set-up both regarding the initial conditions and the simulation parameters.
On the other hand, the results by Allanson and coworkers [45,46] as well as the (preliminary) results for the quadratic transformation presented in this paper seem to suggest that there are equilibrium distribution functions for the force-free Harris sheet in the regime with 1 pl b < that do not deviate massively from shifted Maxwellian distribution functions. This could explain why a set-up using shifted Maxwellians might be justified, although this needs to be investigated in more detail in the future.
One other interesting question that would be worth addressing in the future concerns the non-uniqueness of the 'inverse' problem. In particular, in cases where several distribution functions are known for the same magnetic field profile, it would be interesting to know which of those distribution functions is in any way 'preferred'. A first step would be a linear stability analysis and one would intuitively assume that, for example, multi-peaked distribution functions would be subject to micro-instabilities. However, it would also be of value to think about criteria which would help to select equilibrium distribution functions according to, for example, a (nonlinear) energy principle.