Perpendicular Transport of Energetic Particles in Magnetic Turbulence

Scientists have explored how energetic particles such as solar energetic particles and cosmic rays move through a magnetized plasma such as the interplanetary and interstellar medium since more than five decades. From a theoretical point of view, this topic is difficult because the particles experience complicated interactions with turbulent magnetic fields. Besides turbulent fields, there are also large scale or mean magnetic fields breaking the symmetry in such systems and one has to distinguish between transport of particles parallel and perpendicular with respect to such mean fields. In standard descriptions of transport phenomena, one often assumes that the transport in both directions is normal diffusive but non-diffusive transport was found in more recent work. This is in particular true for early and intermediate times where the diffusive regime is not yet reached. In recent years researchers employed advanced numerical tools in order to simulate the motion of those particles through the aforementioned systems. Nevertheless, the analytical description of the problem discussed here is of utmost importance since analytical forms of particle transport parameters need to be known in several applications such as solar modulation studies or investigations of shock acceleration. The latter process is directly linked to the question of what the sources of high energy cosmic rays are, a problem which is considered to be one of the most important problems of the sciences of the 21st century. The present review article discusses analytical theories developed for describing particle transport across a large scale magnetic field as well as field line random walk. A heuristic approach explaining the basic physics of perpendicular transport is also presented. Simple analytical forms for the perpendicular diffusion coefficient are proposed which can easily be incorporated in numerical codes for solar modulation or shock acceleration studies. Test-particle simulations are also discussed together with a comparison with analytical results. Several applications such as cosmic ray propagation and diffusive shock acceleration are also part of this review.


Introduction
A fundamental problem in the sciences of the 20th and 21st centuries is to understand the physics of cosmic rays. This includes both, the acceleration mechanism leading to the creation of such particles and their propagation through space. The exploration of cosmic rays goes back to famous experimentalists such as Theodor Wulf and Victor Hess who measured the rate of ion production via an electrometer on top of the Eiffel Tower and via balloon flights. Such early measurements lead to the assumption that there is radiation coming from space rather than originating in Earth's crust. Victor Hess received the Nobel Prize of physics in the year 1936 for the discovery of cosmic rays.
Cosmic rays consist mostly of protons, but also of alpha particles, electrons, and heavy ions and their energy ranges from some keV up to about 10 21 eV. This radiation can indeed be dangerous for humanity. For air crews and astronauts, for instance, the radiation can cause damage to their tissue even leading to cancer. Furthermore, cosmic particles can damage electronic devices leading to their malfunction and even some crashes of satellites are believed to be related to the influence of these particles. The danger of cosmic radiation would be a significant problem during a manned mission to Mars.
After the discovery of cosmic radiation by Wulf and Hess, the question arose quickly where these particles are coming from. It is believed that the mechanism leading to the creation of cosmic particles is diffusive shock acceleration sometimes called centrifugal mechanism of acceleration or Fermi acceleration. Candidates for particle accelerators are coronal mass ejection driven shocks, supernova shocks, and active galactic nuclei. In such acceleration sites one finds propagating shock waves similar in comparison to the shock wave associated with an atomic explosion. In astrophysics, however, the shock waves have moving magnetic inhomogeneities. If an electrically charged particle crosses the shock front and encounters a moving change in the magnetic field, it can get reflected back at increased velocity. The same process can occur while the reflected particle crosses the front again. Multiple reflections of this type will increase the particle's energy significantly.
After leaving the acceleration site, cosmic particles propagate through space. However, the space between planets, stars, and even galaxies is filled with a magnetized plasma. Therefore, one finds turbulent magnetic fields which influence the motion of the particles. If those fields would be absent, the electrically charged particles would perform a perfect helical motion due to the interaction with large scale or mean magnetic fields. The turbulent fields lead to diffusion along and across that background field making it very difficult to predict the motion of particles through space. Therefore, the motion of cosmic particles as well as their acceleration are described by complicated transport equations (see, e.g., Schlickeiser 2002;Zank 2014).
In magnetized plasmas one finds chaotic electric and magnetic fields associated with the turbulent motion of the fluid. Due to the turbulent behavior of magnetic fields, for instance, magnetic field lines are stochastic curves. This means that if a field line goes through a particular point, one can only estimate the probability to find that field line somewhere else. To study the stochastic behavior of magnetic field lines in turbulence is done in the theory of field line random walk (FLRW). Energetic particles such as protons, electrons, and heavy ions are electrically charged. Therefore, they interact with turbulent magnetic fields as described by the Newton-Lorentz equation. Originally the interaction of energetic particles and magnetic turbulence was described in the pioneering work of Jokipii (1966). It was shown there that one needs to distinguish between diffusion of particles along and across the mean magnetic field, as both processes behave very differently. Parallel diffusion, for instance, is caused by pitch-angle scattering over an extended period of time. Perpendicular diffusion, on the other hand, is often assumed to be associated with the random walk of chaotic magnetic field lines. However, this random walk alone cannot explain all aspects of cross field diffusion. If particles follow random walking magnetic field lines, this would either lead to an energy independent perpendicular mean free path (if the parallel motion is assumed to be unperturbed) or to sub-diffusive transport (if the parallel motion is assumed to be diffusive). It was shown via test-particle simulations that such results are unrealistic since often perpendicular transport is diffusive rather than sub-diffusive. Furthermore, in general the perpendicular mean free path shows a non-trivial energy-dependence.
The understanding of the motion of energetic particles through a plasma itself is an interesting and fundamental problem of plasma physics, space science, as well as astronomy and astrophysics. However, there are several applications of results obtained from studies of particle transport. Analytical forms of the different diffusion coefficients are frequently used to explore the following problems: 1. In the theoretical investigation of the acceleration of particles at shock waves one usually solves a diffusive transport equation in order to compute the cosmic ray spectrum. Such transport equations contain diffusion coefficients in the different directions of space. The corresponding equation is then either solved analytically or numerically (see, e.g., Zank et al. 2004;Dosch and Shalchi 2010;Li et al. 2012;Ferrand et al. 2014;Hu et al. 2017). 2. Energetic particles such as solar energetic particles or cosmic rays propagate through the solar system, the Milky Way, or external galaxies. Like shock acceleration scenarios their motion is described via diffusive transport equations containing again the different diffusion parameters. In order to understand how the charged particles propagate in the different systems, one needs to understand their diffusion process (see, e.g., Bieber et al. 1994;Shalchi and Schlickeiser 2005;Buffie et al. 2013). 3. In solar modulation studies one explores how the cosmic ray flux is altered by solar activity. This is related because the flux of energetic particles strongly depends on the solar wind. The interaction is either described by the pitch-angle dependent Fokker-Planck equation or the aforementioned diffusive transport equation. Therefore, analytical forms of diffusion coefficients are also highly relevant in investigations of solar modulation (see, e.g. Hitge and Burger 2010;Burger 2010, 2014;Wawrzynczak and Alania 2010;Alania et al. 2011;Potgieter and Nndanganeni 2013;Potgieter 2013;Manuel et al. 2014;Ahluwalia and Ygbuhay 2015;Shen and Qin 2018).
In general, the equation of particle motion in the different scenarios described above can be very complicated. Often the transport can be well described by Parker's famous transport equation (see Parker 1965) (1.1) Here we have used the momentum of the particle p, the plasma bulk velocity u, and the components of the spatial diffusion tensor κ ij . On the left hand side we find variation in time and convection whereas on the right hand side we have the terms of spatial diffusion, energy gains, as well as a source term. A more general version of this equation with momentum diffusion as well as adiabatic deceleration can be found in Schlickeiser (2002, Eq. (12.3.23)).
In the rest frame of the moving plasma we have u = 0 and Eq. (1.1) turns into a usual diffusion equation. Assuming axi-symmetry with respect to the mean magnetic field allows us to use a diffusion tensor of the following form where we have used the parallel diffusion coefficient κ , the perpendicular diffusion coefficient κ ⊥ , and the drift coefficient κ A . In the current review article we focus on the perpendicular diffusion coefficient and treat the parallel diffusion coefficient as a variable entering non-linear theories for perpendicular transport. The drift coefficient is not further discussed in this review but it was explored in Engelbrecht et al. (2017) and its importance in solar modulation studies was emphasized in Engelbrecht and Burger (2015). Problematic in analytical studies of particle transport is that turbulent magnetic fields have to be evaluated along the true particle trajectory. Due to the chaotic nature of the particle motion, however, this is not possible without employing approximations. The original work of Jokipii (1966) employed a perturbation theory approach in which the magnetic fields are evaluated along the unperturbed particle orbit corresponding to a perfect helical motion. This approach is usually called quasi-linear theory (QLT). From a more mathematical point of view this approach is questionable because true particle orbits decorrelate from unperturbed orbits if time passes and strictly speaking the quasi-linear approximation cannot be justified. However, quasi-linear theory was somewhat successful with reproducing observed parallel diffusion coefficients in interplanetary space as shown in the famous work of Bieber et al. (1994). However, even in the context of parallel transport, there are problems associated with this theory. For instance, quasi-linear theory cannot explain pitch-angle scattering at pitch-angles close to 90 • (see, e.g., Shalchi 2009a for a review). Furthermore, there seems to be a non-linear contribution to pitch-angle scattering related to perpendicular scattering and the transverse structure of the turbulence (see again Shalchi 2009a).
While quasi-linear theory is somewhat successful with describing transport of particles along the mean magnetic field, it fails for perpendicular transport in almost all cases. 2 Over years several alternative theories were developed such as non-linear closure approximation of Owens (1974) and the Bieber and Matthaeus (1997) model. Over time the performance of computers was improved drastically leading to an alternative tool to compute particle orbits and transport parameters of energetic particles interacting with turbulence (see, e.g., Giacalone andJokipii 1994, 1999;Michałek and Ostrowski 1996;Mace et al. 2000, and Qin et al. 2002a, 2002b. In such simulations one is still required to specify a certain turbulence configuration as it is necessary in analytical theories but apart from this, the obtained results are exact. However, the simulations do not provide analytical forms for particle diffusion coefficients as needed in the various applications listed above. As test-particle simulations became available, they were used to test the validity and accuracy of analytical theories derived previously. It was shown that in particular for perpendicular transport previous approaches fail and in some cases there was complete disagreement between simulations and analytical results. A major step forward in the theory of perpendicular transport was the development of the so-called non-linear guiding center (NLGC) theory of . The latter theory is based on a series of assumptions and approximations leading to a non-linear integral equation for the perpendicular diffusion coefficient. The theory was the first analytical approach which showed good agreement with simulations. However, there were at least two problems associated with this theory: 1. As NLGC theory was developed, it was already known from simulations that perpendicular transport is sub-diffusive if the turbulence lacks transverse complexity. However, NLGC theory provides a finite perpendicular diffusion coefficient for this type of turbulence which contradicts the sub-diffusive result. Later it was shown that NLGC theory also fails for three-dimensional turbulence and small and intermediate Kubo numbers. 2. The agreement between NLGC theory and simulations was achieved only by incorporating a factor a 2 = 1/3. The physics and the reason for this value remained unclear.
Therefore, it was concluded that NLGC theory was a major step forward in the analytical description of perpendicular transport, but this theory cannot be the final solution of the problem.
In order to further improve our understanding of perpendicular transport, several alternative theory have been developed in more recent years. The unified non-linear transport (UNLT) theory provided an integral equation similarly compared to NLGC theory (see Shalchi 2010). However, UNLT theory has several strengths compared to previous approaches. Most important is that UNLT theory provides a vanishing perpendicular diffusion coefficient for turbulence without transverse structure as seen in test-particle simulations. Furthermore, the theory contains the Matthaeus et al. (1995) theory of random walking magnetic fields lines. Therefore, UNLT theory can be understood as unified theory for magnetic field lines and perpendicular particle transport. Furthermore, the theory works for threedimensional turbulence in the small Kubo number regime (see . Then, on the other hand, the theory still requires the correction parameter a 2 = 1/3 if turbulence with large Kubo numbers is considered. A further improvement of our understanding was achieved via the time-dependent version of UNLT theory developed in Shalchi (2017) and Lasuik and Shalchi (2017). This theory is mathematically more complicated because it provides an integro-differential equation for perpendicular transport. This improved version of UNLT theory describes timedependent phenomena such as ballistic and sub-diffusive transport, but it was also the first theory which explained the recovery of diffusion and that this effect is entirely due to the transverse complexity of turbulence. Unfortunately, time-dependent UNLT theory still requires the parameter a 2 in some cases. Therefore, further improvements of non-linear analytical theories of perpendicular transport are required.
Besides the development of systematic analytical theories such as NLGC and UNLT theories, there were also attempts to describe perpendicular transport based on heuristic arguments. The most famous work in this context is the article published by Rechester and Rosenbluth (1978) with the focus on energetic particle transport in fusion plasmas. In the latter paper the importance of magnetic fields lines as well as the suppression of perpendicular transport to a sub-diffusive level due to parallel diffusion was discussed. However, due to Coulomb collisions and exponential field line separation, perpendicular transport can be restored. However, in space and astrophysical plasmas, Coulomb collision should be irrelevant and magnetic field lines should not separate exponentially. Therefore, it was often stated that Rechester and Rosenbluth (1978) does not apply in astrophysical scenarios (see, e.g.,  for a brief discussion of this matter). However, in Shalchi (2015a) a Rechester & Rosenbluth type of diffusion coefficient was derived from UNLT theory indicating that there is at least some connection between Rechester and Rosenbluth (1978) and energetic particle transport in space plasmas. In Shalchi (2019a) a heuristic approach for perpendicular transport was developed providing different formulas for the perpendicular diffusion coefficient. This finally provided an explanation of the parameter a 2 used in previous analytical theories. Furthermore, this heuristic approach explained how systematic theories could be improved in the future which could lead to a complete understanding of perpendicular transport.
It is the purpose of this review article to discuss developments in the analytical theory of perpendicular diffusion over the past 50 years. This also includes a brief review of heuristic approaches and test-particle simulations. It will be shown that perpendicular diffusion depends on the properties of the turbulent magnetic fields but also on parallel diffusion. Therefore, this review will start with a discussion of various turbulence models which were proposed in the literature over the past view decades (Sect. 2) followed by a review of theories developed for field line random walk (Sect. 3) a process that often controls perpendicular transport. Thereafter, the reader will find a short discussion of parallel particle transport (Sect. 4). However, parallel diffusion itself is complicated and still subject of current research. The main focus of this review is on perpendicular diffusion of energetic particles (Sect. 5) with the emphasis on the unified non-linear transport theory including a discussion of different transport regimes (Sect. 6), time-dependent transport (Sect. 7), simple analytical forms (Sect. 8), and a recently developed heuristic approach (Sect. 9). Thereafter, there is a discussion of numerical tools used in transport theory as well as a comparison between simulations and analytical theory (Sect. 10). Although not the central point of this review, the reader can also find some applications of the results discussed in this review (Sect. 11) such as particle propagation through interplanetary and interstellar spaces as well as the role of perpendicular diffusion in the theory of diffusive shock acceleration. At the end of this article there will be a summary, a conclusion, and a short outlook (Sect. 12) discussing unsolved problems and possible future projects.

Analytic Models for Magnetic Turbulence
In analytical theories for perpendicular diffusion the components of the so-called spectral tensor are required as input as shown in Sect. 5. In the following we discuss different models which were proposed in the past. It needs to be emphasized that the theoretical study of turbulence is an ongoing field of research. Therefore, the models discussed in the following are not supposed to be the final solution to the problem. Instead they should be understood as examples sometimes motivated by solar wind observations or theoretical work. After presenting these models, fundamental turbulence scales such as integral scales and the ultrascale are discussed.

Correlation and Spectral Tensors
Especially in astrophysics and space science we deal with magnetic turbulence. The knowledge of the properties of these stochastic magnetic fields is important in several applications such as the theory of field line random walk and cosmic ray propagation. We consider a physical system where the total magnetic field is a position of a mean field B 0 and a turbulent component δB n B(x, t) = B 0 e z + δB (x, t). (2.1) Whereas the mean field is assumed to be constant, the turbulent field is a stochastic quantity depending on space and time. We can describe magnetic turbulence via the magnetic correlation tensor in the configuration space whose components are defined as R nm (x, x 0 , t, t 0 ) = δB n (x, t)δB * m (x 0 , t 0 ) . (2. 2) The brackets used here denote the ensemble average and δB * m corresponds to the complex conjugate of the turbulent magnetic field component δB m . The functions R nm describe how the magnetic field at position x and time t is related to the magnetic field at position x 0 and time t 0 . Therefore, the functions R nm are also called the two-point-two-time correlations. Alternatively, we can describe turbulence in the Fourier space. The magnetic fields in configuration and wave vector space are linked to each other via a usual Fourier transformation δB n (x, t) = d 3 kδB n (k, t)e ix·k . (2.3) Therewith, the two-point-two-time correlations defined via Eq. (2.2) can be written as It has to be noted that here x and x 0 are well-defined positions in space. An example for a well-defined orbit would be a space probe measuring magnetic fields via a magnetometer. In this case the trajectory of the probe is uncorrelated to the magnetic field. The situation is different in the theory of field line random walk and energetic particle transport where the position vectors are stochastic quantities somehow related to the magnetic fields. This will require special techniques as discussed later in this review. If the turbulence is homogeneous, the correlations discussed here must only depend on the distance between the two considered points. Furthermore, if the turbulence is stationary, we can apply the same argument on time meaning that R nm (x, x 0 , t, t 0 ) = R nm (x − x 0 , t − t 0 ).
(2.5) Therefore, we require δB n (k, t)δB * m k , t 0 = P nm (k, t)δ k − k (2.6) where we have used the Dirac delta and the time-difference t = t − t 0 . Eq. (2.6) defines the components of the spectral tensor P nm . This tensor is a very fundamental quantity in the theory of turbulence and energetic particle transport as we shall see later. By employing Eq. (2.6), we derive from Eq. (2.4) R nm (x, t) = d 3 kP nm (k, t)e ik·(x−x 0 ) (2.7) clearly showing that now the turbulence is indeed homogeneous and stationary. Usually we set x 0 = 0 and t 0 = 0 for convenience so that we find R nm (x, t) = d 3 kP nm (k, t)e ix·k . (2.8) Here we still allow the spectral tensor to be time-dependent. Table 1 The dynamical correlation function (k, t) controls the time-dependence of magnetic turbulence (see, e.g., Eq. (2.9) of the current review). This table shows the different models commonly used in the literature. The damping model of dynamical turbulence as well as the random sweeping model were discussed in Bieber et al. (1994). Various plasma wave propagation models including damping, on the other hand, were discussed in Schlickeiser (2002). In this table we have used the Alfvén speed v A , the plasma wave frequency ω, the damping parameter γ , as well as the free parameter α

Model Dynamical correlation function
Magnetostatic turbulence (k, t) = 1 Undamped shear Alfvén waves (k, t) = e ±iv A k t General wave propagation model with damping (k, t) = e iωt−γ t Damping model of dynamical turbulence (k, t) = e −αv A kt Random sweeping model (k, t) = e −(αv A kt) 2 Often we assume the same temporal behavior of all tensor components meaning that we write (see, e.g., Edwards 1964;Bieber et al. 1994 for more details) P nm (k, t) = P nm (k) (k, t) (2.9) where we used the static tensor components P nm (k) and the so-called dynamical correlation function (k, t). The latter function contains all information about the turbulence dynamics. This includes, for instance, wave propagation effects. Over the past three decades scientists developed a more complete understanding of the time scales of turbulence (see, e.g., Matthaeus et al. 1990;Tu and Marsch 1993;Zhou et al. 2004;Oughton et al. 2006). Based on this improved understanding, Shalchi et al. (2006) have formulated the so-called non-linear anisotropic dynamical turbulence (NADT) model for the dynamical correlation function. Some examples of models for dynamical turbulence are listed in Table 1. In some articles perpendicular diffusion was explored for dynamical turbulence (see, e.g., Shalchi et al. 2004bShalchi et al. , 2006Shalchi 2014;Hussein and Shalchi 2016) but it seems that dynamical turbulence effects are less important for perpendicular transport. 3 Therefore, we focus on static turbulence where, by definition (k, t) = 1, throughout this review. However, it has to be emphasized that dynamical turbulence can be very important for parallel diffusion (see, e.g., Bieber et al. 1994;Shalchi et al. 2006). As shown later in this review, the perpendicular diffusion coefficient depends upon the parallel diffusion coefficient. Therefore, dynamical turbulence effects do indeed have an influence on perpendicular diffusion but this influence occurs only indirectly via the parallel diffusion coefficient. 4 In most analytical considerations presented in this review, the parallel diffusion coefficient is only an input parameter in theories for perpendicular transport and thus, we do not have to worry too much about dynamical turbulence effects. In the context of parallel transport the magnetostatic approximation can often be justified by considering particles moving much faster than the Alfvén speed. An example are galactic cosmic rays propagating through the interplanetary space.
In studies of particle transport one also often assumes that the turbulence is axisymmetric with respect to the mean magnetic field B 0 = B 0 e z . However, there are several studies who dealt with cases where turbulence is no longer axi-symmetric. Ruffolo et al. (2008), for instance, studied perpendicular diffusion of energetic particles in two-component asymmetric turbulence. Furthermore, we often neglect magnetic helicity but this effect and its importance in particle transport theory was also explored in some previous work (see, e.g., Dung and Schlickeiser 1990;Tautz and Lerche 2011). As shown in Matthaeus and Smith (1981), the components of the spectral tensor have the following form P nm (k) = g(k , k ⊥ ) δ nm − k n k m k 2 ⊥ (2.10) for axi-symmetry and vanishing magnetic helicity. In Eq. (2.10) we have used the Kronecker delta δ nm and n, m = x, y. Furthermore, we set P nm = 0 if n = z and/or m = z meaning that we consider the reduced case δB z = 0 corresponding to incompressible turbulence. The theories for field line random walk and perpendicular particle transport discussed in this review were derived for this case but it is usually assumed that those theories are valid as long as δB z < B 0 . Since we usually deal with axi-symmetric cases, we employ cylindrical coordinates (k ⊥ , , k ) rather than Cartesian coordinates (k x , k y , k z ) in the spectral function g(k , k ⊥ ). The relations between these two sets of coordinates are The form given by Eq. (2.10) automatically satisfies the solenoidal constraint ∇ · δB = 0 corresponding to k · δB = 0. In order to show this we consider n k n P nm = g(k , k ⊥ ) n k n δ nm − n k 2 n k m k 2 ⊥ = 0. (2.12) In the next few subsections, we discuss different models for the function g(k , k ⊥ ). This will include models with reduced dimensionality, superpositions thereof, but also full threedimensional turbulence models.

Slab Turbulence
One of the simplest models for magnetic turbulence is the slab model. Physically this model is motivated in terms of shear Alfvén waves propagating along the mean field and it was often employed in early studies of particle diffusion (see, e.g., Jokipii 1966). In slab turbulence the stochastic magnetic field depends only on the coordinate along the mean field, i.e., δB(x) = δB(z). Therefore, the components of the magnetic correlation tensor are given by for n, m = x, y. Clearly this is a special case of the form given by Eq. (2.10). Due to the Dirac delta δ(k ⊥ ) in Eq. (2.13), all wave vectors are aligned parallel with respect to the mean magnetic field. Figure 1 shows flux surfaces for slab turbulence. In this case there is no transverse structure of the turbulence. Transverse complexity, however, is an important feature in the theory of perpendicular diffusion as we shall see later.

Reprinted with permission from
The American Astronomical Society-  Different models have been discussed in the literature for the model spectrum g slab (k ). Figure 2 shows a measured spectrum as obtained in the heliosphere via magnetometer measurements performed by the Helios 2 space probe. Figure 2 shows the spectrum as a function of frequency. A comprehensive discussion of the relation between frequency and wave number is given in . Assuming that one has pure slab turbulence, for instance, this relation is v sw k = 2πf (see equation (11) of  for the case = 0) where v sw is the solar wind speed. In general one expects to find three different ranges of spectral scales as explained nicely in Zhou and Matthaeus (2005): 1. The energy range: The energy containing scales control the overall dynamics of turbulence and are responsible for energy transport. Those are the largest scales corresponding to small wave numbers or frequencies. Figure 2 indicates that in the energy range the spectrum scales like ∝ k −1 . 2. The inertial range: It is assumed that the dynamics of the scales of the inertial range are not influenced by the low-frequency scales of external energy sources nor are they influenced by the high-frequency dissipation and magnetic diffusivity scales. This is the fundamental reason for the universal nature of the inertial range. According to Fig. 2 we find that the spectrum scales like ∝ k −5/3 as predicted in the famous work of Kolmogorov (1941). 3. The dissipation range: The energy of the turbulence is transferred through the inertial range to the small scales where dissipation takes place. The turnover from the inertial range to the dissipation range is characterized by a sudden change of the spectral index. In Fig. 2, for instance, we can clearly see a change from a Kolmogorov behavior to a steeper spectrum with ∝ k −2.85 .
As also discussed in Zhou and Matthaeus (2005), there could also be what they call a fardissipation range where the spectrum should decrease exponentially with wave number. Matthaeus et al. (2007), on the other hand, discussed in more detail the largest scales of magnetic turbulence and discussed the possibility of having different sub-regimes in the energy range. Important for perpendicular diffusion are the large scales of the energy range (see, e.g., Qin and Shalchi 2012). Therefore, this part of the spectrum cannot be neglected. Furthermore, the intermediate scales of the inertial range are at least somewhat Fig. 2 The measured turbulence spectrum in the heliosphere as obtained by the Helios 2 mission. As discussed in , frequency and wave number k are directly proportional to each other. Clearly, we can see the three ranges of the spectrum: for small wave numbers we find ∝ k −1 (energy range), for intermediate wave numbers we find a Kolmogorov spectrum ∝ k −1.7 (inertial range), and for large wave numbers, we find a steep spectrum ∝ k −2.85 (dissipation range). Reprinted with permission from NASA Conference Publication- Denskat and Neubauer (1982) important. They are not directly relevant for perpendicular transport but essential for parallel diffusion. As shown later in this review perpendicular diffusion depends on parallel transport meaning that there in an indirect influence of the inertial range. A slab spectrum which takes into account energy and inertial ranges has been presented in Bieber et al. (1994) where the following model has been used (2.14) In this form it is assumed that the spectrum is perfectly flat for the small wave numbers of the energy range. The characteristic scale at which the turnover from the energy range to the inertial range occurs is called the bendover scale. For larger wave numbers the Bieber et al. spectrum decreases with ∝ k −s where the inertial range spectral index s is often set to s = 5/3 as originally obtained in Kolmogorov (1941). It has to be pointed out that the original work of Kolmogorov dealt with a purely hydrodynamic system and the considered turbulent field is the velocity field. Here we consider magnetic turbulence which can be different. However, as shown in Fig. 2, there are indications that one can indeed find a Kolmogorov spectrum in real astrophysical scenarios. This is at least true in the context of solar wind turbulence. It is convenient to introduce the notation where δB n is the total turbulent magnetic field strength in the n-direction. Any type of turbulence model described by the corresponding spectral tensor has to be normalized so that where we used the total turbulent magnetic field δB. The quantity δB 2 corresponds to the total magnetic energy density except some numerical factors. 5 Often in this review, the slab model is just an example, we consider the incompressible case where δB z = 0 and, thus, P zz = 0. If the slab model (2.13) and spectrum (2.14) are used in condition (2.16), it follows that 6 where we solved the k -integral with the help of Gradshteyn and Ryzhik (2000). Therefore, we find for the normalization function where we have used gamma functions. For a Kolmogorov spectrum, for instance, we have C(s = 5/3) ≈ 0.1188.

Two-Dimensional Turbulence
Another popular model for magnetic turbulence is the two-dimensional (2D) model. This model can be understood as the opposite compared to the slab model discussed in the previous subsection. By definition we assume here that δB(x) = δB(x, y). A detailed motivation 5 The energy density of a magnetic field is given by w = B 2 /(2μ 0 ) in SI units and w = B 2 /(8π) in cgs or Gaussian units. 6 In this review we do not delve too much into the theory of delta distributions but it should be noted for clarity that in cylindrical coordinates we have ∞ 0 dk ⊥ δ(k ⊥ ) = 1. If the Dirac delta involving the parallel wave number is integrated, on the other hand, we have +∞ −∞ dk δ(k ) = 1. Both relations are frequently used in this review.
for two-dimensional modes is given by . Although the focus in the latter paper as well as in this review is on space and astrophysical plasmas, it should be noted that two-dimensional modes are also well-known in studies of tokamak plasmas (see, e.g., Zweben et al. 1979).
In the case of two-dimensional turbulence, the components of the spectral tensor are if n, m = x, y and P nz = P zm = P zz = 0. In this particular model the magnetic field vector as well as the spatial dependence are two-dimensional meaning that they are contained in the plane perpendicular with respect to the mean field. Equation (2.19) is a special case of Eq. (2.10) meaning that the solenoidal constraint is automatically satisfied. Again we need to specify the model spectrum g 2D (k ⊥ ). As for the slab model, we employ a form which contains energy and inertial ranges. The following form was proposed in Shalchi and Weinhorst (2009) (2.20) The latter spectrum contains the characteristic scale ⊥ denoting again the turnover from the energy to the inertial range. In the inertial range the spectrum scales like ∝ k −s ⊥ whereas in the energy range it scales like ∝ k q ⊥ . In order to determine D(s, q), we employ the normalization condition (2.16). Using Eqs. (2.19) and (2.20) in that condition allows us to derive (see again Gradshteyn and Ryzhik 2000 for the solution of the occurring integral) where we assumed that q > −1 and s > 1. Otherwise the k ⊥ -integral would not be convergent. This means that spectrum (2.20) is only correctly normalized if these restrictions for the two spectral indices are satisfied. Furthermore, it follows that the normalization function in the spectrum has the following form  Fig. 3 for different values of the energy range spectral index q. For most calculations involving two-dimensional modes, we set s = 5/3 and q = 3. The possible values of the energy range spectral index were discussed in the work of Matthaeus et al. (2007). The latter authors computed different scales of turbulence such as integral scales and the ultra-scale for different spectra. Requiring that those scales are finite, the work of Matthaeus et al. (2007) leads to the assumption that q > 1 (see Sect. 2.10 for the calculation of integral and ultra-scales) in the model spectrum given by Eq. (2.20). It has to be mentioned that for q = 0, the two-dimensional spectrum of Eq. (2.20) becomes identical to the slab spectrum (2.14). This type of spectrum was very popular a few decades ago (see, e.g., Bieber et al. 1994;Shalchi et al. 2004a) but motivated by Matthaeus et al. 2007, spectra with q > 1 are now favored. It has to be emphasized that spectrum (2.20) is a particular case satisfying the various physical requirements discussed in general terms by Matthaeus et al. (2007). In the latter article one can find alternative forms of spectra satisfying those requirements.

Two-Component Turbulence
Above we have considered the slab and the two-dimensional models as examples. It is often assumed that turbulence in the solar wind can be approximated by a two-component model in which a superposition of slab and two-dimensional modes is considered. This type of turbulence approximation is motivated and supported by analytical investigations, simulations, and solar wind observations (see, e.g., Matthaeus et al. 1990Zank and Matthaeus 1993;Oughton et al. 1994;Shaikh and Zank 2007;Hunana and Zank 2010;). In Fig. 4 we show magnetic correlations obtained via solar wind observations. The data shown there was published in Matthaeus et al. (1990) and is based on magnetometer measurements made aboard the International Sun/Earth Explorer 3 (ISEE-3) spacecraft also known as International Cometary Explorer (ICE). We can clearly see the cross structure supporting the idea of having a superposition of slab and two-dimensional modes. Therefore, we superpose the magnetic field associated with slab modes and the field associated with two-dimensional modes and assume that these two are uncorrelated. Therefore, within the slab/2D composite model, the components of the spectral tensor are written as where the slab and two-dimensional tensors are given by Eqs. (2.13) and (2.19), respectively. For the spectra we can still employ Eqs. (2.14) and (2.20). The quantities δB 2 slab and δB 2 2D therein are the magnetic energy densities associated with slab and two-dimensional modes. Often we replace them by the so-called slab fraction δB 2 slab /δB 2 and the two-dimensional fraction δB 2 2D /δB 2 , respectively. Here we have used the total turbulent magnetic field δB = δB 2 slab + δB 2 2D . Realistic values in the solar wind at 1 AU heliocentric distance should be δB 2 slab /δB 2 = 0.2 and δB 2 2D /δB 2 = 0.8 as discussed in Bieber et al. (1994) as well as Bieber   et al. (1996). Compared to the slab model we now generated turbulence with transverse structure because the total turbulent magnetic field is now (2.24) Figure 5 shows flux surfaces for two-component turbulence.

Noisy Slab Turbulence
So far we have considered models with reduced dimensionality and superpositions thereof. One can also find a variety of full three-dimensional models in the literature. One example is the so-called noisy slab model. The idea is to consider the slab model and include some noise so that the model becomes three-dimensional. This idea was originally developed in Weinhorst and Shalchi (2010) but was also used in Ruffolo and Matthaeus (2013) in the context of two-dimensional turbulence (see next subsection). The noisy slab model, as discussed in the following, was formulated the first time in Shalchi (2015a). In this case the components of the spectral tensor are given by where ⊥ denotes a characteristic scale describing the correlation of magnetic fields across the mean field. In Eq. (2.25) we have employed the Heaviside step function (x) which is defined so that For the model spectrum g slab (k ) one can employ the same form as used above for pure slab turbulence (see, e.g., Eq. (2.14) of this review). Due to the step function in Eq. (2.25), the wave vectors are no longer aligned perfectly parallel with respect to the mean field. However, the idea here is that the parameter ⊥ is large so that the width of the distribution of wave numbers in the perpendicular direction is small. The noisy slab model is useful to study the importance of transverse complexity in the theory of perpendicular diffusion as shown throughout this paper. In the limit ⊥ → ∞ one can obtain the usual slab model from the more general noisy slab model. Mathematically this is due to the function ⊥ (1 − k ⊥ ⊥ ) behaving like a Dirac delta for ⊥ → ∞. Furthermore, the weak transverse complexity allows one to study field line random walk and field line separation by using higher-order perturbation theory and to drop some approximations usually required in this field (see Shalchi 2019b).

Noisy Reduced MHD Turbulence
The idea of broadening a model with reduced dimensionality can also be used for twodimensional turbulence. Ruffolo and Matthaeus (2013) called this model the Noisy Reduced MHD (NRMHD) model. In this case the components of the magnetic correlation tensor have the following form where we have employed the Heaviside step function (x) as defined in Eq. (2.26). The function g 2D (k ⊥ ) is the usual two-dimensional model spectrum discussed before in this review. In transport theory the NRMHD model is useful in order to study particle transport in large Kubo number turbulence. Considering the limit → ∞ would restore the twodimensional model discussed in Sect. 2.3.

The Gaussian Correlation Model
Above we discussed different models for magnetic turbulence. Those models were either models with reduced dimensionality or extensions thereof. In the past full three-dimensional models were discussed as well. The simplest model would be the isotropic model but there are also more complicated models such as the anisotropic models used by Lerche and Schlickeiser (2001) and Zimbardo et al. (2006). An alternative model for three-dimensional turbulence is the Gaussian model. This is not very realistic in the context of astrophysics and space science but is often used in theoretical studies of fusion plasmas (see, e.g., Spatschek 2008 for a review). Special care is required since the Gaussian model must not be used to study parallel transport. As will be demonstrated in Sect. 4 of this review, the inertial range, which is neglected in the Gaussian model, controls pitch-angle scattering and parallel spatial diffusion for intermediate particle energies. Therefore, the Gaussian model is unrealistic in the context of particle transport in space plasmas. However, due to its three-dimensional structure, it can be used for the analytical exploration of field line random walk and it can be used to determine the relation between the parallel and perpendicular diffusion coefficients if non-linear theories are considered. In the following we use a generalized Gaussian decorrelation model which contains a general behavior in the energy range of the spectrum in the perpendicular direction. For the incompressible case the components of the spectral tensor are given by Eq. (2.10). For the model spectrum we now consider the form with the normalization function The function E(q) can easily be obtained by combining Eq. (2.28) with the normalization condition (2.16). In Eq. (2.28) we have used the energy range spectral index q, the total magnetic field strength in the x-direction δB x , as well as the bendover scales in the parallel and perpendicular direction and ⊥ , respectively. For the case q = 3 we can easily recover the Gaussian correlation model which was used elsewhere (see, e.g., Neuer and Spatschek 2006). In this case we have for the normalization function.

Goldreich-Sridhar Turbulence
In the famous paper Goldreich and Sridhar (1995), the authors explored strong Alfvénic turbulence. Their investigations focused on the case where oppositely directed Alfvén waves carry equal energy fluxes. As explicitly stated in Goldreich and Sridhar (1995), this precludes the application to the solar wind in which the outward flux significantly exceeds the ingoing one. Therefore, we use the ideas of Goldreich and Sridhar to approximate interstellar turbulence in the current review whereas interplanetary turbulence is usually modeled via two-component turbulence. Furthermore, the work of Goldreich & Sridhar is based on incompressible magnetohydrodynamics. In the context of field line and particle transport, this corresponds to the case δB z = 0 considered before in this review. Therefore, our considerations are based on the spectral tensor given by Eq. (2.10). In the following our aim is to develop a spectral function g(k , k ⊥ ) which exhibits all the aspects of turbulence described in Goldreich and Sridhar (1995). In order to do this we perform the following steps: 1. The inertial range of the spectrum exhibits a so-called critical balance between linear wave periods and non-linear turnover time scales. In terms of parallel and perpendicular wave numbers, this critical balance condition can be written as where is the outer scale of turbulence. In order to estimate the spectrum g(k , k ⊥ ) we take into account the critical balance condition via (2.32) The exponential here was not used in the original work of Goldreich and Sridhar (1995) but it is one of the simplest ways to incorporate the critical balance condition and was originally used in Cho et al. (2002). 2. The one-dimensional spectrum is assumed to be proportional to k −5/3 ⊥ corresponding to Kolmogorov (1941) but in the perpendicular direction. Therefore, we assume and determine the exponent ν by requiring (2.34) After some straightforward algebra we derive ν = −10/3. 3. The work of Goldreich & Sridhar focused on the inertial range of the spectrum. The parameter used above can be used as a cut-off in the spectrum. This can be incorporated via a Heaviside step function so that with the normalization constant C. 4. The remaining step is the normalization of the spectrum. For the incompressible case the normalization condition yields Finally we arrive at a spectral function in accordance with the work of Goldreich and Sridhar (1995), namely Although this might be a useful form for certain applications, the main problem here is that the critical balance condition was derived for the inertial range of the spectrum. Field line and perpendicular particle transport, however, are controlled by the large scales of the energy range. Therefore, it is questionable whether the critical balance condition really matters in this field of research. Nevertheless, we still consider the spectrum given by Eq. (2.38) as an example throughout this review. Cho et al. (2002) have proposed a form for the spectral tensor which is almost identical compared to the form discussed above. However, they considered the compressible case (δB z = 0) and, thus, the components of the spectral tensor are now given by where n, m = x, y, z. Since we now consider the compressible case, this is not longer a special case of Eq. (2.10). Furthermore, Cho et al. (2002) altered the spectrum slightly and proposed The factors in this model were chosen to satisfy the normalization constraint (2.16) which includes, in this case, a non-vanishing δB z . Furthermore, we have used the ratio of the turbulent and the mean field E B ≡ δB/B 0 also known as the Alfvénic Mach number. Problematic in Eqs. (2.38) and (2.40) is that the energy containing range is neglected because there in no turbulent magnetic field for k ⊥ < 1/ ⊥ . The large scales of the energy range are crucial in the theory of field line random walk and perpendicular transport of energetic particles and should not be neglected. However, one can combine Eq. (2.40) with some ideas discussed above in the context of two-dimensional turbulence. Thus, in Shalchi (2013a) the following generalization was presented which is only valid for s = 5/3. Compared to Eq. (2.40), the latter spectrum takes into account a more general behavior at large scales corresponding to the energy range. The normalization function D(s, q) is defined in Eq. (2.22) and the parameter q corresponds again to the energy range spectral index. Sun and Jokipii (2011) used a very similar spectrum in their test-particle simulations but set q = 0.
In order to understand the relation between Eq. (2.41) and other spectra used in the literature, one can integrate the spectrum over all parallel wave numbers and we find the reduced perpendicular spectrum (2.42) Using spectrum (2.41) in the latter equation, and after some straightforward algebra, we deduce in agreement with the Shalchi & Weinhorst spectrum 7 given by Eq. (2.20) if this is also integrated over all parallel wave numbers.

Further Turbulence Effects
Above, we have focused on incompressible cases where δB z = 0. To explore compressible turbulence can become relevant in scenarios where the turbulent field exceeds the mean field. In particular in strong isotropic turbulence, the perpendicular motion of particles can be different and perpendicular and parallel diffusion coefficients can become identical in some limits (see, e.g., Shalchi and Dosch 2009;Plotnikov et al. 2011, as well as Subedi et al. 2017 for theoretical approaches to describe particle transport in isotropic turbulence). A further discussion of this matter can be found in Sect. 12.2.3 of this review. Furthermore, the effect of intermittency (effects of dynamically produced coherent structures) can be included as it was done in the papers by le Roux (2011) and Pucci et al. (2016). It was shown that intermittency can have an influence on parallel and perpendicular diffusion coefficients of energetic particles.

Magnetic Correlations in Configuration Space
So far we described magnetic turbulence based on spectral tensors because the components of those tensors enter transport theories as we shall see later in this review. Of course one can go back to Sect. 2.1 and can try to compute the two-point correlation functions in configuration space. As an illustrative example we consider the slab model. By combining Eq. (2.13) with Eq. (2.8) we find (2.44) With spectrum (2.14) this becomes where we solved the corresponding integral with the help of Gradshteyn and Ryzhik (2000) and used modified Bessel functions of second kind K ν (x). A strong simplification can be achieved if we consider s = 2 which is close to the more realistic value s = 5/3. In this case the parallel wave number integral in the first line of Eq. (2.45) yields corresponding to an exponential decay of the correlation function. Although this is a special case, the model (2.46) can be useful in order to make simple estimations and to understand the meaning of certain quantities such as the integral scale discussed in the next subsection. An analytical treatment of correlation functions was presented in Shalchi (2008) where the reader can also find results for other turbulence models. Detailed studies of two-point correlations in the solar wind via measurements can also be found in the literature (see, e.g., ).

Integral Scales and the Ultra-Scale
There are several important length scales in the theory of magnetic turbulence. We have already used two of them, namely the parallel and perpendicular bendover scales and ⊥ , respectively. In the following we discuss so-called integral scales in the different directions of space as well as the ultra-scale.
If we deal with axi-symmetric and magnetostatic turbulence, we can define the parallel integral scale L via where we used z = 0 for the initial position. We can understand the integral scale as a characteristic length scale over which the turbulent magnetic field decorrelates. Therefore, the parameter L is often called the parallel correlation length. This becomes clearer if we consider Eq. (2.46) as an example. In this case the magnetic field decorrelates over the characteristic scale . Using this in Eq.
(2.47) shows that L = and that L is indeed a characteristic scale for the decorrelation of the magnetic field. However, in the general case, these two scales are not equal as shown below. Using the Fourier representation (2.8) and assuming static turbulence allows us to derive from Eq. (2.47) (2.48) A very useful relation involving Dirac's delta is given by (see, e.g., Zwillinger 2012) which is frequently used throughout this review. Using this in Eq. (2.48) allows us to derive (2.50) For slab turbulence, for instance, we need to employ Eq. (2.13) to derive δB 2 x L = 2π 2 g slab (k = 0). (2.51) For the spectrum given by Eq. (2.14), the parallel integral scale is then meaning that the integral scale and the bendover scale are directly proportional to each other but not identical. For s = 2, however, we have C(s = 2) = 1/(2π) and the two scales becomes equal. In Table 2 we show the relations between the scales and L for the turbulence models considered in this review together with other scales. These relations will be useful because the scales listed there control the diffusion coefficients of field lines and energetic particles in some limits as will be shown later in this review.
Similar compared to the previous paragraph we now consider the perpendicular integral scale L ⊥ which can be defined via 8 for the axi-symmetric case. Using again the Fourier representation (2.8) yields Since our aim is to rewrite the latter equation by using cylindrical coordinates we rename k y → k ⊥ to find (2.57) As an example we now evaluate this for the two-dimensional model defined via Eq. (2.19). We derive (2.58) By using spectrum (2.20) and the integral transformation x = k ⊥ ⊥ , we obtain (2.59) The latter integral can be solved for q > 0 (see, e.g., Gradshteyn and Ryzhik 2000) . (2.60) Using this in Eq. (2.59) and after employing Eq. (2.22), we find for the perpendicular integral scale As for the slab model, the integral scale is directly proportional to the bendover scale. However, the integral scale defined here is only finite for q > 0. One could argue that this excludes spectra with q ≤ 0. Alternatively, one could allow such spectra but then one needs to introduce a cut-off in the two-dimensional spectrum. A detailed discussion of this matter can be found in Matthaeus et al. (2007).
In the theory of wandering magnetic field lines, another important scale occurs, namely the so-called ultra-scale L U . The latter scale is defined via (see, e.g., Matthaeus et al. 2007) In this review the ultra-scale will be important in the non-linear/Bohmian transport regime of magnetic field line random walk (see Sects. 3.7-3.9). For the two-dimensional turbulence model, defined via Eq. (2.19), the ultra-scale becomes For the spectrum (2.20), with δB 2 2D = 2δB 2 x , and by using the integral transformation x = k ⊥ ⊥ , we deduce (2.64) Fig. 6 The fundamental scales of magnetic turbulence for the two-dimensional model as given by Table 2. Shown are the ratios L U /L ⊥ (dotted line), L U / ⊥ (solid line), and L ⊥ / ⊥ (dashed line) versus the energy range spectral index q. Here we have used the bendover scale ⊥ , the integral scale L ⊥ as well as the ultra-scale L U The latter integral can be solved for q > 1 (see, e.g., Gradshteyn and Ryzhik 2000) and we obtain where we have used Eq. (2.22) and the relation (see, e.g., Abramowitz and Stegun 1974) to get rid of the gamma functions. Obviously the ultra-scale is directly proportional to the bendover scale ⊥ but also depends on the inertial range spectral index s and the energy range spectral index q. It has to be emphasized that the ultra-scale is only finite for the spectrum used here if q > 1 corresponding to an increasing spectrum in the energy range (see, e.g., Fig. 3 of this review). Often one uses the values s = 5/3 and q = 3 leading to L U = ⊥ / √ 3. In this case the ultra-scale is slightly smaller than the perpendicular bendover scale (see Fig. 6 for other values of q).

The Random Walk of Magnetic Field Lines
In magnetized plasmas magnetic field lines are stochastic curves due to the turbulent aspects of magnetic fields described in the previous section. In the theory of field line random walk (FLRW) we study the statistics of magnetic field lines by computing the mean square displacement of different field line realizations. In principle on can also obtain the field line distribution function indicating the probability to find a magnetic field line at a certain position in space. The concepts of field line random walk and field line separation are depicted in Fig. 7 of this review. To study the random walk of magnetic field lines is important because the properties of magnetic field lines in turbulence are often controlling the perpendicular Fig. 7 The concepts of field line random walk (FLRW) and field line separation. In the theory of FLRW we study the statistics of a single field line described by the distance of the unperturbed field line and the turbulent field line In the theory of field line separation, on the other hand, we study the distance s between two stochastic magnetic field lines. The latter effects occurs only in turbulence with transverse structure whereas the FLRW occurs also in slab turbulence transport of energetic particles as will be demonstrated later in this review article. Most of the following presentations are based on the work of Taylor and McNamara (1971), Kadomtsev and Pogutse (1974), Salu and Montgomery (1977), and in particular on Matthaeus et al. (1995).

Fundamental Equations
Magnetic field lines in turbulence are stochastic curves which need to be described by methods of statistical physics. Therefore, we are not able to compute individual field lines but rather work with field line mean square displacements ( x) 2 where x = x(z) − x(0). The mean square displacement increases with increasing distance z from the initial position since there is an increase of the uncertainty to find the field line at a certain position of space. Furthermore, a running field line diffusion coefficient can be defined via It has to be emphasized that in the theory of random walking magnetic field lines, the variable is position z and not time as in particle diffusion theory. Furthermore, the field line diffusion coefficient has length units whereas usual diffusion coefficients have the dimension [L] 2 / [T ]. Especially in static turbulence the concept of field line random walk (FLRW) could be confusing because there is no motion of field lines as the term diffusion suggests. Instead diffusion really means an increase of the uncertainty to find the field line at a certain point of space. This also means that for the axi-symmetric case there is only one diffusion coefficient namely the one defined in Eq. (3.1). In particular there is no parallel or perpendicular diffusion coefficient as in particle transport theory. The whole concept of FLRW suggests that it would be more appropriate to talk about field line meandering rather than field line diffusion. The quantity defined in Eq. (3.1) is called a running diffusion coefficient because, in general, it is a function of position z. If we call the considered process normal or Markovian diffusive. In this case the mean square displacement is directly proportional to |z|. In general we can assume

2)
and the field lines can be characterized as listed in Table 3. Non-diffusive transport is also called anomalous diffusion and corresponds to processes with α = 1. In the context of field line and energetic particle transport anomalous diffusion has received attention in the literature (see, e.g., Zimbardo et al. 1995Zimbardo et al. , 2006Zimbardo et al. , 2012Zimbardo 2005;Zimbardo 2007, 2009a,b;Shalchi and Kourakis 2007a;Perri et al. 2015). In the following we review different analytical theories which were developed in the past in order to compute field line diffusion coefficients and mean square displacements. We consider a scenario where the turbulent magnetic field is perpendicular with respect to the mean field. In the solar wind, for instance, magnetic fluctuations have been found to be about 90% transverse in terms of fluctuation energy (see Belcher and Davis 1971). Based on MHD equations  found theoretically that large scale magnetic fields can organize turbulence fluctuations to be mostly transverse. Therefore, the total magnetic field vector is given by where we assumed that the mean field B 0 is constant and points into the z-direction. For the case considered here, the field line equations have the form dx = δB x B 0 dz and dy = δB y B 0 dz (3.5) and the spectral tensor is given by Eq. (2.10). After integrating the first field line equation, the displacement in the x-direction can be written as the following integral equation where we used z = 0 as initial position. Thus, we find for the field line mean square displacement where δB * x is the complex conjugate magnetic field. By employing the Fourier representation (2.3), we derive We can clearly see that a complicated correlation function occurs on the right hand side of this formula. This correlation function is not trivial since the vectors x(z ) and x(z ) denote the stochastic magnetic field line as a function of position. The field lines depend on the magnetic field and, thus, field lines and field vectors are correlated in the general case. To proceed, we distinguish between two different cases: 1. The slab model: here the exponential function in the brackets . . . depends only on the variable z. In this case, we have k · x(z) = k z, and the ensemble average in Eq. (3.8) can be evaluated without using further assumptions or approximations. 2. Models with transverse complexity: for non-slab models the exponential function in the brackets . . . can depend on all components of the field line vector x(z) = (x(z), y(z), z).
In this case, the field line equation (3.8) is a non-linear equation. Thus, further assumptions and approximations have to be used in order to simplify Eq. (3.8).
Before we investigate these two cases in more detail, we discuss the importance of the socalled Kubo number and explore the limit z → 0.

The Role of the Kubo Number in the Theory of Field Line Diffusion
The behavior of magnetic field lines is linked to the spectral tensor discussed before in this review. However, one can make some general estimations without specifying this tensor. In the following we assume that the magnetic field is given by where ⊥ and are characteristic length scales for the decorrelation of the magnetic field such as the bendover scales. Therewith the field line equation (3.5) can be written as where we have used the Kubo number (see Kubo 1963) If we express all positions along the mean field in terms of and all positions across the mean field in terms of ⊥ , the Kubo number K is the only quantity controlling the FLRW. Therefore, we conclude that the Kubo number is crucial in the quantitative description of field line transport.

The Initial Free-Streaming Regime
Equation (3.8) provides a formula for the field line mean square displacement which can be employed for arbitrary position z. In the limit z → 0 we also expect x, y → 0 and, thus, we can derive from Eq. (3.8) (3.12) By using Eq. (2.6), corresponding to the assumption of homogeneous turbulence, this can be simplified to where we have used Eq. (2.15) as well. If we compare this with the general form given by Eq. (3.3) we find α = 2. According to Table 3 this corresponds to ballistic field line random walk. In the following we refer to this limit as the initial free-streaming regime. This behavior is found for short distances regardless what the form of the spectral tensor is. The initial ballistic regime can also be seen in Fig. 8 of this review.

Field Line Random Walk for Slab Turbulence
The slab model is a very special model of turbulence and it is not sufficient to approximate real turbulence in most cases. In the following it is used because it allows for an exact description of FLRW. This is not the case for other turbulence models. For pure slab geometry Eq. (3.8) becomes (3.14) By employing again Eq. (2.6) and by using the slab model given by Eq. (2.13), we find for the mean square displacement of the field lines The two integrals over z and z can easily be evaluated and we find (3.16) Taking the second derivative of the latter equation yields (3.17) Now we compute the running diffusion coefficient defined in Eq. (3.1). However, we concentrate on the limit z → ∞. Therefore, we integrate Eq. (3.17) over all z to find where d FL (z = 0) = 0 due to Eq. (3.13). To continue we employ the relation (2.49) to find the following result It is important to interpret this result in the right way. In reality there should not be any turbulence for k = 0. However, in reality the limit z → ∞ is not permitted either. Therefore, the limit z → ∞ really means that we consider large but not too large distances whereas k = 0 corresponds to the largest scales where we find turbulence. A detailed discussion of this matter, in the context of particle transport, can be found in Shalchi (2005b).
For the model spectrum of Eq. (2.14) the field line diffusion coefficient becomes (3.20) Characteristic here is the scaling with δB 2 slab /B 2 0 . Alternatively, we can go back to Eq. (3.19) and combine this more general form with the parallel integral scale derived in Eq. (2.51). This allows us to write which does not require to specify the slab spectrum. It has to be emphasized that the calculations presented here for slab turbulence did not require any assumptions or approximations and, thus, the obtained results are exact. Equation (3.21) shows the importance of the parallel integral scale in the theory of random walking magnetic field lines.

Quasi-Linear Theory of Field Line Random Walk
For non-slab models, Eq. (3.8) is a non-linear equation. One way of simplifying this is the application of perturbation theory. In the theories of field lines and energetic particles, this approach is usually called quasi-linear theory (QLT) which was developed by Jokipii (1966) and Jokipii and Parker (1969). Within QLT, we replace the field lines on the right hand side of Eq. (3.8) by the unperturbed field lines x(z) = y(z) = 0. Therewith, Eq. (3.8) becomes where we assumed again homogeneous turbulence. This form is very similar compared to Eq. (3.15) and, thus, we employ the same steps as before. The integrals over z and z can easily be solved and we find (3.23) Again we consider the second derivative of the latter equation to obtain The running diffusion coefficient is defined in Eq. (3.1). Therefore, to obtain this quantity, we integrate Eq. (3.24). After employing again relation (2.49) we find in the limit z → ∞ the following result Quasi-linear theory can be understood as lowest order perturbation theory (see, e.g., Shalchi 2019b for a higher order perturbation approach). Of course the question arises what the perturbing parameter is in this case. According to Eq. (3.10) it is the Kubo number. If the latter number is small, the field lines are close to the unperturbed field line x(z) = y(z) = 0. This means that the smaller the Kubo number is, the more accurate QLT should be. For slab turbulence, for instance, we have ⊥ → ∞ and, thus K = 0 meaning that quasi-linear theory is exact. For two-dimensional turbulence, on the other hand, we have → ∞ and, thus K = ∞. One would expect that in this case QLT fails completely. Therefore, non-linear tools need to be developed in order to describe the random walk of magnetic field lines for cases where the Kubo number is not small.

Corrsin's Independence Hypothesis
Equation (3.8) is exact and, thus, we use it as starting point for developing a non-linear approach. Problematic here is the complicated correlation function involving magnetic field components as well as the field line x(z) itself. A strong simplification can be achieved by employing Corrsin's independence hypothesis (see Corrsin 1959) which can be written as (3.26) The widespread use of the Corrsin hypothesis in transport theory is discussed in detail in Tautz and Shalchi (2010). In the theory of FLRW this approximation was already used by Lerche (1973) but at this point it was not called the Corrsin hypothesis. With the Corrsin approximation, Eq. (3.8) becomes corresponding to a significant simplification. Assuming homogeneous turbulence enables us to use Eq. (2.6) and, thus, we obtain (3.28) To further simplify this, we compute the second derivative of the latter equation. To determine the first derivative, we employ the Leibniz integral rule which can be written as (3.31) To continue we assume that the characteristic function depends only on the positiondifference meaning that (3.32) We now employ the integral transformation ξ = z − z to derive (3.33) Using the latter relation in Eq. (3.31) yields where we have used x(z ) = x(z ) − x(0). Now it is straightforward to calculate the second derivative. We find For a Gaussian distribution of field lines, we have 9 where we assumed again axi-symmetry. Finally we find the following second-order ordinary differential equation for the field line mean square displacement This equation corresponds to the integral equation derived in Shalchi and Kourakis (2007b) by using a slightly different approach based on the same set of assumptions. In Shalchi and Kourakis (2007c) it was shown that the FLRW is mainly controlled by the large scales of the energy range and in Shalchi and Qin (2010) the equation derived here was successfully tested by comparing results with simulations performed for two-component turbulence. In Figs. 8 and 9 we show examples for a comparison of analytical theory and simulations confirming the validity of the non-linear tools discussed above. In the following Eq. (3.37) is discussed in more detailed by considering special cases and limits.

FLRW in Two-Dimensional Turbulence
Equation (3.37) can be used for an arbitrary spectral tensor P nm (k). In the following we consider the two-dimensional model as an example and follow the considerations presented in Shalchi (2011c). In this case we employ Eq. (2.19) and therewith Eq. (3.37) becomes where we have used σ = ( x) 2 for convenience. Now we multiply the latter equation by the derivative σ to find (3.39) 9 If the field lines obey a Gaussian statistics, the field line distribution function is for the axi-symmetric case. The characteristic function is then computed via exp (ixk x + iyk y ) = dx dy exp (ixk x + iyk y )f (x, y; z) yielding Eq. (3.36).

Fig. 8
The running field line diffusion coefficient normalized with respect to the slab bendover scale for two-component turbulence versus the parallel position z/ . Shown is the analytical result in the limit z → ∞ (dotted line), the result obtained by solving Eq. (3.37) numerically (dashed line), and the simulations performed by Shalchi and Qin (2010) which is represented by the solid line. All results shown here were obtained for an energy range spectral index of q = 1.5, ⊥ / = 0.1, δB 2 slab /B 2 0 = 0.2 and δB 2 2D /B 2 0 = 0.8

Fig. 9
The field line diffusion coefficient for NRMHD turbulence (see Sect. 2.6). Shown is κ F L / ⊥ versus the magnetic field ratio δB/B 0 for ⊥ = .
The solid line represents the analytical result obtained by solving Eq. (3.47) numerically and the dots represent the simulations performed by Shalchi and Hussein (2014) This can be written as Next we integrate this and use σ (z = 0) = 0 to obtain For large distances the exponential will damp out and we find normal diffusion with the field line diffusion coefficient given by With the ultra-scale L U defined via Eq. (2.63), this result can be written as showing the relevance of the ultra-scale in the theory of FLRW. For spectrum (2.20) the ultra-scale is given by Eq. (2.66) leading to the following field line diffusion coefficient The latter formula provides the field line diffusion coefficient for pure two-dimensional turbulence. It has to be emphasized that a diffusion approximation was not used in order to derive this formula. An alternative approach, where the diffusion approximation is used, will be presented in the next subsection.

The Diffusion Approximation
Equation (3.37) allows for a z-dependent description of FLRW by solving this second-order differential equation either analytically or numerically. However, this can be difficult or at least time-consuming. A strong simplification can be achieved by combining Eq. (3.37) with a diffusion approximation. In the context of FLRW diffusion approximation means that we assume that field lines behave diffusively for all distances z. This means in particular that the initial free-streaming regime is ignored. We also loose the ability to describe FLRW for nondiffusive cases. Mathematically the diffusion approximation corresponds to the assumption that Using this on the right hand side of Eq. (3.37) yields After integrating this over all z, we obtain the following equation for the field line diffusion coefficient as derived in Matthaeus et al. (1995). This results can also be written as which is very similar compared to equation (13) of Kadomtsev and Pogutse (1979). Equation (3.47) can still be evaluated for any given spectral tensor P nm but it is much simpler compared to Eq. (3.37). However, it is less accurate because it is based on the diffusion approximation and less general because it only works for cases where field line random walk is indeed diffusive. Furthermore, Eq. (3.47) can also not describe the initial free-streaming regime as well as the turnover to the diffusive regime.
In the following we employ the slab/2D model as an example. With Eqs. (2.13), (2.19), and (2.23), Eq. (3.47) can easily be evaluated and we find (3.49) To continue we use 10 (3.51) With the field line diffusion coefficient obtained for pure slab turbulence κ slab given by Eq. (3.19), and with the parameter this becomes corresponding to a quadratic formula for the parameter κ FL . It can easily be solved by where we excluded the negative solution. For pure slab turbulence, we have κ 2D = 0 and we get the expected result κ FL = κ slab . For pure two-dimensional geometry, on the other hand, we have κ slab = 0 and we find κ FL = κ 2D . Thus, the parameter κ 2D , defined in Eq. (3.52), can be identified with the diffusion coefficient for pure two-dimensional turbulence. Comparing Eq. (3.52) with Eq. (3.42) leads to the conclusion that there is a factor √ 2 between both diffusion coefficients. Equation (3.42) was derived without employing the diffusion approximation and, thus, it is more accurate. Keeping in mind that √ 2 ≈ 1.4 leads to the assumption that the diffusion approximation should work well in the theory of FLRW but small discrepancies can occur. The factor √ 2 resulting from the diffusion approximation, was originally discovered in Shalchi (2011c) for two-dimensional turbulence but it was also found in Shalchi (2019b) for small Kubo number turbulence. 10 It is not possible to prove this relation without delving deeply into the theory of delta distributions. However, If the ultra-scale, given by Eq. (2.63), is used in Eq. (3.52), we find For spectrum (2.20) the ultra-scale is provided by Eq. (2.66) leading to the following field line diffusion coefficient First we note that the field line diffusion coefficient is directly proportional to the bendover scale ⊥ . Furthermore, there is also a direct proportionality to the magnetic field ratio which is different compared to the quasi-linear scaling.
More theoretical work has been done in the recent years. Ruffolo et al. (2006), for instance, dropped the assumption of axi-symmetry and developed a non-linear theory for FLRW based on the ideas of Matthaeus et al. (1995). Furthermore, one can study magnetic field lines in different turbulence configurations. Besides the often discussed slab, twodimensional, and two-component models, one can study FLRW in full three-dimensional turbulence. Ruffolo and Matthaeus (2013), for instance, studied field lines in noisy reduced magnetohydrodynamic turbulence whereas Shalchi and Kolly (2013) explored FLRW in Goldreich-Sridhar turbulence and Sonsrettee et al. (2016) studied the FLRW in isotropic turbulence with varying mean field. Ruffolo et al. (2004) and Shalchi (2019b) have investigated the separation of two magnetic field lines showing that this process has to be described by non-linear tools in the general case.

The Non-Linear Regime of FLRW
From Eq. (3.47) one can easily recover quasi-linear theory by considering the limiting process κ FL k 2 ⊥ → 0 together with the relation (see, e.g., Zwillinger 2012) (3.57) As discussed before the resulting quasi-linear scaling should be valid in the small Kubo number limit and the field line diffusion coefficient is then given by Eq. (3.21). The question naturally arises what one would obtain for large Kubo numbers K. To explore this further, we now consider the case κ FL k 2 ⊥ → ∞ which allows us to simplify Eq. (3.47) to We like to emphasize that at this point we did not specify the spectral tensor P nm (k) but we can easily see that Eq. (3.58) is identical compared to Eq. (3.52) if the two-dimensional model is used. This is no surprise since the two-dimensional model corresponds to the case of K = ∞. In general one can use the ultra-scale defined via Eq. (2.62) to write Eq. (3.58) as The latter result is valid for arbitrary turbulence but requires that we consider large Kubo numbers. We can see that the diffusion coefficient is directly proportional to the ultra-scale and the magnetic field ratio. In this non-linear regime the field line diffusion coefficient is, therefore, directly proportional to the Kubo number K whereas in the quasi-linear regime it is proportional to K 2 . The regime obtained here for large Kubo numbers is also known as the Bohmian regime of FLRW or as the Kadomtsev and Pogutse (1979) limit. The limit discussed here is sometimes described as being incorrect. The percolation theory developed by Gruzinov et al. (1990) was applied to the FLRW by Isichenko (1991) and a different scaling compared to Eq. (3.59) was obtained by these authors (see Sect. 3.10 for more details). Shalchi and Qin (2010) performed simulations for field line random walk and they confirmed the validity of the non-linear or bohmian regime described above (see, e.g., Fig. 8 of this review for an example). Furthermore, Fig. 9 shows a comparison between analytical theory and simulations. As demonstrated, there is very good agreement for small, intermediate, and large values of the Kubo number. It seems that the Corrsin approximation, which is the only approximation used during the derivation of Eq. (3.37), works well. We conclude that the percolation model probably applies to the limit of extremely high Kubo numbers as already stated in Ghilea et al. (2011).

Scaling Laws in the Theory of FLRW
The results derived so far can also be obtained by using simple scaling arguments (see, e.g., Hauff et al. 2010). First we consider the definition of a field line diffusion coefficient via Eq. (3.1). The latter diffusion coefficient consists of a perpendicular position squared divided by a parallel position. Therefore, a scaling independent diffusion coefficient can be defined as κ FL / 2 ⊥ . Motivated by the considerations presented above, we assume that this diffusion coefficient depends only on one single parameter and that is the Kubo number K. Assuming a power-law dependence leads to the following proportionality (3.60) Using Eq. (3.11) to replace the Kubo number therein yields the following form For certain values of γ we can recover the formulas derived before. For γ = 2, for instance, we obtain the quasi-linear scaling and for γ = 1 the Bohmian regime. If we consider models with reduced dimensionality, this can be simplified further. For slab turbulence, for instance, the field line diffusion coefficient must not depend on the perpendicular scale ⊥ . Therefore, we need to set γ = 2 automatically leading to the correct quasi-linear scaling κ FL ∝ δB 2 x /B 2 0 . For two-dimensional turbulence, on the other hand, the field line diffusion coefficient must not depend on the parallel scale . Thus, we need to set γ = 1 in Eq. (3.61) automatically leading to κ FL ∝ ⊥ δB x /B 0 corresponding to Bohm diffusion. Note that L U ∝ ⊥ in most cases. It is interesting that we automatically arrive at the correct results just by performing simple considerations based on dimensional arguments.
The percolation theory discussed in the previous subsection can also be obtained from Eq. (3.61) by setting γ = 0.7 leading to (3.62) Therefore, we conclude that the scaling argument applied here does indeed contain most cases of FLRW discussed in the literature. Of course such dimensional arguments do not provide the numerical factors in the equations for the field line diffusion coefficient nor do they work for anomalous transport.

Parallel Transport of Energetic Particles
In the current review article we focus on the transport of energetic particles across the mean magnetic field. However, perpendicular transport strongly depends on the parallel particle motion as will be shown in Sects. 5-9. Therefore, it is necessary to discuss parallel diffusion but we only focus on basic ideas because parallel transport itself is very complicated (see Shalchi 2009a).

Equations of Motion and Unperturbed Orbits
As a first step we briefly discuss the equation which describes the motion of energetic particles through turbulence. Thereafter, we solve this equation for the unperturbed case. This is a simple classical mechanics problem but the solution is briefly summarized in order to explain the used notation. Unperturbed orbits and the quantities used below are frequently used throughout this review. The fundamental equation describing the motion of charged particles through purely magnetic turbulence is the Newton-Lorentz equation where the total magnetic field is given by B(x) = δB(x) + B 0 e z as in previous sections. Furthermore, we have used the electric charge of the energetic particle q, the particle velocity v, the relativistic momentum p, as well as the speed of light c. In Eq. (4.1), as well as during the whole review paper, we employ cgs or Gaussian units rather than SI units. First we consider the case that there is no turbulence and that the mean magnetic field is constant and points into the z-direction. In this case Eq. (4.1) simplifies drastically to where we have replaced the relativistic momentum by p = mγ v and used the fact that the speed of the particle v and therewith the Lorentz factor γ are constant if there are no electric fields. Furthermore, we used the gyrofrequency given by 11 Equation (4.2) can be written as the following partly coupled system of linear differential equationsv (4.4) The latter system can easily be solved by where we have used the pitch-angle cosine μ = v z /v which is constant in the unperturbed case. The pitch-angle indicates the angle between velocity vector v and the mean magnetic field B 0 = B 0 e z . Furthermore, we have used the gyrophase (t) = 0 − t and the initial gyrophase 0 . Obviously the parallel motion occurs with constant velocity. The particle trajectory can easily be obtained by integrating Eq. (4.6) over time and we derive It is convenient to define the parameters as well as (4.8) Therewith, the first two lines of Eq. (4.7) can be written as corresponding to the equation of a circle with center (x c , y c ) and radius R L 1 − μ 2 . The orbits found here correspond to a helical motion across the mean magnetic field. Therefore, R L 1 − μ 2 is the so-called gyroradius or Larmor radius whereas the parameter R L is this radius at μ = 0. However, since R L is a characteristic quantity in transport theory, even if there is turbulence, we just call R L the Larmor radius although this is not strictly correct.
Of course we are more interested in the case that there is turbulence. From Eq. (4.1) we can easily derive the corresponding equations of motion. Since we focus on parallel transport in the current section, we only consider the z-component of this equation which can be written asμ (4.10) Therefore, the pitch-angle cosine μ is no longer conserved if the particle moves through turbulence. Since the pitch-angle depends on stochastic magnetic fields, we find pitch-angle scattering which needs to be described by statistical methods and transport equations. Furthermore, the process of pitch-angle scattering is eventually leading to a diffusive motion along the mean magnetic field as discussed in the following.

Transport Equations
The motion of energetic particles is described by transport equations. As explained in detail in the book Schlickeiser (2002), there is the following hierarchy of transport equations: 1. The relativistic Vlasov equation describing the particle motion in phase-space, 2. The Fokker-Planck equation which is an equation for the ensemble averaged particle distribution function, 3. The diffusive transport equation which is obtained from the Fokker-Planck equation by averaging over pitch-angle.
The general Fokker-Planck equation of cosmic ray transport is complicated and contains several terms such as stochastic acceleration and perpendicular diffusion (see, e.g., Skilling 1975;Schlickeiser 2002;Zank 2014). For some applications one can neglect such terms and consider the two-dimensional version of the Fokker-Planck equation where we have used time t , the particle position along the mean magnetic field z, the pitchangle cosine μ, the particle speed v, and the pitch-angle Fokker-Planck coefficient D μμ . The right hand side of this equation corresponds to a pitch-angle scattering term which is important in the theory of parallel diffusion as discussed below. Equation (4.11) has the form specified by Sturm-Liouville theory and the pitch-angle Fokker-Planck coefficient can be computed via A justification of this form will be provided in Sect. 5.4. The time-derivativeμ is, in principle, given by Eq. (4.10). However,μ depends on stochastic fields and stochastic particle properties such as velocities. Therefore, the calculation of D μμ via Eqs. (4.10) and (4.12) is not straightforward and requires special techniques such as perturbation theory. The solution of Eq. (4.11) provides the particle distribution function f = f (μ, z, t). Furthermore, it describes a pitch-angle isotropization process meaning that the solution f (μ, z, t) becomes independent of μ for t → ∞ as shown, for instance, in the numerical work of Lasuik and Shalchi (2017). For certain applications, one could be interested in the pitch-angle averaged particle distribution indicating the probability to find the particle at position z at time t . In the following we derive a differential equation for M(z, t). First we average Eq. (4.11) over all μ to find the exact relation where we have used for the pitch-angle Fokker-Planck coefficient. 12 To continue we define the current density allowing us to derive from Eq. (4.14) corresponding to a one-dimensional continuity equation. Now we integrate the Fokker-Planck equation (4.11) from −1 to μ to find Dividing this by D μμ and using this in the current density (4.16) yields Now we consider the late-time limit t → ∞. Since we expect a pitch-angle isotropization process, we can assume that (4.20) In this limit the current density, given by Eq. (4.19), becomes where in the last step we have employed the continuity equation (4.17). In the late-time limit the current density J tends to zero faster than the distribution M. Therefore, we can neglect the first term on the right hand side of Eq. (4.21) yielding We now define the parallel spatial diffusion coefficient via corresponding to Fick's second law. This is a usual diffusion equation or heat transport equation. Obviously, due to pitch-angle scattering, the particle motion becomes a diffusive motion along the mean magnetic field. Equation (4.23) provides a relation between the parallel spatial diffusion coefficient and the pitch-angle Fokker-Planck coefficient. This formula is one of the most fundamental and important relations in the transport theory of energetic particles. It allows for the calculation of the parallel diffusion coefficient κ if the pitch-angle Fokker-Planck coefficient D μμ is known. This relation was derived the first time in Jokipii (1966) and later, in a more systematic way, in Earl (1974).

The Pitch-Angle Fokker-Planck Coefficient
The calculation of the pitch-angle Fokker-Planck coefficient is problematic. One way of computing the parameter D μμ is the application of QLT (see again Jokipii 1966 or Schlickeiser 2002 for a review) similar compared to the QLT of random walking magnetic field lines. For magnetostatic slab turbulence, for instance, the quasi-linear pitch-angle Fokker-Planck coefficient is given by (see, e.g., Shalchi 2009a) This equation can be obtained by combining Eqs. (4.10) and (4.12). Furthermore, the velocities v x and v y as well as the magnetic fields were evaluated along the unperturbed orbit given by Eqs. (4.6) and (4.7), respectively. The Dirac delta in Eq. (4.26) describes the interaction between the turbulence and particles. Very clearly we can identify the so-called gyro-resonance condition meaning that only if the particles hit the resonance one finds scattering. As a consequence, high energy particles interact with the larger turbulence scales whereas low energy particles interact with smaller scales. This means, for instance, that the small scales of the dissipation range, which are neglected in spectrum (2.14), are essential for the scattering of low energy particles such as solar energetic particles. Of course, the form given by Eq. (4.26) is only valid for the very special case of magnetostatic slab turbulence. More complete forms for the pitch-angle Fokker-Planck coefficient were derived in the past. In Bieber et al. (1994) and Teufel and Schlickeiser (2003), for instance, the influence of dynamical turbulence and dissipation effects was explored. The book of Schlickeiser (2002) contains detailed derivations of the quasi-linear pitch-angle Fokker-Planck coefficient incorporating plasma wave propagation effects.
For the spectrum given by Eq. (2.14) the Fokker-Planck coefficient (4.26) becomes where we have used the relation (see, e.g., Zwillinger 2012) for Dirac's delta. To continue it is convenient to define the parameter sometimes called the dimensionless rigidity. Using this in Eq. (4.28) allows us to write corresponding to the quasi-linear pitch-angle Fokker-Planck coefficient for magnetostatic slab turbulence. However, like for field lines, the quasi-linear approximation is problematic and does not work in the general case. In particular QLT does not work for pitch-angles close to 90 • corresponding to μ = 0. The latter problem is known as the 90 • problem of scattering theory and it was attempted to be solved by incorporating non-linear effects (see, e.g., Völk 1973Völk , 1975Jones et al. 1973Jones et al. , 1978Owens 1974;Goldstein 1976). A more systematic approach is provided by the so-called second-order quasi-linear theory (SOQLT) developed in Shalchi (2005a). In this theory orbit fluctuations are estimated by employing the quasilinear approximation itself. The corrected orbits are then used in the formula for D μμ leading to resonance broadening. The second-order theory agrees well with simulations for cases where QLT fails completely such as parallel diffusion in isotropic turbulence (see, e.g., Tautz et al. 2008). Furthermore, perpendicular diffusion can lead to resonance broadening and can therefore influence the pitch-angle Fokker-Planck coefficient (see Shalchi et al. 2004c). In Qin and Shalchi (2009) the pitch-angle Fokker-Planck coefficient was extracted from test-particle simulations. This approach is difficult because in analytical theories one needs to know the parameter D μμ as a function of μ corresponding to the pitch-angle cosine. However, due to pitch-angle scattering itself, the parameter μ is changing if time passes. Therefore, one needs to find a parameter regime where pitch-angle scattering is so weak that μ is almost constant over an extended period of time. The little changes of μ can then be used in order to compute D μμ numerically for that given μ. In the considered parameter regime it was shown that pitch-angle scattering around μ ≈ 0 is indeed strong and that the aforementioned second-order theory works well. In the current review we abstain from a detailed discussion of non-linear pitch-angle scattering and parallel diffusion. A detailed review of this matter can be found in Shalchi (2009a). It has to be emphasized, however, that parallel transport is far from being straightforward and the quasi-linear approximation as well as the simple gyro-resonance picture fail in the general case.

A Simple Formula for the Parallel Diffusion Coefficient
As noted in the previous subsection, the calculation of the pitch-angle Fokker-Planck coefficient itself is difficult due to the importance of non-linear and dynamical turbulence effects. To obtain a simple formula for the parallel mean free path we perform the following steps: 1. We employ quasi-linear theory but keep in mind that the theory is not valid in the general case. 2. We employ the slab model and keep in mind that two-dimensional modes do exist and can contribute to the parallel diffusion coefficient via non-linear interactions. Alternatively, the slab/2D model can be replaced by a full three-dimensional turbulence model. In the latter case the 90 • problem becomes significant. 3. We consider magnetostatic turbulence. This approximation should be valid if we consider particles moving much faster than the Alfvén speed. 4. We employ Eq. (4.23) to obtain the parallel spatial diffusion coefficient.
Using the first three assumptions enables us to employ Eq. (4.31) for the pitch-angle Fokker-Planck coefficient. Since s > 1 we can clearly see that D μμ (μ = 0) = 0 which is, of course, an unrealistic result. It would even contradict the original assumption of a pitchangle isotropization process which requires that the particles can be scattered through 90 • pitch-angles. Fortunately, this behavior does not cause a singularity if Eq. (4.31) is used in Eq. (4.23) and even leads to a useful result for the parallel diffusion coefficient. Furthermore, it is slightly more convenient to compute the parallel mean free path defined via λ = 3κ /v. An explanation for this relation is given in Sect. 5.4 of this review. Using Eq. (4.31) in Eq. (4.23), and solving the corresponding μ-integral, yields (see, e.g., Gradshteyn and Ryzhik 2000 for the solution of the occurring integrals) where we have used the Gaussian hypergeometric function 2 F 1 (a, b; c; z). Formula (4.32) has the following asymptotic properties Clearly we can see how the parallel mean free path increases with rigidity. A simple formula which has the same asymptotic properties was given by Shalchi (2009a), namely which is visualized in Fig. 10 together with Eq. (4.32) as well as the asymptotic limits (4.33). Equation (4.34) can be used for simple estimations of the parallel mean free path. It needs to be emphasized, however, that Eq. (4.32) and its asymptotic forms, are a very special case. Fig. 10 The quasi-linear parallel mean free path for magnetostatic slab turbulence. Shown is the exact formula as given by Eq. (4.32) which is represented by the solid line. Also shown are the asymptotic limits as given by Eq. (4.33) represented by the dotted lines, as well as approximation (4.34) represented by the dashed line As pointed out before in this review, the slab model with the spectrum used here as well as the magnetostatic approximation does not work in the general case (see again Bieber et al. 1994;Schlickeiser 2002;Shalchi et al. 2006). Furthermore, quasi-linear theory fails in a lot of cases such as three-dimensional turbulence (see, e.g., Shalchi 2009a; Tautz et al. 2008).

A Subspace Approximation to the Solution of the Fokker-Planck Equation
It will be shown later in this review that one needs to know certain moments of the Fokker-Planck equation in analytical theories for perpendicular transport. Therefore, this matter will be discussed in the following. Equation (4.11) is difficult to solve analytically for the general case. Therefore, often an isotropic pitch-angle Fokker-Planck coefficient is used (4.35) A detailed discussion of the analytical form of D μμ and the validity of the isotropic regime can be found in  based on the aforementioned second-order quasi-linear theory. It was shown in the latter article that Eq. (4.31) is valid for |μ| δB/B 0 whereas the isotropic form is accurate for |μ| δB/B 0 . In Eq. (4.35) the parameter D does not depend on μ but it is a complicated function of turbulence and particle properties (see again ). For the isotropic form the parallel spatial diffusion coefficient can very easily be computed by using Eq. (4.23). Then this diffusion coefficient and the corresponding parallel mean free path are given by where we have used again λ = 3κ /v. In order to solve Eq. (4.11), we employ the Fourier transform so that the Fokker-Planck equation (4.11) becomes To continue we expand the solution of the latter equation in a series of Legendre polynomials where the coefficients C n are functions of time. Using Eq. (4.39) in differential equation (4.38) yields whereĊ n denotes the time-derivative of the coefficient C n . Two well-known relations for Legendre polynomials are given by (see, e.g., Abramowitz and Stegun 1974) as well as Therewith, Eq. (4.40) can be written as We now multiply this equation by the Legendre polynomial P m , integrate over μ, and use the orthogonality relation of Legendre polynomials (see again Abramowitz and Stegun 1974) Equation (4.45) corresponds to an infinite set of coupled differential equations. As an example we consider the case m = 0 and we finḋ For m = 1, on the other hand, we havė In the following we employ a two-dimensional subspace approximation. This means we set so that only the coefficients C 0 and C 1 are used. In Lasuik and Shalchi (2019) one can find the solution obtained by employing one-and three-dimensional subspace approxima-tions. It was shown there that the one-dimensional solution is not useful whereas the threedimensional solution is too complicated for most applications. Therefore, the two-dimensional solution provides a good compromise in terms of accuracy as well as analytical tractability. Within the two-dimensional subspace approximation expansion (4.39) reduces to In this case Eqs. (4.46) and (4.47) can be combined to find the second-order differential equationC Using the ansatz leads to the quadratic equation Thus, the two eigenvalues ω ± can be determined to be The coefficient C 0 can now be written as the linear combination and from Eq. (4.46) we derive (4.55) The solution F (μ, t) can, therefore, be expressed as In order to find the coefficients b ± , we can employ the initial condition This form simply means that the particle has its initial position at z = 0 and the initial pitch-angle cosine μ 0 . Using this in the inverse Fourier transform of Eq. (4.37) yields Together with the expansion (4.39) we obtain Multiplying this by P m and integrating over μ leads to where we have used again orthogonality relation (4.44). For m = 0 and m = 1 this gives as well as Using the latter two formulas in Eqs. (4.54) and (4.55) yields the following system of equations (4.64) Using this and Eq. (4.53) in Eq. (4.56) provides the two-dimensional subspace approximation to the solution F (μ, t). Of course, our solution is based on expansion (4.49). One can easily demonstrate by using Eqs. (4.13) and (4.16) together with the orthogonality relation (4.44) that the function C 0 (t) corresponds to the Fourier transform of the pitch-angle averaged distribution function M(z, t) and C 1 (t) corresponds to the Fourier transform of the current density or diffusion flux J (z, t).
In the theory of perpendicular diffusion one needs to know the quantity μ 0 μe −ik z as shown in Sect. 7 of this review. First we need to define the meaning of the ensemble average in particle transport theory. Since the particle distribution function derived above is a function of μ 0 , μ, and z, we understand the ensemble average used here as (4.65) Therefore, we need to evaluate the following quantity (4.66) Replacing f (z, μ, t) by Eq. (4.37) therein leads to With the help of Eq. (2.49) this becomes (4.68) We now replace F (k, μ, t) by Eq. (4.39), and use P 1 (μ) = μ to show that which, due to the orthogonality of Legendre polynomials (4.44), reduces to To solve the remaining integral we employ Eq. (4.55) with Eq. (4.64) to derive 13 (4.71) The latter formula contains the two parameters ω ± given by Eq. (4.53). We shall see later how Eq. (4.71) is used in analytical theories for perpendicular diffusion. In such theories, depending on how they are derived, one needs to know the complex conjugate of Eq. (4.71).
If ω + and ω − are real, taking the complex conjugate does not change the result. If they are complex, on the other hand, we have ω + = ω * − and ω − = ω * + . Therefore, we can easily see that This result was originally derived in  by using a slightly different approach.
The same calculations performed above can easily be repeated to compute (see Lasuik and Shalchi 2019 for details) corresponding to the characteristic function. This result can be further simplified by considering limits. We find (4.74) The first line therein corresponds to the characteristic function of a diffusion equation. However, this form is only valid for smaller wave numbers. The second line, on the other hand, corresponds to a damped unperturbed orbit. The conditions found here can be rewritten in terms of the parallel mean free path. We find the diffusive result for λ 2 k 2 3/4 and the damped unperturbed result for λ 2 k 2 3/4. We conclude that the shorter the parallel mean free path is, the more important the diffusive solution is.
The pitch-angle averaged function (see, e.g., Eq. (4.13)) can be obtained via As an example we consider the limit λ → 0 so that we can employ the first line of Eq. (4.74). We can easily derive corresponding to a Gaussian solution. The result obtained here is the diffusive solution that one would expect in this case. However, it has to be emphasized that it is only correct in the limit λ → 0. We can also compute the velocity correlation function based on this method. We find after some algebra It has to be emphasized that this form was derived from the Fokker-Planck equation but an isotropic pitch-angle Fokker-Planck coefficient was used. For other forms different velocity correlation functions can be obtained (see Shalchi 2011a).
Also useful is to determine the second moment z 2 because it controls the parallel spatial diffusion coefficient. After the same calculations performed above for the other quantities, we find (see Lasuik and Shalchi 2019 for a detailed derivation and Shalchi 2006a for an alternative derivation of this result) We can easily see that asymptotically This means that for early times parallel transport is ballistic (unperturbed) whereas for late times diffusion is restored. Therefore, Eq. (4.78) contains the important features of parallel transport. Furthermore, we can easily learn that the turnover from the ballistic to the diffusive regime occurs around the characteristic time vt ≈ λ meaning that the particles travel about a parallel mean free path before they become diffusive. This fact will be useful for a heuristic discussion of perpendicular transport as we shall see later.

Modeling Parallel Diffusion
In test-particle simulations one usually obtains parallel and perpendicular diffusion coefficients as a function of time (see Sect. 10 of this review). In the current review we consider theories in which the parallel diffusion coefficient occurs as input parameter. In the following we discuss a simple model for time-dependent parallel transport which has only the parallel mean free path and magnetic rigidity as input parameters. In the previous subsection the second moment was computed and Eq. (4.78) was obtained. In the following we define the running parallel diffusion coefficient via Usually a time-dependent or running diffusion coefficient is defined via the time-derivative of the corresponding mean square displacement. However, our aim is to compare the parallel diffusion coefficient calculated analytically with test-particle simulations. In principle it is not difficult to compute a running diffusion coefficient defined via time-derivatives from such numerical data. However, as well known, calculating derivatives of noisy data leads to an even stronger noise and this would make the aforementioned comparison difficult. Thus, we define the running parallel diffusion coefficient via Eq. (4.80) and keep in mind that for diffusive parallel transport, both definitions become equal anyway. By using the parallel spatial diffusion coefficient in the late-time limit and combining this with Eqs. (4.78) and (4.80), yields In Sect. 10 the form (4.82) will be compared with test-particle simulations. In order to do this, it is convenient to define the following dimensionless quantities: where is a characteristic length scale of turbulence (see Sect. 2 for more details), R is the dimensionless rigidity used before in Eq. (4.30), and is the gyrofrequency defined via Eq. (4.3). With this set of parameters, Eq. (4.82) can be written as If the dimensionless rigidity R and parallel diffusion coefficient K are specified, Eq. (4.84) can easily be used in order to plot D versus time T . Equation (4.84) was originally derived in Arendt and Shalchi (2018) and is compared with simulations in Figs. 16 and 17 of this review. Alternatively, one can express the parallel position in terms of the parallel bendover scale and the dimensionless time τ = vt/ so that Eq. (4.78) becomes The latter form has the advantage that only one parameter controls the solution, namely λ / .

Perpendicular Transport of Energetic Particles
In the previous sections we reviewed turbulence models, analytical theories for the random walk of magnetic field lines, and discussed parallel transport. These are three important ingredients in the theory of particle transport across the mean magnetic field. In the following we discuss quasi-linear transport, compound sub-diffusion, and the non-linear guiding center theory.

The FLRW Limit
One of the first descriptions of perpendicular transport was entirely based on the effect of field line random walk (see Jokipii 1966). This approach is simple but does not work in the general case as shown later. As an approximation we assume that FLRW is diffusive for all distances meaning that we neglect the initial ballistic regime entirely. Thus, we can employ Eq. (3.45). As a simple model for the perpendicular transport of energetic particles we assume that the particle follows a single magnetic field line while it moves with constant velocity in the parallel direction meaning that where we have used the pitch-angle cosine μ as before. Furthermore, if particles follow field lines, the perpendicular mean square displacement of particle trajectories and field lines are the same. Therefore, we can combine Eqs. (3.45) and (5.1) to find This corresponds to a diffusive motion of the particle where the diffusion coefficient is given by Compared to the parameter κ ⊥ used so far in this review, the Fokker-Planck coefficient D ⊥ (μ) allows for a pitch-angle dependent description of perpendicular transport. The two parameters are related to each other via (see, e.g., Schlickeiser 2002 for a detailed derivation) corresponding to a simple pitch-angle average. Solving the latter integral for the form given by Eq. (5.3) yields for the perpendicular diffusion coefficient We can clearly see that the perpendicular diffusion coefficient of the energetic particles is directly proportional to the diffusion coefficient of magnetic field lines. The only property of the particles entering this simple equation is their speed v. Alternatively, one can compute the perpendicular mean free path via Obviously the latter quantity does not depend on any particle property. This means that the perpendicular mean free path does even not depend on particle energy, momentum, or rigidity. Later in this review paper we will argue that this behavior can typically be found for high particle energies. For slab turbulence the field line diffusion coefficient is given by Eq. (3.21) or even Eq. (3.20) if spectrum (2.14) is used. Thus we have allowing for a straightforward calculation of the perpendicular mean free path. The result obtained here was obtained by combining the simple model that particles follow diffusive magnetic field lines. One can develop and apply a systematic quasi-linear approach as it was originally done in Jokipii (1966). For slab turbulence this would lead to exactly the same result, namely Eq. (5.7). The original work of Jokipii (1966), however, also allows the case δB z = 0 leading to another term in the quasi-linear perpendicular diffusion coefficient. This term is basically a gyro-resonant term which is negligible in most cases.

Compound Sub-Diffusion
Real particles do not perform an unperturbed motion in the parallel direction as assumed in the previous paragraph. Therefore, we now assume a diffusive parallel motion. However, we still assume that the particle is tied to a single magnetic field line as before. First we start with some simple considerations. If a typical particle moves diffusively, we can use instead of Eq. (5.1). Here we have used the parallel diffusion coefficient of the particle κ . Combining Eqs. (3.45) and (5.8) yields corresponding to a sub-diffusive motion since the mean square displacement is directly proportional to √ t . This type of transport is called compound sub-diffusion and is relevant in turbulence without transverse complexity (see, e.g., Mace et al. 2000;Qin et al. 2002a). According to Eq. (5.9) parallel diffusion suppresses perpendicular transport to a sub-diffusive level.
A more complete and accurate description of this type of transport was presented in Webb et al. (2006). The latter authors employed an approach based on the so-called Chapman-Kolmogorov equation (see, e.g., Gardiner 1985) where the particle distribution in the perpendicular direction f ⊥ (x, y; t) is given as convolution integral of the parallel distribution function f (z; t) and the field line distribution function f FL (x, y; z). Whereas Webb et al. (2006) computed the function f ⊥ (x, y; t) by solving Eq. (5.10), we follow Shalchi and Kourakis (2007a) and only consider the second moment of f ⊥ (x, y; t) leading to giving us a relation between particle and field line mean square displacements. In the following we employ a Gaussian particle distribution and assume that parallel transport is diffusive for all times. Thus, the parallel distribution function has the form Using this and Eq. (3.45) in Eq. (5.11) leads to This integral can easily be solved yielding similar compared to Eq. (5.9). The running diffusion coefficient is then corresponding to sub-diffusion. However, the Chapman-Kolmogorov equation based approach allows us to compute not just second moments but also the perpendicular particle distribution function f ⊥ (x, y; t). This was done in Webb et al. (2006) where a result depending on a so-called Fox function was found. Therefore, for compound sub-diffusion, one expects to find a non-Gaussian distribution function. This form, and therewith the result of Webb et al. (2006), was confirmed via test-particle simulations by Arendt and Shalchi (2018) and can also be seen in Fig. 16 of this review. A velocity correlation function based approach to describe compound sub-diffusion was presented in Kóta and Jokipii (2000).
Compound sub-diffusion is highly relevant for particle transport in slab turbulence. According to the so-called theorem on reduced dimensionality, proposed in Jokipii et al. (1993) and Jones et al. (1998), particles are tied to magnetic field lines if the turbulence has reduced dimensionality. Undoubtedly this is the case for slab turbulence. Numerically it was confirmed by Qin and Shalchi (2015) that for slab turbulence particles follow indeed field lines. However, if one includes wave propagation effects, for instance, perpendicular transport can be diffusive even if particles follow field lines (see Shalchi et al. 2007). It remains unclear whether the theorem on reduced dimensionality is valid for pure two-dimensional turbulence and how perpendicular transport behaves in this specific case (see again Qin and Shalchi 2015) but for slab/2D and full three-dimensional turbulence it is usually assumed that perpendicular transport is diffusive except some transient regime where one can still find sub-diffusion because normal diffusion is not yet restored. In Sect. 9 we shall provide a more complete description of the physics of perpendicular transport. Before doing this we focus on systematic theories.

Fundamental Equations for Perpendicular Diffusion
So far we discussed two simple models for perpendicular transport. Although perpendicular transport is by far more complicated, both types of transport can be seen in test-particle simulations as discussed later in this review. The FLRW limit should be valid in the limit of high particle energies corresponding to very long parallel mean free paths. Compound subdiffusion, on the other hand, should describe particle transport in slab turbulence correctly. Furthermore, there are intermediate times where we find sub-diffusion (see Fig. 14 of this review or Zimbardo et al. 2006). In the next few sections we discuss more complete theories for perpendicular transport which contain the aforementioned limits. Those considerations will lead to the so-called unified non-linear transport (UNLT) theory and its time-dependent generalization (see Sects. 6 and 7).
The fundamental equation describing the motion of charged particles through purely magnetic turbulence is the Newton-Lorentz equation (4.1). We now replace the particle position x therein by the guiding center coordinates (see, e.g., Schlickeiser 2002) where the gyrofrequency is given by Eq. (4.3). The Cartesian components of Eq. (5.16) have the form (5.17) In the unperturbed case, the coordinates (X, Y, Z) correspond to the coordinates of the guiding center. In reality those coordinates have a different meaning due to the chaotic motion of the particles. Using those coordinates, however, is still useful as it allows for a strong simplification of fundamental equations. Even in the case of particle motion through turbulence, we call those coordinates guiding center coordinates but we keep in mind that if there is turbulence the particle motion is no longer a perfectly helical motion. The guiding center velocity can be obtained by considering the time-derivative of Eq. (5.16). We can easily derive where we have employed the Newton-Lorentz equation (4.1) also. Using the standard formula for the double cross product (sometimes called vector triple product or Grassmann identity) allows us to write Considering Cartesian components yields (5.21) For the incompressible case δB z = 0 we find the following equations of motion Even if real turbulence is compressible, these equations should still provide a good approximation assuming that δB z < B 0 . In the current review we assume that the diffusion coefficients based on particle and guiding center coordinates are the same. This assumption is motivated by the fact that the distance between the coordinates (X, Y, Z) and the coordinates (x, y, z) is limited meaning that it does not increase with time. Furthermore, diffusion coefficients were calculated for the different coordinates based on test-particle simulations in Qin and Shalchi (2016) showing that there is no difference between the diffusion coefficients regardless whether guiding center or particle coordinates are used (see Sect. 10.5 of the current review for more details). Therefore, our investigations of perpendicular diffusion are based on Eq. (5.22) as equations of motion. Due to the assumption of axi-symmetry we only consider the x-component. Using again a Fourier representation as given by Eq. (2.3), allows us to write this equation of motion as ( 5.24) The latter equation will be the starting point of the so-called non-linear guiding center theory as well as the unified non-linear transport theory. Both approaches will be discussed in the following. Before doing so, we discuss another important tool frequently used in transport theory, namely the so-called Taylor-Green-Kubo formulation.

The Taylor-Green-Kubo Formula
In the section about magnetic field line random walk, we computed diffusion coefficients via mean square displacements (see, e.g., Eq. (3.1) of the current review). A powerful tool in transport theory is the so-called TGK (Taylor-Green-Kubo) formula (see Taylor 1922;Green 1951;Kubo 1957) which relates diffusion coefficients to so-called velocity auto-correlation functions. In the following this relation is derived and discussed. The difference between final and initial position of the particle can be written as where the coordinate x can be replaced by any other Cartesian coordinate. Therefore, the mean square displacement of possible particle trajectories can be written as where we have used the velocity auto-correlation function v x (t 1 )v x (t 2 ) . A time-dependent or running diffusion coefficient can be defined via 14 which is more common than definition (4.80). In order to find the time-derivative of the mean square displacement given by Eq. (5.26), we can use the Leibniz integral rule (3.29) to derive (5.28) Therewith the running diffusion coefficient (5.27) becomes The latter equation is exact and its derivation did not require any assumptions or approximations. Alternatively, we can use again Eq. (5.25) to write the running diffusion coefficient as To continue we assume that the transport is stationary, meaning that the velocity correlation function in Eq. (5.29) depends only on the time difference t − τ . Then we can write After employing a simple integral transformation we obtain In general the latter function can depend on time. This is the case for anomalous transport such as compound sub-diffusion. If we compute the time-derivative of Eq. (5.32), on the other hand, we find which is another useful relation as we shall see later in this review. If we consider the latetime limit t → ∞ and assume that the transport eventually becomes diffusive, Eq. (5.32) turns into which is the standard version of the TGK formula. This formula can also be used in order to understand the relation between the particle's mean free path and its diffusion coefficient. For simplicity we assume that the velocity correlation function is given by the exponential form where we have assumed isotropic initial conditions v 2 x (0) = v 2 /3. Using this velocity correlation function in the TGK formula (5.34) gives us and, therefore, This relation is frequently used throughout this review. We understand the mean free path, as defined via Eq. (5.35), as a characteristic length scale for the decorrelation of velocities. This means that after the particle has traveled the distance vt = λ ⊥ , the velocity correlation function basically decayed to zero. All considerations presented above were done for the x-components of position and velocity. All these relations can be used for any direction of space. However, special care is required if it comes to the drift coefficient (see, e.g., Shalchi 2011b). Then, on the other hand, the TGK formula (5.34) can be used in order to describe pitch-angle scattering (see, e.g., Eq. (4.12)) or even momentum diffusion (see, e.g., Schlickeiser 2002).

The Non-Linear Guiding Center Theory
The first systematic theory for perpendicular diffusion which is discussed here is the socalled non-linear guiding center (NLGC) theory. After some simple model-based approximations such as the non-linear closure approximation of Owens (1974) and the Matthaeus (1997) model, Matthaeus et al. (2003) developed a more systematic approach in order to describe perpendicular diffusion. The equations of motion (5.22) can be combined with the TGK formula (5.34) to obtain Besides quantities used before in this review, we have introduced the parameter a 2 . In the original paper by , the parameter a was introduced directly in Eq. (5.22) but it was shown in Qin and Shalchi (2016) that Eq. (5.22) is correct as it is without containing any type of additional parameter. In the following we just keep this parameter in our equations and assume that this parameter is there in order to balance out some inaccuracies of the theory related to some of the used approximations. To proceed,  assumed that the fourth-order correlation in Eq. (5.38) can be replaced by a product of two second-order correlations Approximation (5.39) is problematic and does not work in the general case as shown in the numerical work of Qin and Shalchi (2016). A more detailed discussion of this matter can be found in Sect. 10.5 of this review. It will be shown there that the parameter a 2 occurring in Eq. (5.40) is there to balance out the inaccuracy of approximation (5.39) in slab/2D turbulence. Furthermore, it was already shown in Shalchi (2005b) that Eq. (5.39) is the reason why NLGC theory provides a diffusive result for pure slab turbulence where one expects to find sub-diffusive transport. Therefore, approximation (5.39) is dropped in the derivation of UNLT theory.
Furthermore,  modeled the parallel velocity correlation function by the isotropic exponential model where we used the parallel mean free path λ . This corresponds to the velocity correlation function derived above (see Eq. (4.77) of this review). It was shown in Shalchi (2011a) that this form is only exact for an isotropic pitch-angle Fokker-Planck coefficient (see Eq. (4.35) of the current article). However, for other forms of D μμ , Eq. (5.41) should still provide an accurate approximation. To proceed, the magnetic correlation function δB x (t)δB * x (0) has to be replaced by the Fourier representation (2.3) leading to where we used x(t = 0) = 0. Now we employ Corrsin's independence hypothesis as in the theory of FLRW (see, e.g., Eq. (3.26) of this review) and use Eq. (2.6), corresponding to the assumption of homogeneous turbulence, to derive Next, the characteristic function e ik·x must be approximated. Matthaeus et al. assumed a Gaussian distribution of the particles and that the particle motion is diffusive for all times. Therefore, one can use corresponding to the characteristic function of a usual diffusion equation for the axisymmetric case. 15 It has to be noted that here it is assumed that perpendicular and parallel transport are diffusive for all times. For pure slab turbulence, however, we can set k ⊥ = 0 and the assumption of diffusive perpendicular transport is no longer part of the theory in this special case. Combining Eq. (5.43) with Eq. (5.44) and evaluating the time-integral, yields This is the non-linear integral equation of the NLGC theory as originally derived by . This integral equation can be used for turbulence models with purely magnetic fluctuations where δB z < 0. Furthermore, we assumed static turbulence here but the theory can easily be generalized to allow for dynamical turbulence (see Shalchi et al. 2004c). As shown in  and Bieber et al. (2004), Eq. (5.45) agrees with some simulations as well as solar wind observations as long as a two-component turbulence model 15 A diffusion equation for the axi-symmetric case has the form ∂f/∂t = κ ∂ 2 f/∂z 2 + κ ⊥ (∂ 2 f/∂x 2 + ∂ 2 f/∂y 2 ). In Fourier space this turns into ∂F /∂t = −κ k 2 F − κ ⊥ k 2 ⊥ F having the solution shown in Eq. (5.44). The solution of the diffusion equation in Fourier space corresponds to the characteristic function apart from factors depending on how the Fourier transform is defined. is used and one sets a 2 = 1/3. This was undoubtedly a breakthrough in the theory of perpendicular diffusion. However, Eq. (5.45) cannot be seen as the final solution to this problem. One of the major problems is that NLGC theory does not work for slab turbulence where perpendicular transport should be sub-diffusive. This problem was also pointed out in the review of Giacalone (2013). Since NLGC theory was presented, several authors developed alternative theories. Examples are the weakly non-linear theory of Shalchi et al. (2004b), the random ballistic interpretation of NLGC theory of , or the approaches presented by le Roux and Webb (2007) and Qin and Zhang (2014). These alternative theories definitely provided some interesting insight and contain useful ideas. However, they did not really lead to a theory which provides a correct description of perpendicular transport in slab turbulence, nor do they contain the Matthaeus et al. (1995) theory of random walking magnetic field lines. An exception is the extended non-linear guiding center theory presented in Shalchi (2006b) which describes correctly sub-diffusion for slab turbulence but the theory was specifically designed for two-component turbulence and cannot be used for three-dimensional models.
Another problem is that the correction factor a 2 is needed in order to achieve agreement with simulations. Since the theory does not contain an explanation of this parameter nor is it clear why a 2 = 1/3 is needed, NLGC theory is incomplete. Unfortunately, the correction factor a 2 will be kept in the next three sections although in several cases a 2 is no longer needed in more advanced theories. The physics of a 2 and its value will eventually be explained in Sect. 9.

The Unified Non-Linear Transport Theory
The NLGC theory discussed before is based on approximation (5.39). It was shown numerically and analytically that this approximation is not valid in the general case (see Shalchi 2005b;Qin and Shalchi 2016, and Sect. 10.5 of the current review). Therefore, UNLT theory was developed which is not based on this approximation. Before we discuss solutions of the latter theory, we derive UNLT theory by employing the approach based on the Fokker-Planck equation as presented in Shalchi (2010).

Derivation from the Fokker-Planck Equation
In order to determine the correlations between the velocity component v z (t) and the particle trajectory x(t) in the magnetic field, one needs a velocity-dependent transport equation. This equation is given by the pitch-angle dependent Fokker-Planck equation similar compared to Eq. (4.11) of this review. However, in order to derive diffusive UNLT theory, we need to take into account perpendicular diffusion so that the Fokker-Planck equation becomes The solution of this equation provides the pitch-angle dependent particle distribution function f (x, μ, t). Before using this, we go back to Eq. (5.38) and use the Fourier representation (2.3) to find As before we employ Corrsin's approximation to obtain where we included the correction factor a 2 as in NLGC theory. Ideally we would have a 2 = 1 but a more detailed discussion of this factor can be found in Sects. 9 and 10 of this review. Now we assume homogeneous turbulence meaning that we can employ Eq. (2.6). Thus, we can write the perpendicular diffusion coefficient as with the wave vector dependent function (6.5) Above we assumed magnetostatic turbulence but diffusive UNLT theory can also be formulated for dynamical turbulence (see discussion at the end of this subsection). The ensemble average used here is defined like in Eq. (4.65) but with additional x-and y-integrals. Thus, the correlations in Eq. (6.5) can be written as with the pitch-angle dependent characteristic function (k, μ, t) := d 3 xe ik·x f (x, μ, μ 0 , t). (6.7) In the notation used here we emphasize that the distribution function also depends on the initial pitch-angle cosine μ 0 . To continue, we multiply the Fokker-Planck equation (6.1) by exp (ik · x) and thereafter we integrate over the whole space to obtain corresponding to a partial differential equation for . To proceed it is useful to define the pitch-angle dependent function and with Eqs. (6.5) and (6.6) we obtain where we have used integration by parts. In order to derive an ordinary differential equation for the function S(μ, k), we multiply Eq. (6.8) by μ 0 , integrate over time t , and average over the initial pitch-angle cosine μ 0 to find Due to the assumption of sharp initial conditions we can use f (x, μ, t) = 2δ(x)δ(μ − μ 0 ) and, thus, we have Furthermore, due to the fact that we find a pitch-angle isotropization process (see, e.g., Sect. 4 of this review), the particle distribution function f and therewith the function become μ 0 -and μ-independent for late times. Thus, if we multiply the pitch-angle independent function (k, t → ∞) by μ 0 and integrate over μ 0 , we find zero. Therewith, Eq. (6.11) becomes If the function S would be known, we can easily compute the function T via Eq. (6.10) and thereafter the perpendicular diffusion coefficient via Eq. (6.4). Equation (6.13) provides an ordinary differential equation for the function S which depends on μ as well as the wave vector k. Unfortunately, it is not possible to solve Eq. (6.13) for the general case. To continue we assume that the two Fokker-Planck coefficients D μμ and D ⊥ are even functions in μ. Furthermore, we split S in an even contribution S + and an odd contribution S − . Then we can derive from Eq. (6.13) a system of two ordinary differential equations for S + and S − , namely which is even in μ and which is odd in μ. Equation (6.14) can be averaged over μ. By using Eqs. (4.15) and (6.5), we can easily derive which is an exact relation. The left hand side of this equation corresponds to the function T defined via Eq. (6.10) so that we can write Here we can already see that for slab turbulence, where we have k ⊥ = 0, we find T = 0 and, thus, κ ⊥ = 0 corresponding to sub-diffusion. To proceed we assume which is a simple model based on the assumption that D ⊥ (μ) is symmetric in μ and increases with |μ|. The assumption that D ⊥ ∼ |μ| is supported by the test-particle simulations performed by Qin and Shalchi (2014). Furthermore, this pitch-angle dependence is also provided by the FLRW limit (see, e.g., Eq. (5.3) of the current article). The factor 2 in Eq. (6.18) has been chosen so, that Eq. (5.4) is satisfied. With the form given by Eq. (6.18), the exact relation (6.17) turns into Now we integrate Eq. (6.15) over the pitch-angle cosine from 0 to 1 and employ Eq. (6.19) to obtain where we have used S − ≡ ∂S − /∂μ and D ≡ D μμ (μ = 0). Using Eqs. (6.18) and (6.19) yields where we have employed Eq. (6.10) again. Due to the factor (1 − μ 2 ) in Eq. (6.10), we approximate Therefore, Eq. (6.21) becomes .

(6.24)
This result is valid for arbitrary but finite D ≡ D μμ (μ = 0). A further simplification can be achieved by employing the isotropic pitch-angle Fokker-Planck coefficient as given by Eq. (4.35). This enables us to use Eq. (4.36) and Eq. (6.4) becomes where we have used the function This can easily be integrated for arbitrary D μμ so that (6.28) Using this in Eq. (6.10) yields where we have used Eq. (4.23). One can easily see that this is in agreement with Eq. (6.25) if we set k = 0 therein, suggesting that Eq. (6.25) is also valid if D μμ is not isotropic. One can repeat the calculations performed above for the more general case of dynamical turbulence. In this case one would obtain the formula derived in Shalchi (2011d). Diffusive UNLT with plasma wave propagation effects was derived in Hussein and Shalchi (2014b). Furthermore, the theory was extended in Wang and Qin (2018) by including the effect of adiabatic focusing, which occurs if the mean magnetic field is curved. We also note that Eq. (6.25) has some similarity with the NLGC theory represented by Eq. (5.45). However, the two integral equations are not identical but they can provide similar results in certain limits. This is in particular the case for pure two-dimensional turbulence. For slab and three-dimensional turbulence, on the other hand, the solutions of these two equations are very different.
Diffusive UNLT theory provides the non-linear integral equation for the perpendicular diffusion coefficient given by Eq. (6.25). In the following we discuss asymptotic limits which are contained in this theory. This will link Eq. (6.25) to previously derived and discussed theories such as QLT. Before we derive these limits systematically, we show the relevance of the Kubo number in the theory of perpendicular diffusion.

The Importance of the Kubo Number
Usually, in the theory of energetic particles, one needs to specify the spectral tensor describing the properties of magnetic turbulence. In the current subsection we employ a different approach. We assume that the components of the spectral tensor have the general form where we have used the dimensionless function f (x, y) which depends only on x = k and y = k ⊥ ⊥ . Furthermore, this function decays with increasing parallel and perpendicular wave numbers. Therefore, and ⊥ are characteristic length scales for the decorrelation of the turbulence namely the bendover scales as discussed in Sect. 2 of this review. The reason why Eq. (6.30) is valid in the general case is the following. The components of the spectral tensor need to satisfy the normalization constraint (2.16). Therefore, we have  Table 4. As demonstrated, all incompressible models discussed in this review have spectral tensors of the form given by Eq. (6.30). Furthermore, we define the diffusion ratio via and use the Kubo number K defined via Eq. (3.11). Therewith, and with Eq. (6.30), Eq. (6.25) can be written as (6.34) It follows from Eq. (6.33) that within the units used here, we measure all distances in the perpendicular direction with respect to ⊥ and all distances in the parallel direction with respect to . According to Eq. (6.34), there are only two parameters controlling the diffusion

Fig. 11
Diffusive UNLT theory contains four asymptotic limits. We obtain these limits for short or long parallel mean free paths and extreme values of the Kubo number defined in Eq. (3.11). For long parallel mean free paths and small Kubo numbers, for instance, we find the quasi-linear limit. For a Kubo number equal to zero, corresponding to slab turbulence, and a finite parallel mean free path we find a vanishing perpendicular diffusion coefficient usually interpreted as sub-diffusive transport ratio D, namely the parallel mean free path normalized with respect to the parallel scale λ / and the Kubo number K. If we focus on small and large values of those parameters, there are four different asymptotic limits and, therewith, four different transport regimes. All limits are compared with each other in Table 5 and they are visualized in Fig. 11. The validity of these limits and further regimes are discussed in Sects. 7 and 9 of this review. In the following we explore these four limits. This could either be done by starting with Eq. (6.25) or Eq. (6.34). In this review we use the first option but keep in mind Eq. (6.34) and the relevance of the Kubo number. This will help us to understand in which regime the considered limit is valid. Also important here is to note that in some cases one has equal bendover scales in the two directions of space = ⊥ . In this particular case the Kubo number is nothing else than the magnetic field ratio δB x /B 0 sometimes called the Alfvénic Mach number. Examples for turbulence models where this is the case are, of course, isotropic turbulence but also models based on Goldreich & Sridhar's critical balance.

The Quasi-Linear Regime
First we consider the limit λ / → ∞ corresponding to the case that pitch-angle scattering and therewith parallel diffusion are suppressed. From a more practical point of view this corresponds to the case of higher particle rigidities or energies. In the case considered here we can set v/λ = 0 in Eq. (6.25) to derive the simple equation (6.35) To evaluate this further, we use the notation Using this in Eq. (6.35) yields for the parameter κ FL the following equation where we have usedB 0 := B 0 /a. Equation (6.36) agrees with Eq. (5.5) corresponding to the FLRW limit. Equation (6.37), on the other hand, is in perfect agreement with Eq. (3.47) if we set a 2 = 1. Therefore, it is easy to understand the result obtained here. For λ / → ∞, corresponding to suppressed pitch-angle scattering, perpendicular diffusion is caused by FLRW and no other effect is controlling the transport. Furthermore, the Matthaeus et al. (1995) theory for field line wandering is contained in Eq. (6.25). Therefore, UNLT theory can be seen as a unified transport theory for magnetic field lines and energetic particles, hence the name unified non-linear transport theory.
The quasi-linear limit can be obtained if we additionally consider κ FL → 0 on the right hand side of Eq. (6.37). This corresponds to the limit of small Kubo numbers. Using relation (3.57) in Eq. (6.37) yields corresponding to the well-known quasi-linear result (see Eq. (3.25) of the current paper). As shown here, QLT for perpendicular diffusion can be obtained in the limit of long parallel mean free paths and small Kubo numbers. Using the general form (6.30) in Eq. (6.38) as well as Eqs. (6.36) and (5.37) allows us to derive the scaling law It is usually assumed that the parallel mean free path increases with increasing magnetic rigidity (see, e.g., Fig. 10 of the current review). Therefore, the quasi-linear scaling should become relevant if high energy particles propagating through small Kubo number turbulence are considered.

The Non-Linear FLRW Limit
Here we still consider the case of λ / → ∞ and, thus, Eqs. (6.36) and (6.37) are still valid. However, we now assume a large field line diffusion coefficient corresponding to large values of the Kubo number. Considering κ FL → ∞ on the right hand side of Eq. (6.37) yields in agreement with Eq. (3.58). The latter limit is sometimes called the non-linear regime or the Bohm limit of field line diffusion and was originally obtained by Kadomtsev and Pogutse (1979). In this case particles still follow magnetic field lines but the field lines are highly non-linear (see Sect. 3.9 of the current review). For the general spectral tensor (6.30) we find after some simple algebra the scaling law which is different compared to the quasi-linear scaling given by Eq. (6.39). This result should be valid for high particle energies and large Kubo number turbulence. It has to be noted that the scale in Eq. (6.41) is really the ultra-scale L U but in most cases we expect L U ∝ ⊥ as indicated by Table 5.

The Fluid Limit
We now consider the limit λ / → 0 corresponding to strong pitch-angle scattering. This case is important for smaller rigidities. If pitch-angle scattering is strong we expect that the perpendicular diffusion coefficient is somehow related to the parallel diffusion coefficient.
In the case considered here, Eq. (6.25) becomes Here we kept the function F (k , k ⊥ ) in the equation because for a small parallel diffusion coefficient, we also expect that the perpendicular diffusion coefficient is small. This means that we can neglect the term directly proportional to κ ⊥ k 2 ⊥ in Eq. (6.25) but not the term F (k , k ⊥ ). Replacing the parallel mean free path via λ = 3κ /v and using Eq. (6.26) allows us to write (6.43) In order to simplify this further, we consider two cases for the ratio κ ⊥ /κ . Since this ratio depends on the Kubo number, small and large values of κ ⊥ /κ correspond to small and large Kubo numbers, respectively. First we assume that this ratio is large and we find Important here is to note that this result does not depend on the details of turbulence. Only the ratio of the turbulent field with respect to the mean field enters the latter equation. The result obtained here was derived before by Owens (1974) as well as Zybin and Istomin (1985). Krommes et al. (1983) called Eq. (6.45) with a 2 = 1 the fluid limit. Characteristic here is that the ratio κ ⊥ /κ does not depend on any particle properties such as energy or rigidity. However, this ratio is typically in the order of unity in the considered case. The physical reason for the parameter a 2 and its value will be discussed in Sect. 9.

The Collisionless Rechester & Rosenbluth Scaling
Missing is the case relevant for short parallel mean free paths, corresponding to strong pitchangle scattering, and small Kubo numbers. Therefore, we still consider the case λ / → 0 and, thus, Eq. (6.43) is still valid. In order to simplify the latter equation in the limit considered here, it is useful to write Eq. (6.43) as Now we assume that the ratio κ ⊥ /κ is small. Using Eq. (3.57) yields The latter formula can easily be combined with any turbulence model as long as the occurring integrals are convergent. With the general form (6.30) we now find the scaling law which sensitively depends on the scale ratio / ⊥ as well as the magnetic field ratio δB x /B 0 . Characteristic here is that the ratio κ ⊥ /κ does not depend on magnetic rigidity. Furthermore, we expect to find small values of that ratio in this limit (see next subsection for some examples). As explained in detail in Sect. 9, this result can also be obtained by assuming a Rechester and Rosenbluth (1978) type of transport. Since Rechester & Rosenbluth considered collisions whereas we explore the collisionless case, we call the result obtained here the collisionless Rechester & Rosenbluth (CLRR) limit as originally suggested in Shalchi (2015a).

Illustrative Examples
The CLRR limit as given by Eqs. (6.47) or (6.48) provides a very sensitive dependence on turbulence parameters. Therefore, we focus now on the CLRR limit and consider two examples.
The first example is the noisy slab model defined via Eq. (2.25). This example is perfect because it corresponds to the small Kubo number limit. In this case Eq. (6.47) becomes To continue we employ the model spectrum given by Eq. (2.14) to obtain which is clearly a special case of the scaling law given by Eq. (6.48). To estimate a number for the ratio κ ⊥ /κ in the limit considered here, we use δB 2 /B 2 0 ≈ 0.5, / ⊥ ≈ 0.75, and a 2 ≈ 1. This choice for a 2 can be justified by comparing diffusive UNLT theory with simulations performed for noisy slab turbulence (see Sect. 10.6 of this review). For s = 5/3 we find C(s = 5/3) ≈ 0.12 and, thus, Eq. (6.50) provides κ ⊥ /κ ≈ 0.005 corresponding to a small and rigidity independent ratio.
As a second example we employ the Gaussian correlation model defined via Eq. (2.28). For simplicity we consider the special case q = 3. After evaluating Eq. (6.47) one finds for this case again in agreement with the scaling law given by Eq. (6.48). Typical for the CLRR regime is that the ratio κ ⊥ /κ does not depend on particle rigidity or energy. Furthermore, this ratio is much smaller than one for small Kubo number turbulence.

Scaling Laws in the Theory of Perpendicular Diffusion
Above we have derived four different scaling laws from diffusive UNLT theory and they are summarized in Table 5 as well as Fig. 11. One could try to find general scaling laws for perpendicular transport in magnetic turbulence. This corresponds to the considerations performed in Sect. 3.10 for magnetic field lines but is now done for the particles. This idea was used in Hauff et al. (2010) and is reviewed and then generalized in the current subsection. First we focus on the case of ballistic parallel transport so that v = vμ = const. (6.52) This corresponds to the case that pitch-angle scattering is suppressed. We now try to find a general form of the pitch-angle Fokker-Planck coefficient of perpendicular diffusion D ⊥ (μ). First we note that the diffusion coefficient D ⊥ has the dimensions [L] 2 /[T ]. Therefore, it is natural to assume that D ⊥ ∝ v|μ|. Furthermore, we measure all distances along the mean field with respect to the parallel scale and all distances across the mean field with respect to the perpendicular scale ⊥ . Therefore, one can assume that D ⊥ / 2 ⊥ ∝ v|μ|/ . The missing constant of proportionality does not have any dimensions but should depend on the properties of turbulence. We assume that this constant depends only on the Kubo number K and that this dependence is a simple power-law with the exponent γ . Therefore, it is appropriate to assume Replacing the Kubo number by employing Eq. (3.11) yields the following scaling law We can easily see that the quasi-linear scaling given by Eq. (6.39) can be obtained by setting γ = 2 and the non-linear FLRW limit, given by Eq. (6.41), can be obtained for γ = 1. The percolation regime (see discussion in Sects. 3.9 and 3.10), on the other hand, corresponds to γ = 0.7 so that as, for instance, discussed in Hauff et al. (2010). The quasi-linear regime can be found for small Kubo numbers K 1. The non-linear regime, on the other hand, arises from diffusive UNLT theory for the opposite case K 1. It is still unclear, whether the percolation regime exists in cases relevant in space and astrophysics. However, it is unlikely that it can be derived from UNLT theory because this theory is based on Corrsin's hypothesis whereas the percolation regime is beyond this approximation. The considerations presented above are based on the assumption that pitch-angle scattering is suppressed meaning that v|μ| is a constant. In general this is unrealistic because in reality particles experience strong pitch-angle scattering leading to a diffusive parallel motion. However, we can repeat the considerations presented above for this case as well. First we assume that the perpendicular diffusion coefficient is directly proportional to the parallel diffusion coefficient κ ⊥ ∝ κ because both parameters have the dimensions [L] Again we can recover two cases obtained before in this review. By setting either γ = 2 or γ = 4 we obtain the fluid limit given by Eq. (6.45) or the CLRR scaling given by Eq. (6.48). However, there could be exotic cases where one finds different scaling laws. In Shalchi (2015b), for instance, perpendicular transport for the NRMHD model defined via Eq. (2.27) was explored based on diffusive UNLT theory. The following scaling was obtained corresponding to the case γ = 3 in Eq. (6.58).
Interesting here is that more detailed formulas can be obtained if models with reduced dimensionality are considered. For slab turbulence, for instance, the perpendicular scale ⊥ must not occur. Therefore, Eq. (6.54) can only be valid for γ = 2 automatically leading to the quasi-linear scaling as the only solution. For the case of non-vanishing pitch-angle scattering, one has to use Eq. (6.58). This would lead to a diffusive result contradicting the sub-diffusive behavior. The reason for this is that the simple dimension-based considerations performed here only work if perpendicular transport is diffusive. However, we can generalize the considerations presented above, to cover anomalous transport as well. Based on dimensions, we assume that which can be written as Since we consider transport in pure slab turbulence, the scale ⊥ must not occur in the latter formula. This means that we need to set γ = 2 leading to Assuming that the perpendicular diffusion coefficient is directly proportional to the field line diffusion coefficient and using that κ FL ∝ allows us to set α = 1/2 leading to in agreement with compound sub-diffusion discussed before (see, e.g., Eq. (5.15) of this review). For two-dimensional turbulence we can perform the same considerations but in this case the corresponding scaling law must not depend on the parallel scale . Therefore, Eq. (6.54) is only valid for γ = 1 leading to These are exactly the limits one finds for two-dimensional and two-component turbulence as we shall see in Sect. 8 of this review.

Time-Dependent Perpendicular Transport
For certain applications, a time-dependent description of perpendicular transport could be desired. Therefore, we now focus on time-dependent UNLT theory where the diffusion approximation is no longer used. Ideally a time-dependent theory should contain sub-diffusive transport for slab turbulence but the theory should also explain why and when we get normal diffusion. The time-dependent UNLT theory was presented in Shalchi (2017) as well as Lasuik and Shalchi (2017). In the following we discuss this theory and show how diffusive UNLT theory can be derived from the time-dependent approach. Thereafter, we discuss special cases.

Time-Dependent UNLT Theory
Perpendicular transport is described by the auto-correlation function V x (t)V x (0) where the guiding center velocity is given by Eq. (5.24). It is useful to write the latter equation as where we have used the notation x · k = zk + x ⊥ · k ⊥ with the two-dimensional vectors x ⊥ = (x, y) and k ⊥ = (k x , k y ). Therefore, we find for the velocity correlation function where we have set x(t = 0) = 0. The central idea of time-dependent UNLT theory is to use the following approximation meaning that we group together magnetic fields, all particle properties associated with their parallel motion, and all particle properties associated with their perpendicular motion. This can be understood as an extension of the Corrsin approximation used before in this review. Clearly this is different compared to approximation (5.39) used during the derivation of NLGC theory. In UNLT theory we no longer assume that the parallel position and the parallel velocity are uncorrelated. However, as in NLGC theory, we introduce the parameter a 2 which balances out inaccuracies arising from the approximation used here. Using Eq. (7.3) in (7.2) and employing Eq. (2.6) allows us to write the velocity auto-correlation function as where we have used the parallel correlation function

ξ(k , t) = v z (t)v z (0)e iz(t)k
= v 2 μμ 0 e izk (7.5) corresponding to the moment defined and calculated in Eqs. (4.71) and (4.72). For the perpendicular characteristic function in Eq. (7.4) we set e ix ⊥ ·k ⊥ = e − 1 2 ( x) 2 k 2 ⊥ (7.6) corresponding to a Gaussian distribution with vanishing mean. In Lasuik and Shalchi (2018) other forms of distribution functions, such as kappa distributions, were employed but it was shown that the assumed statistics has only a minor influence on the resulting perpendicular diffusion parameter. In Eq. (7.6) we do not distinguish between particle and guiding center coordinates as before. However, finite gyroradius effects could be important for high particle energies as discussed in Sect. 7.5. We now have all ingredients necessary to formulate a time-dependent theory for perpendicular transport. After combining Eqs. (7.4) and (7.6), we can write the velocity correlation function as and after employing Eq. (5.33) we find The latter formula is a second-order differential equation for the mean square displacement ( x) 2 . What the solution of this equation is, depends on the spectral tensor P xx (k, t), describing magnetic turbulence, and the parallel correlation function ξ(k , t). It has to be emphasized that at no point have we assumed that perpendicular transport is diffusive. Within a two-dimensional subspace approximation, it was derived in Eqs. (4.71) and (4.72) of the current review that ξ(k , t) can be approximated by with the parameters ω ± given by Eq. (4.53). After using Eq. (4.36) this can be expressed by the parallel mean free path via Since the function ξ(k , t) is known, the ordinary differential equation (7.8) can be evaluated numerically for any given turbulence model described by the tensor component P xx .
We like to emphasize that Eq. (7.8) can easily be evaluated for dynamical turbulence as well. Equation (7.8) is very powerful because it allows for a time-dependent description of the transport and covers also cases where the transport never becomes diffusive such as compound sub-diffusion in slab turbulence.

Diffusive UNLT Theory
In the following we combine time-dependent UNLT theory with a diffusion approximation in order to explore how Eq. (7.8) is related to diffusive UNLT theory represented by Eq. (6.25). Using a diffusion approximation means that we set ( x) 2 = 2κ ⊥ t ∀t (7.11) on the right hand side of Eq. (7.7). This approximation was already used in the theory of FLRW and during the derivation of NLGC theory. In reality one expects a ballistic motion (see next sub-section) and thereafter there could be a sub-diffusive regime (see, e.g., Lasuik and Shalchi 2017). Eventually diffusion is recovered if there is transverse structure (see, e.g., Shalchi 2017; Lasuik and Shalchi 2017). A detailed discussion of these effects can be found in Sects. 7.3, 7.5, 7.6, 9, and 10 of this review. If a diffusion approximation is used, however, all these effects are ignored. Using Eqs. (7.9) and (7.11) in (7.7) yields According to the Taylor-Green-Kubo formula (5.34) we find for the perpendicular diffusion coefficient in the case of static turbulence (7.13) From Eq. (7.10) it follows that ω + + ω − = −v/λ (7.14) and Using the latter two relations in Eq. (7.13) yields where F (k , k ⊥ ) is still given by Eq. (6.26). This result agrees with Eq. (6.25) apart from a factor 4/3 in the second term of the denominator. The derivation presented here has the advantage that one can see it as a special case of the time-dependent theory. The only difference between Eqs. (7.16) and (6.25) is the factor 4/3. Although this factor is a minor difference, strictly speaking Eq. (6.25) should be more accurate than Eq. (7.16) because this factor arises from a pitch-angle averaging procedure. The pitch-angle dependence of ( x) 2 is neglected during the derivation of time-dependent UNLT theory. Therefore, we employ Eq. (6.25) if diffusive UNLT theory is used.

The Initial Free-Streaming Regime
Before we solve Eq. (7.8) for specific turbulence models we focus on early times. From Eq. (7.9) it follows for t → 0 that Here we have set a 2 = 1 because for early times, this correction factor is not needed. Also assuming ( x) 2 → 0 allows us to approximate Eq. (7.8) by where we have employed the normalization condition (2.16). Integrating over time twice yields corresponding to ballistic perpendicular transport at early times. The running diffusion coefficient is in this case corresponding to a linear increase with time.

Suppressed Pitch-Angle Scattering
Equation (7.8) can also be simplified significantly if we assume that pitch-angle scattering is suppressed. In this case the function ξ(k , t) defined in Eq. (7.5) becomes ξ(k , t) = v 2 μ 2 e ivμk t (7.21) corresponding to an unperturbed parallel motion. Using this in Eq. (7.8) yields Although we consider the special case of suppressed pitch-angle scattering, the latter result is interesting because it allows for a time-and pitch-angle dependent description of perpendicular transport. We can combine this result with a diffusion approximation similar compared to Eq. (7.11) ( X) 2 = 2D ⊥ t (7.23) but now the corresponding diffusion coefficient D ⊥ is pitch-angle dependent. Using this diffusion approximation in Eq. (7.22), assuming static turbulence, and integrating over all times, yields This can be understood as special case of diffusive UNLT theory because of the assumption of suppressed pitch-angle scattering. Then, on the other hand, it is also more general because it allows for a μ-dependent description of the transport. Equation (7.24) can easily be written as where the parameter κ FL is again the solution of Eq. (6.37). The solution (7.25) corresponds to the pitch-angle dependent version of the field line random walk limit (see, e.g., Eq. (5.3) of this review).

Finite Gyroradius Effects
So far we have neglected finite gyroradius effects completely. However, these effects can be important for high particle energies as discussed in Neuer and Spatschek (2006), Shalchi (2015cShalchi ( , 2016b, and Qin and Shalchi (2016). In all those previous articles it was found that at high energies there is a reduction of the perpendicular diffusion coefficient if finite gyroradius effects are taken into account.
In the current subsection we give a brief introduction of finite gyroradius effects and explain how they can be incorporated into UNLT theory. Since gyroradius effects are assumed to be important in the high rigidity limit, we concentrate on the case of suppressed pitchangle and gyrophase scattering as it was done in the previous subsection. In this case we can write Eq. (7.2) as The perpendicular position x ⊥ (t) therein is the particle position and not the guiding center position. In the previous sections we did not differentiate between particle and guiding center coordinates. However, this is done in the following. To continue we write the first two lines of Eq. (5.17) as where we have used the velocity components given by Eq. (4.6). To continue it is useful to write Eq. (7.26) as where the vector x g (t) describes the gyro-rotation. We now evaluate this term further. First we combine Eq. (7.27) with (2.12) to write where we have used 7.30) and the corresponding addition theorem for trigonometric functions. In order to rewrite the exponential function involving x g in Eq. (7.28), we can use Eq. (7.29) to derive where we have used the so-called Jacobi-Anger expansion (see, e.g., Cuyt et al. 2008). Therefore, we obtain To continue we average over the initial gyrophase 0 and employ (see, e.g., Zwillinger 2012) Using this in Eq. (7.28), assuming magnetostatic turbulence, using the diffusion approximation (7.23), and integrating over time yields with the resonance function (7.36) Please note that the latter resonance function corresponds to the one used in Shalchi et al. (2004b) if we set D μμ = 0 therein. Furthermore we can recover quasi-linear theory by considering the limit D ⊥ → 0 at the right hand side. Equation (7.35) with (7.36) describes perpendicular transport for suppressed velocity diffusion but with finite gyorradius effects. We can recover the zero gyroradius limit by considering W → 0. In this case we can use (see, e.g., Abramowitz and Stegun 1974) (7.37) and Eq. (7.35) becomes (7.39) This result agrees with the UNLT theory for the case of suppressed pitch-angle scattering (see Eq. (7.24) of this review). For arbitrary Larmor radii R L we can compute the perpendicular mean free path by solving Eq. (7.35) numerically. For two-dimensional turbulence and spectrum (2.20) the results are visualized in Fig. 12. We have computed λ ⊥ / ⊥ versus R L / for different values of the ratio ⊥ / . In all considered cases we find that finite gyroradius corrections reduce the Fig. 12 The perpendicular mean free path λ ⊥ versus the unperturbed gyroradius R L . The first parameter is normalized with respect to the perpendicular bendover scale ⊥ and the second parameter to the parallel bendover scale . All curves shown here are obtained for two-dimensional turbulence and q = 3. The solid line corresponds to the zero Larmor radius limit and the other curves are obtained by taking into account the non-vanishing gyroradius. Shown are results for ⊥ / = 10 (dotted line), ⊥ / = 1 (dashed line), and ⊥ / = 0.1 (dash-dotted line) perpendicular diffusion coefficient. For ⊥ = 0.1 and high rigidities, for instance, the perpendicular mean free path is about a factor 3 smaller compared to the zero gyroradius limit. If the parameter a 2 is understood as finite gyroradius correction, this means that a 2 ≈ 1/3 for the case ⊥ = 0.1 . That is in agreement with what was obtained by  and other papers. However, this only explains the meaning and the value of a 2 in the high rigidity regime. The meaning of a 2 for small rigidities is discussed in Sect. 9.

The Ballistic Approximation
Previous results were directly obtained from time-dependent UNLT theory or via a diffusion approximation. One can combine Eq. (7.22) with the ballistic model ( x) 2 = v 2 μ 2 t 2 δB 2 x B 2 0 (7.40) to find after time integration for the pitch-angle dependent diffusion coefficient For small Kubo number turbulence this corresponds to the FLRW limit with the quasi-linear field line diffusion coefficient. For large Kubo numbers, on the other hand, we can set k ≈ 0 and solve the time integral to derive The wave number integral therein corresponds to the perpendicular integral scale (see, e.g., Eq. (2.57) of this review) so that (7.43) Fig. 13 The running perpendicular diffusion coefficient versus time for slab turbulence, a Kubo number of K = 0.75, and a parallel mean free path of λ = 0.2 . The results were obtained by solving Eq. (7.45) numerically (solid lines). The dotted lines describe compound sub-diffusion as obtained from Eq. (7.50) and the dashed lines the initial ballistic regime as given by Eq. (7.20). In the lower panel the graphs are shown as double-logarithmic plot in order to emphasize the turnover from the ballistic to the sub-diffusive regime. Perpendicular transport in slab turbulence does not depend on the perpendicular bendover scale ⊥ but the units have been chosen so that the results shown in this plot can easily be compared with the results for noisy slab turbulence visualized in Fig. 14 This result is relevant if the particles follow a ballistic motion at the time they experience transverse complexity (see Sect. 9 for more details).

Pure Slab Turbulence
Now we consider pure slab turbulence where we expect compound sub-diffusion. This has to come out of time-dependent UNLT theory as well. For slab turbulence, defined via Eq. (2.13), Eq. (7.7) becomes where we set a 2 = 1. Using Eqs. (5.32) and (7.9) yields for the running perpendicular diffusion coefficient (7.45) Figure 13 shows the numerical solution of this equation for spectrum (2.14) and how this agrees perfectly with compound sub-diffusion. In the following we shall simplify this equation by considering the late-time limit. From Eq. (7.10) it follows that ω − has always a non-vanishing negative real part. Therefore, the exponential exp (ω − t) in Eq. (7.45) damps out for t → ∞. The parameter ω + , on the other hand, can become zero for very small wave numbers k . In this case we can easily derive from Eq. (7.10) as well as Therewith, Eq. (7.45) becomes 16 We can further simplify this by approximating for late times. For the turbulence spectrum we employ Eq. (2.14) and therewith, we derive from Eq. (7.49) corresponding to compound sub-diffusion discussed in Sect. 5.2. Alternatively, we can use the field line diffusion coefficient for slab turbulence as given by Eq. (3.20) to find agreement with Eq. (5.14). Clearly we can see how compound sub-diffusion is correctly described by time-dependent UNLT theory.

Diffusive Perpendicular Transport in Noisy Slab Turbulence
As demonstrated above, time-dependent UNLT theory provides compound sub-diffusion for slab turbulence. In the following we try to explore what ingredients are needed in order to restore normal diffusion. Therefore, we employ a turbulence model which is very close to the slab model but contains some transverse structure. This model is the noisy slab model defined via Eq. (2.25). Combining this with Eq. (7.8) yields where we set a 2 = 1. The k ⊥ -integral can be expressed by an error function 16 Note that Eq. (7.48) was in principle derived the first time in Shalchi (2005b) but in the latter paper a slightly different definition of the diffusion coefficient was used.

Fig. 14
The running perpendicular diffusion coefficient versus time for noisy slab turbulence, a Kubo number of K = 0.75, and a parallel mean free path of λ = 0.2 . The results were obtained by solving Eq. (7.53) numerically (solid lines). The dotted lines describe compound sub-diffusion as obtained from Eq. (7.50) and the dashed lines the CLRR limit represented by Eq. (6.50) for a 2 = 1. In the lower panel the graphs are shown as double-logarithmic plot in order to emphasize the turnover from the sub-diffusive to the normal diffusive regime. The vertical gray line in the bottom plot separates the regimes with ( x) 2 ≤ 2 2 ⊥ and ( x) 2 ≥ 2 2 ⊥ , respectively so that (7.53) The error function has the following asymptotic properties (see, e.g., Abramowitz and Stegun 1974) Erf ≈ 2 √ π x for x 1 1 forx 1.
(7.54) Therewith, it follows from Eq. (7.53) that we find compound sub-diffusion as long as the condition ( x) 2 2 2 ⊥ is satisfied. A numerical solution of Eq. (7.53) is visualized in Fig. 14. Clearly we can see that diffusion is restored for later times. From Fig. 14 but also from Eq. (7.53) we can learn that in order to recover normal diffusion, we need to satisfy the condition ( x) 2 ≥ 2 2 ⊥ (diffusion condition). (7.55) For slab turbulence we have ⊥ = ∞ and we never satisfy this. As soon as there is transverse complexity, corresponding to a finite ⊥ , diffusion is restored as soon as the perpendicular mean square displacement exceeds the perpendicular bendover scale squared. Also shown in Fig. 14 is that the obtained diffusion coefficient is close to the CLRR limit as given by Eq. (6.50). The discrepancy between the latter limit and the numerical solution of Eq. (7.53) is due to the diffusion approximation not being used in the time-dependent theory but it was used during the derivation of diffusive UNLT theory and, thus, it is part of the CLRR limit. It has to be emphasized that the result obtained here is very important. We have seen that everything which is needed in order to get normal diffusion in the perpendicular direction is transverse structure in magnetic turbulence and condition (7.55) needs to be satisfied. The latter condition will be crucial for the heuristic discussion presented in Sect. 9.

Analytical Forms of the Perpendicular Diffusion Coefficient
In several applications in astrophysics and space science simple analytical forms of the perpendicular diffusion coefficient are desired. In the current section we employ the slab/2D model described in Sect. 2.4. This model is supported by observations in the solar wind and theoretical work as explained before. Another strength of this model is that it leads to a strong simplification of fundamental equations. To find analytical solutions for the timedependent case is difficult. However, some analytical solutions of Eq. (7.8) were found in Shalchi (2018). In the following we employ diffusive UNLT theory. First of all we note that the explicit contribution from the slab modes disappears in the late-time limit due to its sub-diffusive behavior. This can, for instance, be obtained from Eq. (6.25). For slab turbulence we have k ⊥ → 0 and, therewith, F (k , k ⊥ → 0) → ∞ and the resulting perpendicular diffusion coefficient is zero. Therefore, the perpendicular diffusion coefficient in slab/2D turbulence is entirely controlled by the two-dimensional modes. However, there could be an implicit contribution of the slab modes as shown in Shalchi (2016a). This contribution is due to the non-linearity of fundamental equations for the perpendicular diffusion coefficient. This effect is automatically taken into account in time-dependent UNLT theory. In the following this effect is neglected in order to obtain simple analytical forms for the perpendicular diffusion coefficient. Furthermore, we neglect finite gyroradius effects.
If we combine Eqs. (6.25) with the spectral tensor for two-dimensional turbulence given by Eq. (2.19), we obtain (8.1) In the following we evaluate the latter equation for two different forms of the spectrum. We like to emphasize that Eq. (8.1) can also be derived from NLGC theory apart from the factor 4/3 in the denominator. This is the case because for pure two-dimensional turbulence, NLGC theory and diffusive UNLT theory are very similar.

The Shalchi & Weinhorst Spectrum
Early investigations of non-linear perpendicular diffusion were based on two-dimensional turbulence and a model spectrum similar compared to Eq. (2.14) which was formulated in the context of slab turbulence (see, e.g., Shalchi et al. 2004a;Zank et al. 2004). In the following we consider the more general spectrum proposed by Shalchi and Weinhorst (2009) (see Eq. (2.20) of this review). The latter spectrum incorporates requirements such as having finite length scales of turbulence originally formulated in Matthaeus et al. (2007). With spectrum (2.20) the perpendicular diffusion coefficient given by Eq. (8.1) becomes ( 8.2) To solve this non-linear integral equation we employ the methods presented in . With the integral transformation x = k ⊥ ⊥ and by using the parameter where we have used Here we have also replaced the perpendicular diffusion coefficient by the perpendicular mean free path. With the integral transformation y = x 2 and by using Gradshteyn and Ryzhik (2000), the integral in Eq. (8.5) can be expressed by a product of a beta function B(x, y) and a Gaussian hypergeometric function 2 F 1 (a, b; c; z) so that Note, the integral in Eq. (8.5) is only convergent for q > −1. Thus, we only consider such values of the energy range spectral index. The perpendicular mean free path λ ⊥ is obtained if Eq. (8.6) is combined with Eq. (8.4). However, the parameter itself depends on λ ⊥ and, thus, we have to deal with a non-linear equation. To simplify the expression found here, we consider three different cases of the two parameters and q.

The Case
1 and −1 < q < 1 First we rewrite the hypergeometric function in Eq. (8.6). In order to do this we employ the transformation (see, e.g., Abramowitz and Stegun 1974) allowing us to write In the case discussed in this paragraph we consider 1 and, therefore, we use (see Abramowitz and Stegun 1974) as well as (1) = 1 to obtain . (8.10) Whether the first or second term in Eq. (8.10) dominates depends on the value of q. For the case that −1 < q < 1, the hypergeometric function can be approximated by the first term in Eq. (8.10). If we additionally replace the beta function by (see again Abramowitz and Stegun 1974) B(x, y) = (x) (y) (x + y) (8.11) Equation (8.6) turns into giving us the perpendicular mean free path as a function of the parallel mean free path. In the case considered here, the parameter is assumed to be small. Going back to Eq. (8.3) tells us that this means 4λ λ ⊥ 9 2 ⊥ . Since both mean free paths increase with rigidity, this result should be important for high particle rigidities. Interesting here is the dependence on the parallel mean free path. For q = 0, for instance, we find λ ⊥ ∝ λ 1/3 as originally obtained in Shalchi et al. (2004a) and Zank et al. (2004).

The Case
1 and q > 1 In the following we still assume that is small as before. However, we now consider the case q > 1 and, thus, the second term in Eq. (8.10) is dominant and we have .
(8.15) Using this term in Eq. (8.6) yields where we have used Eq. (2.67) as well. Using this result in Eq. (8.4) gives us Alternatively this result can be obtained by combining Eq. (5.6) with Eq. (3.56). Equation (8.18) corresponds to the FLRW limit in the case where the field lines are highly non-linear.

The Case 1
So far we considered cases where the parameter is small. In the last out of the three cases we now explore 1 corresponding to small rigidities. After employing the so-called Pfaff transformation (see, e.g., Abramowitz and Stegun 1974) the hypergeometric function in Eq. (8.6) can be written as By additionally using (Abramowitz and Stegun 1974) in Eq. (8.20), we get In the last step we have used Eq. (2.67) as well. The beta function in Eq. (8.6) can also be written as and for the perpendicular mean free path, given by Eq. (8.4), we therefore find This case was also obtained before (see Eq. (6.45) of this review). Above we have shown how different limits come out of diffusive UNLT theory for twodimensional turbulence. Equations (8.18) and (8.25) were derived before in this review by considering the corresponding limits directly. Equation (8.14), however, disagrees with those directly obtained limits. The reason for this is that Eq. (8.14) was obtained for a case where the ultra-scale is not finite. This would also mean, as shown in Shalchi and Kourakis (2007b) and Shalchi (2011c), that the magnetic field lines are super-diffusive for the case q < 1. If one drops the Matthaeus et al. (2007) requirement that the ultra-scale needs to be finite, Eq. (8.14) would still be a valid result. However, during the remainder of this review, we focus on spectra satisfying those requirements and this means that Eqs. (8.18) and (8.25) are the only considered asymptotic limits for perpendicular transport in two-dimensional turbulence.

The Simplified Spectrum
Above we derived simple formulas for the perpendicular diffusion coefficient and a general and realistic spectrum. However, the derived formulas are only valid in asymptotic limits. An alternative approach is to consider a special spectrum with q = 2 and s = 2 (see Shalchi 2013b). For the choice q = 2 the integral scale and the ultra-scale are finite (see, e.g., Table 2 of this review) and the random walk of magnetic field lines is normal diffusive (see Sects. 2 and 3 of this review). For these values of the spectral indices we have D(s = 2, q = 2) = 1/π and the spectrum (2.20) becomes (8.26) In the following we call this the simplified spectrum but we keep in mind that this spectrum contains all features important for describing perpendicular diffusion in two-dimensional turbulence. In order to compute the perpendicular diffusion coefficient we can employ again Eq. (8.1). In combination with the simplified spectrum we find where we have used the integral transformation x = k ⊥ ⊥ and the parameters α = 4κ ⊥ /(3 2 ⊥ ) as well as β = v/λ . We can solve the integral by using (see, e.g., Gradshteyn and Ryzhik 2000) Using this in Eq. (8.27) and by replacing κ ⊥ = vλ ⊥ /3 we find (8.29) By additionally using Since x is defined as positive number, the latter quadratic equation has only one solution.
By replacing x, α, and γ again by physical quantities, this solution is The latter formula can easily be employed to compute the perpendicular mean free path if the parallel diffusion coefficient and turbulence parameters are known. The strength of Eq. (8.33) is that we now have one single formula for the perpendicular mean free path in two-dimensional turbulence. The cases given by Eqs. (8.18) and (8.25) can easily be restored by considering the corresponding limits. Furthermore, Eq. (8.33) is an exact solution of diffusive UNLT theory for two-dimensional turbulence and the simplified spectrum. In Fig. 15 this solution is visualized together with the numerical solution of Eq. (8.1) for the more realistic case s = 5/3 and q = 3 as well as the two asymptotic solutions given by Eqs. (8.18) and (8.25).

The Heuristic Description of Perpendicular Diffusion
In the previous sections we discussed advanced systematic theories for perpendicular transport. As it was shown over the past few years, time-dependent UNLT theory in particular agrees with most test-particle simulations (see, e.g., Arendt and Shalchi 2018). However, there are two remaining questions which need to be explored. First, there are cases where the agreement between theory and simulations can only be achieved if the parameter a 2 is set to 1/3. Strictly speaking the theory fails in such cases because ideally this type of correction parameter should not be needed. Furthermore, in physics one usually desires to achieve a complete understanding of fundamental processes. The question arises how and why are particles experiencing perpendicular diffusion. In the current section we, therefore, review the heuristic approach towards perpendicular transport developed in Shalchi (2019a). This will not only explain the underlying physics of perpendicular diffusion but it will also provide an explanation of why the correction factor a 2 is sometimes needed.

The Three Rules of Perpendicular Diffusion
Before we derive formulas for the perpendicular diffusion coefficient based on heuristic arguments, we discuss the three rules of perpendicular diffusion formulated in Shalchi (2019a). Those are: 1. Perpendicular transport is only controlled by three effects, namely parallel transport, the random walk of magnetic field lines, as well as transverse complexity. The last of these three effects leads to the particles getting scattered away from the original magnetic field lines they were tied to. 2. We assume that the bendover scales and ⊥ , the integral scales L and L ⊥ , the ultrascale L U , as well as the Kolmogorov scale L K are finite and non-zero. Furthermore, the parallel motion is assumed to be ballistic at early times and thereafter turns into a diffusive motion described by the parallel diffusion coefficient κ . The FLRW is initially ballistic and becomes diffusive for larger distances. In this case it is described by the field line diffusion coefficient κ FL which depends on some of the aforementioned scales.  Table 7 The different routes leading to perpendicular diffusion. In the first three cases perpendicular transport starts as ballistic motion which then turns directly into a normal diffusive motion. Only in the fourth case the ballistic motion is followed by a sub-diffusive regime and then, at later times, diffusion is restored. This table is from Shalchi (2019a)

Route
Final state Diffusion coefficient 3. In order to obtain normal diffusion, the particles need to leave the original magnetic field lines they followed. This happens as soon as transverse complexity becomes significant corresponding to What the perpendicular diffusion coefficient is depends solely on the state of parallel and field line transport at the time particles start to satisfy condition (9.1).
It is assumed here that ⊥ is the scale at which transverse complexity becomes significant.
In principle this could be a different scale such as the integral scale L ⊥ . To use the bendover scale, however, is motivated by time-dependent UNLT theory (see Sect. 7 and in particular Eq. (7.55)). In the following we try to understand perpendicular transport based on the three rules formulated above. This will require to consider eight different cases which are summarized in Table 6. The different routes to perpendicular diffusion are listed in Table 7.

The FLRW Limit
First we assume that the random walk of magnetic field lines is diffusive at time scales where perpendicular transport is of interest. This means we assume where we have used the diffusion coefficient of magnetic field lines κ FL . In order to continue we need to model the parallel motion of the particles. If we assume that there are no collisions and no pitch-angle scattering, we can use z(t) = vμt where we used the pitch-angle cosine μ and the particle speed v as before. Combining this with Eq. (9.2) and averaging over μ yields corresponding to the field line random walk limit as given by Eq. (6.36). This limit is stable because if condition (9.1) is met, it does not alter the transport. Sometimes this type of transport is called first diffusion and it is highly relevant in the limit of very long parallel mean free paths. However, it needs to be emphasized that it has to be understood as asymptotic limit which is never really achieved (see, e.g., Figs. 21 and 22 of this review).

Compound Sub-Diffusion
The FLRW limit was obtained by assuming that the parallel motion is unperturbed meaning ballistic. However, if there is strong pitch-angle scattering the parallel motion is diffusive meaning that Assuming that at a certain time field lines are diffusive and particles follow field lines, we can combine Eqs. (9.2) and (9.4) to find The running diffusion coefficient can then be obtained via 17 corresponding to sub-diffusive transport. In other words: if nothing else happens and we just have a combination of particles following random walking magnetic field lines, while they move diffusively in the parallel direction, perpendicular transport will be sub-diffusive forever. However, diffusion will be restored as soon as condition (9.1) is satisfied as discussed in the next paragraph. For slab turbulence, on the other hand, this condition is never satisfied due to ⊥ = ∞. Thus, for slab turbulence, we find compound sub-diffusion as the final state of the transport.

The CLRR Regime
We now assume that diffusion starts as soon as the particles scatter away from their original magnetic field line. This happens as soon as condition (9.1) is satisfied. We also assume that this happens after the particles travel the distance L K in the parallel direction. Therefore, we have (9.7) In order to replace L K we can use the field line diffusion coefficient where The latter formula can be rewritten as Using this in Eq. (9.7) yields Alternatively, Eq. (9.9) can be used in order to replace the perpendicular bendover scale ⊥ in Eq. (9.7) so that κ ⊥ ≈ κ FL κ L K (9.11) in agreement with equation (8) of Rechester and Rosenbluth (1978) as well as equation (4) of Krommes et al. (1983). This formula was also derived in Shalchi (2019b) by using arguments based on field line separation theory. The quantity L K is sometimes called the Kolmogorov length or the Kolmogorov-Lyapunov length (see, e.g., Krommes et al. 1983;Spatschek 2008). However, here L K is not an exponentiation length but a characteristic distance along the mean field at which transverse complexity becomes significant. The work by Rechester and Rosenbluth (1978) emphasized the importance of collisions. In space plasmas collisions are not important but there is strong pitch-angle scattering leading to parallel diffusion as discussed in Sect. 4 of this review. Thus, the importance of parallel diffusion remains regardless whether it is caused by collisions or pitch-angle scattering. Thus, we call Eq. (9.10) the collisionless Rechester & Rosenbluth (CLRR) limit. 18 One can also obtain Eq. (9.10) by using a slightly different derivation. Let's assume that we are in a regime where particles experience compound sub-diffusion. Then we can use Eqs. (9.5) and (9.6) to describe the transport. We now assume that we find compound subdiffusion until the particles satisfy condition (9.1) which is assumed to happen at the time t d so that Eqs. (9.5) and (9.6) turn into as well as Combining the latter two equations in order to eliminate t d yields again Eq. (9.10).
In order to evaluate Eq. (9.10) further, we consider two cases, namely small and large values of the Kubo number, respectively. For small Kubo numbers the field line diffusion coefficient is given by Eq. (3.21) leading to as well as (9.15) Equation (9.14) is in perfect agreement with the scaling obtained from diffusive UNLT theory (see, e.g., Eq. (6.48) of the current review). For large Kubo numbers, on the other hand, we use Eq. (3.59) in Eq. (9.10) yielding Furthermore, the Kolmogorov scale is, in this case, Since we usually expect L U ∝ ⊥ (see, for instance, Sect. 5), Eq. (9.16) looks very similar compared to the fluid limit (see next subsection). However, the way how particles move across the field is different. Below we will evaluate this further and relate it to previous results and test-particle simulations.

The Fluid Limit
Let's now assume that parallel transport is diffusive but magnetic field lines are still ballistic when the particles start to satisfy condition (9.1). Then we can easily derive Consequently, the perpendicular diffusion coefficient is given by corresponding to the fluid limit. We conclude that the fluid limit is obtained if the parallel motion becomes diffusive when the field lines are still ballistic.

Ballistic Perpendicular Transport
The simplest case is obtained for the very early times when parallel transport and field lines are ballistic. In this case (9.21) This means that the running diffusion coefficient increases linearly with time. This type of transport can be seen in test-particle simulations for very early times (see, e.g., Sect. 10 for some examples). However, this is not a stable regime since we only find this type of transport before condition (9.1) is met.

The Double-Ballistic Regime
We now consider the case where perpendicular transport is ballistic as described in Sect. 9.6 when particle experience transverse complexity. As soon as we satisfy condition (9.1), normal diffusion is restored. Therefore we use Eqs. (9.20) and (9.21) to derive as well as Combining the latter two equations in order to replace the diffusion time t d leads to very similar compared to Eq. (7.43) derived from time-dependent UNLT theory by assuming a ballistic motion. However, Eq. (7.43) contains the perpendicular integral scale L ⊥ rather than the perpendicular bendover scale ⊥ . Therefore, we need to keep in mind that the heuristic approach is based on simple arguments and is not necessarily providing a very accurate result.

Time-Scale Arguments
Above we derived several cases for perpendicular transport. Of course the question arises which case is valid for which scenario. In order to answer this question one needs to explore how long it takes for a certain process to take place. If it comes to parallel transport, for instance, one usually assumes that the particles need to travel a parallel mean free path in order to get diffusive (see, e.g., Eqs. (4.78) and (4.79) of this review). Therefore, we can estimate To find a condition for magnetic field line diffusion is more difficult. For small Kubo numbers the field lines usually become diffusive for |z| ≈ (see, e.g., Sect. 3.4). In this case one needs to specify how the particles move when they satisfy this condition. For a very short parallel mean free path, parallel diffusion is reached quickly. Therefore the field lines become diffusive for Then, on the other hand, if we assume that condition (9.1) is satisfied while the field lines are still ballistic, we can estimate (9.28) In order to find the fluid limit as the final state of perpendicular diffusion, we need to satisfy because only then we find that parallel transport becomes diffusive first, then we meet condition (9.1) and we obtain a stable regime for perpendicular transport, namely the fluid limit. If, on the other hand, the field lines become diffusive before condition (9.1) is met. This means that we find compound sub-diffusion first. At even later time condition (9.1) is eventually met and diffusion is restored. The corresponding diffusion coefficient is then the one corresponding to the CLRR limit. Using the formulas for the times discussed above, this means that we find CLRR diffusion for where we have used the Kolmogorov length L K given by Eq. (9.15). This means that for λ we either find the fluid limit or CLRR diffusion. If additionally L K we find the fluid limit but for L K we get CLRR diffusion. We can clearly see that L K / ≈ K −2 1 meaning that for small Kubo numbers we should always find CLRR diffusion rather than the fluid limit. This is in principle in agreement with UNLT theory where we find CLRR diffusion for small λ and small K but the fluid limit for small λ and large K.
The arguments presented above were performed for small Kubo number turbulence. For large Kubo numbers, on the other hand, it is more difficult to estimate the time t FL since it is unclear for which distance the field lines become diffusive. However, we expect that this time will be somewhere between the t FL derived above and t F luid . We can easily see that if the Kubo number is large, the time t F luid given by Eq. (9.28) becomes small. We therefore, conclude that if the Kubo number is large enough one should find the fluid limit but in most applications we expect to find CLRR diffusion.
Also problematic is that one can easily end up in a regime somewhere between CLRR and fluid limits. If we focus on the large Kubo number regime, Eqs. (9.16) and (9.19) are very similar. One would expect that the fluid limit acts as an upper limit meaning that if L U > ⊥ one finds the fluid limit.

Further Comments
Important here is to note that the results obtained here are sometimes not comparable to previous results. However, this is because previously obtained results often violate one of the three rules of perpendicular diffusion listed above. First of all there are cases such as pure slab or pure two-dimensional turbulence. In the former case condition (9.1) is never satisfied leading to compound sub-diffusion as the final state. In the two-dimensional case parallel transport is not diffusive (see, e.g., Arendt and Shalchi 2018) violating the second rule.
In some work (see, e.g., Qin et al. 2002aQin et al. , 2002bor Shalchi et al. 2004a) a spectrum was used for the two-dimensional modes which is perfectly flat at large scales. For this type of spectrum some fundamental scales of turbulence are not finite. This is in particular the case for the ultra-scale and, thus, it violates the second rule of perpendicular diffusion.
Sometimes the FLRW limit as given by Eq. (6.36) together with Eq. (3.59) is called the Kadomtsev & Pogutse limit. However, there is another limit which is sometimes also called the Kadomtsev & Pogutse limit (see, e.g., Table 1 of Krommes et al. 1983) which is where χ ⊥ is the perpendicular diffusion coefficient due to collisions. However, in astrophysical scenarios there are no collisions and perpendicular diffusion is only caused by pitchangle scattering and transverse complexity. Thus we set χ ⊥ = κ ⊥ so that Eq. (9.33) turns into κ ⊥ = κ 2 FL κ / 2 ⊥ in perfect agreement with Eq. (9.10). Therefore, we conclude that in the collisionless case the second Kadomtsev & Pogutse limit and the Rechester & Rosenbluth limit are the same.
In order to determine the form of the field line diffusion coefficient and to distinguish between the different regimes of perpendicular diffusion, we use the Kubo number as given by Eq. (3.11). However, in some cases such as the Goldreich & Sridhar model (see Sect. 2.8) there is only one scale and, thus, = ⊥ . In this case the Kubo number becomes K = δB x /B 0 often called the Alfvénic Mach number. The results derived above are, of course, still valid.

A Composite Formula
A problem of the heuristic approach is that the obtained formulas are only valid in asymptotic limits. Therefore, we now design a formula which contains the two most important cases as asymptotic limits. Those are the CLRR and FLRW limits. Motivated by Eq. (8.33) we start with the ansatz We now choose the constants b and c so that for λ → 0 we obtain Eq. (9.10) and for λ → ∞ we get Eq. (6.36). After straightforward algebra we obtain Of course this formula has its limitations. For instance it does not contains the fluid limit and, thus, it is not valid for too large Kubo numbers. However, it can be used for arbitrary turbulence. In Sect. 10 we compare this formula with UNLT results and test-particle simulations.

Two-Component Turbulence as an Example
As an example we employ a two-component turbulence model where we approximate turbulence by using a superposition of slab and two-dimensional modes. In this case perpendicular transport is mostly controlled by the two-dimensional modes whereas the slab modes ensure that parallel transport behaves well. One still needs to specify the spectrum in order to obtain concrete results.
In the case of two-dimensional turbulence we employ Eq. (2.19) for the spectral tensor and Eq. (2.20) for the spectral function. According to Eq. (2.66) the ultra-scale is, in this case, given by and, thus, it is directly proportional to the perpendicular bendover scale but depends also on the two spectral indices.
With the parameter a 2 included, non-linear theories provide for the perpendicular diffusion coefficient in the limit of short parallel mean free path (see, e.g., Eq. (6.45) of this review) If one is in the fluid regime we have a 2 = 1. According to the heuristic approach described above, however, we expect more a CLRR type of behavior. Comparing Eqs. (9.39) and (9.16) yields a = L U ⊥ (9.40) and using Eq. (9.38) for the ultra-scale gives us Previously, as quantitative theories were compared with simulations, it was often assumed that s = 5/3 and q = 3 (see for instance Arendt and Shalchi 2018). Therefore, we can easily show that The arguments presented here provide a natural explanation of why the parameter a 2 is needed. This parameter balances out the differences of the fluid limit and CLRR diffusion in the large Kubo number regime. This is important because the problem of needing this correction factor a 2 in systematic theories is known since more than 15 years (see . Due to the heuristic approach discussed here we now understand the physics of this parameter as well as its value. It should be emphasized that Eq. (9.42) suggests that a 2 is exactly 1/3. This is probably not true. First of all the heuristic approach is based on rough estimations and is not supposed to be exact. Furthermore, Eq. (9.39) with a 2 = 1/3 is already close to the fluid limit where a 2 = 1. The parameters often used in the context of slab/2D turbulence lead to a scenario where perpendicular transport ends up somewhere between CLRR and fluid limits. One would expect that changing turbulence and particle parameters could easily lead to different values of a 2 and in some cases one would even expect a 2 = 1.
More examples and a comparison between UNLT theory, the heuristic approach, and test-particle simulations is presented in the next section (see Figs. 21 and 22).

Test-Particle Simulations
The work discussed so far is entirely based on analytical considerations. An alternative tool for describing the transport of energetic particles through magnetic turbulence is provided by test-particle simulations. They are also useful in order to test the validity and accuracy of the analytical theories discussed so far. Ideally there is agreement between theory, simulations, and observations. Test-particle simulations consist of the following three steps: 1. First we need to generate the magnetic fields by employing a turbulence model very similar compared to the models used in analytical theories (see Sect. 2 of this review). A problem is that the computational time needed for turbulence creation increases significantly if three-dimensional or time-dependent models are employed. 2. The second step is to solve the Newton-Lorentz equation (4.1) to obtain the particle orbits. Since the Newton-Lorentz equation is a second-order ordinary differential equation, one usually employs standard solvers such as fourth-order Runge-Kutta integrators with adaptive time steps (see, e.g., Press et al. 2007). An alternative is provided by symplectic integrators leading to better energy conservation (see Arendt and Shalchi 2018). 3. The last step is to ensemble average the obtained orbits to compute the diffusion parameters in the different directions of space. Furthermore, one can also obtain particle distribution functions, mean square displacements, and velocity correlation functions. In order to get results with high accuracy and to reduce noise, one has to perform such simulations by using thousands of particles.
In the following we present a more detailed description of these steps and thereafter we show numerical results together with a comparison with analytical results for different turbulence configurations.

Generating the Turbulence
As described above, test-particle simulations are a common tool in the theory of particle transport. In the following we summarize the method and employ the notation used before in Hussein and Shalchi (2014a), , Hussein and Shalchi (2016) as well as Arendt and Shalchi (2018). This approach, however, is based on previous papers such as Giacalone and Jokipii (1999) as well as Tautz (2010b). The most common technique to generate magnetic turbulence is to sum over a large number of plane waves with random phases, polarizations, amplitudes, and orientations. In test-particle simulations the turbulent magnetic field vector is created via where N is the number of wave modes. The quantity A(k n ) denotes the amplitude function and is discussed below. In Eq. (10.1) we have also used the wave vector k n = k n e k,n with the Table 8 The values of the parameters used in the simulations for the different turbulence models. Here we have listed the parameters η n , α n , and n used in Eqs. (10.1)-(10.3) in order to create the turbulent magnetic field. The parameters , ⊥ , and 0 are the corresponding bendover scales. The used values for the energy range spectral index q are also listed. The first three cases are based on a single sum in the field creation whereas the last two cases are based on a double-sum Turbulence model η n α n n Wave numbers q wave number k n and the unit vector (10.2) Furthermore, we have used the random phase β n , restricted by 0 ≤ β n < 2π , and the polarization vector − sin φ n cos α n + η n cos φ n sin α n cos φ n cos α n + η n sin φ n sin α n − 1 − η 2 n sin α n ⎞ ⎟ ⎠ (10.3) with η n = cos θ n . The angles θ n , φ n , and α n can have a specific value or they are random numbers depending on the simulated turbulence model (see Table 8 of the current paper for some examples). If the slab model is simulated, for instance, we set η n = 1 so that e k,n = e z for all n. Furthermore, we set α n = 0 so that ξ n = − sin φ n e x + cos φ n e y for all n. This choice of η n and α n ensures that all wave vectors are aligned parallel with respect to the mean field and that δB z = 0. These are exactly the conditions which need to be satisfied if the slab model is considered. The remaining angles φ n are random numbers satisfying 0 ≤ φ n < 2π emulating the chaotic nature of the field δB.
In a very similar manner we can generate two-dimensional turbulence. In this case we set α n = 0 as before to ensure that δB z = 0. However, we now use η n = 0 so that e k,n = cos φ n e x + sin φ n e y and ξ n = − sin φ n x + cos φ n e y as required for two-dimensional turbulence.
Isotropic turbulence can also be generated via Eqs. (10.1)-(10.3). In this case all angles θ n , φ n , and α n are random numbers leading to a turbulence model where the wave vector as well as the magnetic field vector are isotropic.
One can easily show by combining Eqs. (10.2) and (10.3) that for the general case we have e k,n · ξ n = 0, corresponding to the solenoidal constraint k · δB = 0, as required. The amplitude function A(k n ) used in Eq. (10.1) depends on the spectrum G(k n ) via (10.4) The parameter k n describes the spacing between wave numbers and is discussed below. Note, the wave numbers used in simulations are unitless meaning that the physical wave numbers are multiplied by a characteristic scale of turbulence, usually one of the bendover scales. This means, for instance, that for slab turbulence the parameter k n in Eq. (10.4) is really k n . The amplitude function (10.4) has to be normalized so that For the spectrum G(k n ) we use a form corresponding to the analytical models described in Sect. 2, namely The parameters q and s are energy and inertial range spectral indices as before. The used values for these two parameters and the meaning of k n in the different models are summarized in Table 8. In most simulations a logarithmic spacing in k n is implemented so that which is constant. In the context of two-dimensional turbulence the grid created in wave number space via Eq. (10.8) is visualized in Fig. 3 together with the spectra used in analytical treatments of field line and particle transport. 19 Furthermore, the possible values for the wave numbers are restricted by k min ≤ k n ≤ k max . Typical values for minimum and maximum wave numbers are k min = 10 −5 and k max = 10 3 . There are two constraints that should be taken into account in test-particle simulations. First we know that pitch-angle scattering happens mostly close to the resonance condition μR L k = 1 (see Sect. 4.3). Therefore, the simulations have to be performed so that we hit the resonance. Furthermore, we need to ensure that no particle travels more than the distance L max = k −1 min . This is done via the relation vt max < L max corresponding to a restriction of time. For k min = 10 −5 , for instance, we have L max = 10 5 leading to the condition vt max < 10 5 . Since in the simulations we use T = t for time, this turns into T max R L / < 10 5 . Obviously, this restriction becomes more relevant for high rigidities.

More General Models for Turbulence
As obvious from Eq. (10.1), we only discuss the case of magnetostatic turbulence here. Of course, one can incorporate wave propagation effects in such simulations (see, e.g., Michałek and Ostrowski 1996;Tautz 2010b;Hussein and Shalchi 2014b). Test-particle simulations in dynamical turbulence were performed for the first time in Hussein and Shalchi (2016) based on a more-dimensional Fourier technique. Furthermore, the models discussed above were models with either reduced dimensionality (slab or two-dimensional modes), superpositions thereof (two-component turbulence), or isotropic models. Therefore, the simulated wave numbers are only one-dimensional and thus there is only one sum over the different modes in Eq. (10.1). There are in principle two ways of generalizing this to three-dimensional anisotropic but axi-symmetric turbulence: 1. One can keep Eqs. (10.1)-(10.8) and adjusts the spectrum and the values of the parameters η n , φ n , and α n accordingly. Of course one has to ensure that the whole wave number space of the simulated turbulence model is covered and that the right spectrum is chosen.
Since a more-dimensional wave number space has to be covered with the same spacing k n , significantly more wave modes N are needed. This is the method which was used, for instance, by Tautz (2010aTautz ( , 2010b. 2. One can follow  and create the turbulent magnetic field vector via the double-sum to be modified as well but this is straightforward because one just has to replace the sums in these equations by double-sums. The second method was specifically developed for incompressible turbulence. It is in particular useful for so-called noisy models where one starts with a model with reduced dimensionality such as the slab model and then adds some noise. In the case of the noisy slab model, for instance, the turbulence is very narrow in the perpendicular direction. Therefore, one has to generate only a few wave numbers k n in Eq. (10.9) but a large number of parallel wave numbers k m . The main problem is that the creation of anisotropic turbulence is much more time-consuming computationally. With the help of modern super-computing, however, one can overcome this problem and perform simulations for time-dependent or anisotropic turbulence without too many problems.

Solving the Newton-Lorentz Equation
The Newton-Lorentz equation for a particle in purely magnetic turbulence is given by Eq. (4.1). This equation is solved numerically in test-particle simulations with the initial conditions discussed in the following. Initially all particles have the same z-coordinate and the same rigidity which is conserved in a purely magnetic system. However, their initial x-and y-coordinates as well as their initial pitch-angle cosine μ 0 are different. There are, in principle, two ways of implementing the turbulent magnetic fields in test-particle codes: 1. Using discrete grid: In this case we first create and save magnetic field data at each grid point using Eqs. (10.1)-(10.8) for the whole space and then interpolate for where the particle is moving. This could be done using a two-dimensional grid for the perpendicular directions accompanied with a one-dimensional grid for the parallel direction or rigorously using a three-dimensional grid. The grid method was used by Mace et al. (2000), Casse et al. (2002), Qin et al. (2002a,b), Pommois et al. (2007), and Reville et al. (2008). 2. Creating fields along the trajectory: An alternative is to create fields anew at each time step. The given initial position allows us to create the field initially which is then seeded back to the numerical integrator to solve for position which is then seeded to turbulence creation and so on and so forth. This type of simulations were performed by Giacalone andJokipii (1994, 1999), Tautz (2010aTautz ( , 2010b, Hussein and Shalchi (2014a), and Arendt and Shalchi (2018).
The second method listed above saves time and uses less memory compared to the grid method because it generates magnetic fields only where the particle is actually moving, not on all the provided space as the first method does. On the other hand, when it comes to visualization, the grid-based system allows to visualize the magnetic field lines across the whole space and therefore one can see how particles are moving in the vicinity of field lines. This could be useful in order to develop a deeper understanding of the physics of particle transport.
In order to solve the second-order Newton-Lorentz equation one can use a fourth-order Runge-Kutta solver with an adaptive time step option. Although this can be seen as a standard method in this field, more recently a modified third-order symplectic integration method was used as an alternative (see Arendt and Shalchi 2018). This ensures energy conservation and should provide an important improvement of test-particle simulations if stochastic acceleration due to turbulent electric fields is studied.
In simulations parameters are made to be dimensionless. For instance, all length scales are divided by the turbulence bendover scale . Furthermore, we define the dimensionless running time via T = t and the dimensionless rigidity vector via R := v/( ). With these parameters, we can derive the dimensionless Newton-Lorentz equation where the turbulent field δB(x) is given by Eq. (10.1). The relation between position and velocity v = dx/dt turns into the dimensionless equation Special care is required if there is more than one bendover scale . For two-component turbulence, we usually choose = . Then, however, if the two-dimensional modes are created via Eq. (10.1), positions are measured in terms of the slab bendover scale . In the two-dimensional modes we also have k n = k ⊥ ⊥ , therefore, the ratio ⊥ / appears in front of x · k and this scale ratio controls the transport of particles. Reprinted with permission from Springer- Arendt and Shalchi (2018) In the numerical solution of Eqs. (10.12) and (10.13), one needs to specify several parameters as well. The initial time is usually set to zero but there is also a final time t max . This has to be chosen so that one finds the stable regime which is often the time where the particles have reached diffusive behavior. In Figs. 16 and 17, for instance, the choice was t max = 5000. Furthermore, the constant step size in the symplectic solver was set to t = 10 −3 meaning that the total number of time steps in such runs was 5 × 10 6 . The procedure explained so far has to be performed for a huge amount of particles to obtain results with a high accuracy and in order to reduce the noise in the running diffusion coefficients as much as possible. Often the number of particles is a few thousands but the results visualized in Figs. 16 and 17 were created by using 12000 particles. That is the point where parallel computing becomes a required tool so that the different particles can be distributed among the different processors.
If particle trajectories are obtained numerically, the remaining step is the calculation of the diffusion coefficients via mean square displacements. To do this, a diffusion coefficient is preferably defined as the ratio of the corresponding mean square displacement and time (see, e.g., Eq. (4.80) of this review) rather than the time-derivative (see, e.g., Eq. (5.27)). This is entirely done with the purpose of reducing noise in the diffusion coefficient. If one is more interested in the late-time limit and assuming that one indeed finds diffusive transport, there is no difference between dividing by time and computing the time-derivative. However, special care is required for anomalous transport where these two definitions do not yield the  R,K ,and D are defined in Eq. (4.83). In the upper left panel the solid line represents the test-particle simulations and the dotted line the analytical formula (4.84) for K = 0.60. In the upper right panel the solid line represents the test-particle simulations, the dotted line the numerical solution of Eq. (7.8), corresponding to time-dependent UNLT theory for a 2 = 1, and the dashed line for a 2 = 1/3. The bottom panels show parallel and perpendicular distribution functions for the different times T = 0, T = 2500, and T = 5000. Shown are the simulations (solid lines) and Gaussian overlays (dashed lines). Reprinted with permission from Springer- Arendt and Shalchi (2018) same diffusion coefficients. Of course, if numerical and analytical data are compared with each other, one has to use the same definition of a diffusion coefficient in both cases.
In addition to running diffusion coefficients in both directions of space, Arendt and Shalchi (2018) also obtained particle distribution functions giving more insight into the transport of particles. Such distribution functions are shown in Figs. 16 and 17 together with Gaussian overlays.
In the following we consider test-particle simulations for different turbulence configurations. In all considered cases we compare simulations and analytical results with each other. The results of this comparison are summarized in Table 9.

Slab Turbulence
As a first example we explore the transport in magnetostatic slab turbulence. The aim is to explore whether perpendicular transport in indeed sub-diffusive as described by compound sub-diffusion models. Furthermore, we can test time-dependent UNLT theory for that specific case. Arendt and Shalchi (2018) performed simulations for slab turbulence and not just parallel and perpendicular diffusion coefficients were obtained, but also the corresponding distribution functions of the energetic particles. In the simulations performed for slab turbulence q = 0 and s = 5/3 were used in the spectrum defined via Eq. (10.7). Therefore, the simulated spectrum agrees with model spectrum (2.14). In Fig. 16 we show the simulations Table 9 Diffusive and time-dependent UNLT theories have been tested by comparing perpendicular diffusion parameters with test-particle simulations in several articles. The table lists those comparisons as well as the parameter values for a 2 for which the best agreement was obtained. Ideally we have a 2 = 1 meaning that the correction factor a 2 is not needed at all demonstrating that we find indeed sub-diffusive perpendicular transport. The conclusion was also obtained by Giacalone and Jokipii (1994) and Giacalone and Jokipii (1999) as well as Mace et al. (2000) and Qin et al. (2002a). Furthermore, the simulations shown in Fig. 16 confirm the parallel transport model listed in Eq. (4.84). We can also see that the parallel distribution function is Gaussian. However, the perpendicular distribution function is clearly non-Gaussian. Webb et al. (2006) used the Chapman-Kolmogorov equation (5.10) to compute the distribution function across the mean field and they found a result depending on a so-called Fox function. This distribution was plotted in Fig. 3 of Webb et al. (2006) and is very similar compared to the distribution obtained and visualized in Fig. 16. We can also clearly see that time-dependent UNLT theory, represented by Eq. (7.45), agrees very well with the simulations for all considered times. This includes the initial ballistic regime as well as the peak and the following sub-diffusive regime. It needs to be emphasized that the theory does not include any free parameter since in Eq. (7.45) we set a 2 = 1.

Slab/2D Composite Turbulence
In order to restore Markovian diffusion one needs to include transverse complexity. Starting with the slab model this can be done in two ways. First we can add two-dimensional modes to create a quasi three-dimensional model, or we broaden the slab model leading to the noisy slab model. The first option is discussed in the following whereas the second option is considered in Sect. 10.6.
The two-component model describes the turbulence as superposition of slab and twodimensional modes. For the slab modes we set q = 0 as before and for the two-dimensional modes q = 3. Furthermore, we set δB 2 slab /δB 2 = 0.2 and δB 2 2D /δB 2 = 0.8 as suggested by . Figure 17 shows simulations performed for slab/2D turbulence. Clearly we can see that parallel and perpendicular transport behave diffusively in the late-time limit. Again we can confirm the parallel transport model given by Eq. (4.84). Both, the parallel as well as the perpendicular distribution functions are now Gaussian. We can clearly see that time-dependent UNLT theory agrees well with the simulations if we set a 2 = 1/3 in front of the two-dimensional modes. In front of the slab modes we kept a 2 = 1. As demonstrated, the theory describes the initial ballistic regime as well as the following diffusive regime correctly.
In  test-particle simulations were used to test NLGC and diffusive UNLT theories for two-component turbulence. In Fig. 18 the influence of the slab fraction δB 2 slab /δB 2 on the perpendicular diffusion coefficient is shown. The simulations indicate sub-diffusive transport for pure two-dimensional turbulence as well as for pure slab Fig. 18 The perpendicular mean free path for two-component turbulence versus the slab fraction δB 2 slab /δB 2 . Visualized are the test particle simulations (dots), NLGC theory (dashed line), and diffusive UNLT theory (solid line). The arrows indicate that the simulations provide a sub-diffusive result in the considered case. All analytical results shown here were obtained by setting a 2 = 1/3. The simulations are from  turbulence. Therefore, the perpendicular diffusion coefficient should go to zero for high values of the slab fraction. This is exactly what diffusive and time-dependent UNLT theories provide. NLGC theory, on the other hand, provides a finite perpendicular diffusion coefficient even if the pure slab model is considered. This is clearly not in agreement with simulations.
The common approach in order to test analytical theories is to compute diffusion parameters via test-particle simulations and compare the numerical results with analytical findings. This can either be done in the late-time limit (see, e.g. Fig. 18 of the current review) or by considering the whole time-dependent diffusion process (see Fig. 17). An alternative test-method was proposed by Qin and Shalchi (2016). In the latter work different approximations were artificially incorporated into the simulations in order to test their validity. Then the corresponding diffusion parameter was computed from the simulations and compared to the result which is not based on this approximation. The tests have been performed for twocomponent turbulence and the different diffusion parameters were plotted versus the slab fraction. The following approximations were tested: 1. First the authors tested the validity of using guiding center coordinates instead of particle coordinates. This basically tests the accuracy of using Eq. (5.22) as equation of motion.
In order to do this the following parameter was defined where κ P ⊥ is the perpendicular diffusion coefficient computed by using particle coordinates and κ GC ⊥ is the diffusion coefficient based on guiding center coordinates. The way how the parameter a 2 is defined, corresponds to . However, according to Fig. 19, we find numerically a 2 ≈ 1 confirming that using guiding center coordinates as well as Eq. where the two diffusion parameters D P ⊥ and D GC ⊥ have been computed via mean square displacements whereas the two parameters κ P ⊥ and κ GC ⊥ were based on the TGK formula. Fig. 19 The upper panel shows perpendicular diffusion coefficients for two-component turbulence versus the slab fraction E slab /E total = δB 2 slab /δB 2 . Here we have used δB/B 0 = 1.0 and R L / = 0.1. The lower panel shows the parameters a 2 , b 2 , c 2 , and d 2 which are defined in Eqs. (10.14)-(10.17). The aim of both panels is to test approximations commonly used in the theory of perpendicular transport. If the corresponding parameter is close to one, this means that the tested approximation works perfectly well. Reprinted with permission from The American Astronomical Society- Qin and Shalchi (2016) Again the ratio b 2 allows us to test the validity of using guiding center coordinates. Furthermore, we can also compare the diffusion parameters κ P ⊥ and D P ⊥ with each other to test the TGK formula. According to Fig. 19, we have b 2 ≈ 1 as well as κ P ⊥ ≈ D P ⊥ meaning that using guiding center coordinates as well as the TGK formula works very well in diffusion theory. 3. As pointed out in the previous sections of this review, a problematic approximation is given by Eq. (5.39). The NLGC theory, for instance, is fully based on this approximation whereas UNLT theories are based on a different approximation (see Eq. (7.3) of this review). In order to check Eq. (5.39), Qin and Shalchi (2016) defined where the diffusion coefficient κ 4 ⊥ is based on approximation (5.39) and κ GC ⊥ is not based on this approximation. As shown in Fig. 19, this approximation is indeed not valid. In particular for slab turbulence it breaks down completely explaining why NLGC theory does not work in this case. Interesting here is that even for the typical scenario of 80% two-dimensional and 20% slab turbulence we have c 2 ≈ 1/3 which means that we need this parameter (in the previous section it was called a 2 ) to balance out the inaccuracy of approximation (7.3). There will be more comments about this matter in Sect. 12.2.1. 4. As the last approximation Qin and Shalchi (2016) tested the validity of the guiding center approximation in the exponential of the characteristic function. In the general case, one needs to include the gyro-rotation as it was, for instance, done in Sect. 7.5. This can be called the second guiding center approximation. Therefore, Qin and Shalchi (2016) defined (10.17) According to Fig. 19, this approximation works well for small and intermediate rigidities. Qin and Shalchi (2016) performed more runs to show that for higher and higher rigidities the second guiding center approximation fails more and more.

Fig. 20
The perpendicular mean free path versus the parallel mean free path for noisy slab turbulence. Shown are test-particle simulations (dots) and diffusive UNLT theory (solid line). The dashed line corresponds to the CLRR limit given by Eq. (6.50) and the dotted line corresponds to the FLRW limit. All these results were obtained for δB/B 0 = 1 and / ⊥ = 0.5 to ensure that the Kubo number is smaller than one. The simulations were originally presented in Heusen and Shalchi (2017) As demonstrated above, using Eq. (5.22) as equation of motion but also the TGK formulation work very well in the theory of perpendicular diffusion. Approximation (5.39), on the other hand, fails in almost all cases. This is in particular the case for pure slab turbulence. Special care is also required for the second guiding center approximation which becomes more and more inaccurate for high particle energies and one needs to include finite Larmor radius corrections in such cases as explained in Sect. 7.5 of this review.

Noisy Slab Turbulence
A simple three-dimensional turbulence model is provided by the noisy slab model given by Eq. (2.25). This model allows us to incorporate transverse complexity by using a minimalistic approach. Therefore, one can use this to explore our understanding of perpendicular transport and to test the validity of analytical theories in the small Kubo number regime. Figure 20 shows simulations for a noisy slab model in comparison with UNLT theory for a 2 = 1. First of all we note that diffusion is restored if the pure slab model is replaced by the corresponding noisy model. Furthermore, we can easily see that diffusive UNLT theory agrees well with simulations. The latter theory predicts a turnover from the CLRR regime, as described by Eq. (6.50), for short parallel mean free paths to the FLRW limit for long parallel mean free paths. In the CLRR regime the perpendicular mean free path is directly proportional to the parallel mean free path. For long parallel mean free paths, on the other hand, the perpendicular mean free path becomes constant corresponding to the FLRW limit. According to diffusive UNLT theory this behavior of the perpendicular mean free path is universal and should be observed for almost all turbulence models.

Noisy Reduced MHD Turbulence
In the previous subsection we have employed the noisy slab model in order to test our understanding of particle transport in the small Kubo number regime. In a similar way we can use Fig. 21 The perpendicular mean free path versus the parallel mean free path for the noisy reduced MHD model. Compared are the results obtained from the simulations of Shalchi and Hussein (2014) (dots), diffusive UNLT theory for a 2 = 1/3 and a 2 = 1 (solid lines), the FLRW limit (dotted line), and the CLRR limit (dashed line). Also shown is the composite formula (grey line) as given by Eq. (9.37). In the case considered here we have used δB/B 0 = 1 and the field line diffusion coefficient is in this case κ F L = 0.23 ⊥ . Reprinted with permission from The American Astronomical Society- Shalchi (2019a) the noisy reduced MHD model as given by Eq. (2.27). The latter model allows us to explore the transport in larger Kubo number turbulence. Figure 21 shows the simulations performed by Shalchi and Hussein (2014). Diffusive UNLT theory does agree with the simulations but only for the right choice of a 2 . The heuristic composite formula given by Eq. (9.37), however, agrees perfectly with the simulations. Like in all considered cases, the theory provides a perpendicular diffusion coefficient which is directly proportional to the parallel diffusion coefficient for small λ and for long parallel mean free paths, the perpendicular mean free path approaches asymptotically the FLRW limit. Interesting here is that for the case of short parallel mean free paths we find agreement with simulations for a 2 = 1/3 whereas for long parallel mean free paths we find agreement if a 2 is close to one. The heuristic approach described in Sect. 9 explains why these values of a 2 are needed.

Goldreich-Sridhar Turbulence
In the recent years spectral tensors based on the Goldreich-Sridhar model became more popular and they were used in test-particle simulations by Sun and Jokipii (2011). The latter authors used a model spectrum which corresponds to the one given by Eq. (2.41) for the case q = 0. Diffusive UNLT theory was combined with this spectral tensor to test the theory via a comparison with the aforementioned simulations (see Shalchi 2013a; . This comparison is visualized in Fig. 22. As in the cases before, we find good agreement between analytical theory and simulations. Interesting here is that we find good agreement for a 2 = 1 meaning that the correction factor a 2 is not needed at all. This was also found for slab as well as noisy slab turbulence but not for two-component and NRMHD turbulence where we need to set a 2 = 1/3. This suggests that UNLT theory works very well for small and intermediate Kubo numbers and only for large Kubo numbers one needs the correction factor a 2 . Furthermore, we can clearly see by considering Fig. 22 that λ ⊥ ∝ λ for small λ whereas the perpendicular mean free path becomes independent of the parallel mean free path for Fig. 22 The perpendicular mean free path λ ⊥ versus the parallel mean free path λ for Goldreich-Sridhar turbulence. Both parameters are normalized with respect to the characteristic turbulence scale . Shown are the mean free paths obtained from the simulations (dots) performed by Sun and Jokipii (2011). Also shown are the mean free paths obtained from diffusive UNLT theory for a 2 = 1 (solid line), the CLRR limit (dashed line), and the FLRW limit (dotted line). Also shown is the composite formula (grey line) as given by Eq. (9.37). In the case considered here we have used δB/B 0 = 1 and the field line diffusion coefficient is in this case κ F L = 0.38 . Reprinted with permission from The American Astronomical Society- Shalchi (2019a) large λ . This behavior was also seen in other cases. Thus there is some universality in the transport of particles across the mean magnetic field (see, e.g., Shalchi 2014;. However, this universal behavior of the perpendicular mean free path requires that all fundamental turbulence scales such as the integral scales and the ultra-scale are finite. In Fig. 22 we have also shown the heuristic composite formula. In this case diffusive UNLT theory is closer to the simulations than the heuristic approach. This is what one would naturally expect since systematic theories should be more accurate than simple heuristic arguments.

Applications
The motion of energetic particles in magnetized plasmas is a fundamental problem and there are applications of transport theory in a variety of scenarios ranging from fusion reactors to the space between galaxies. In the following we discuss some examples to emphasize the importance of theories and analytical results for the perpendicular diffusion coefficient.

Transport of Galactic Protons in the Solar System
Different energetic particles, such as galactic cosmic rays and solar energetic particles propagate through the solar system and interact with the interplanetary magnetic field. Therefore, they experience parallel and perpendicular diffusion. It is one aim of the theory of perpendicular transport to reproduce measurements performed by space probes such as Ulysses. In Burger et al. (2000), for instance, the rigidity dependence of the perpendicular diffusion coefficient was obtained for galactic protons. As also discussed in Bieber et al. (2004), it was found that their perpendicular mean free path is about λ ⊥ ≈ 0.006 − 0.008 AU at higher rigidities. Palmer (1982), on the other hand, focused on the parallel spatial diffusion coefficient. A detailed discussion of this analysis and our ability to reproduce such measurements theoretically can be found in Bieber et al. (1994). A comparison of simulations and observed parallel diffusion coefficients can be found in Hussein and Shalchi (2016). Palmer (1982) also proposed an average perpendicular diffusion coefficient at 1 AU of κ ⊥ c/v ≈ 10 21 cm 2 /s and therefore λ ⊥ ≈ 0.0067 AU. However, it was pointed out that the spread around this average is rather large. Bieber et al. (2004) combined the NLGC integral equation with a two-component turbulence model and they have shown that good agreement with the aforementioned Ulysses measurements can be obtained. For this specific case, diffusive UNLT theory provides very similar results except that, compared to NLGC theory, there is no longer an explicit contribution due to the slab modes because this contribution damps out sub-diffusively as discussed before in this review. Furthermore, Bieber et al. (2004) used a flat spectrum for the energy range of the two-dimensional modes which contradicts the Matthaeus et al. (2007) conditions which require either an increasing spectrum at small wave numbers or a cut-off.
The best approach in order to reproduce solar wind observations is to either combine testparticle simulations or an advanced theory, such as UNLT theory, with an accurate turbulence model and an appropriate value for the parallel mean free path. Besides UNLT theory we also use a much simpler approach. First we assume that the slab/2D model provides a good approximation to solar wind turbulence. If we also take into account that the aforementioned observations were obtained for high rigidities, corresponding to long parallel mean free paths, we can employ Eq. (8.18) to approximate the perpendicular mean free path. This means that in the considered case perpendicular diffusion is predominantly caused by random walking magnetic field lines which are, however, in the Bohmian/non-linear regime. In order to compute the perpendicular mean free path via Eq. (8.18), we need to specify the parameters therein. For the energy range spectral index we set q = 3 (in agreement with Matthaeus et al. 2007) and for the inertial range spectral index we use s = 5/3 as usual. For the magnetic fields we assume δB 2 2D /δB 2 = 0.8 (see, e.g.,  as well as δB 2 /B 2 0 = 0.5 and for the two-dimensional bendover scale we use ⊥ = 0.03 AU. Motivated by test-particle simulations performed for two-component turbulence, we set a 2 = 1/3 (see, e.g., Fig. 17). With those values we find from Eq. (8.18) (11.1) which is in perfect agreement with the Ulysses measurements. This result is compared with diffusive UNLT theory and observations in Fig. 23. There are also results for Jovian electrons obtained in Chenette et al. (1977) which are at smaller rigidities but the corresponding values for the perpendicular mean free path are similar. As pointed out in Bieber et al. (2004), the assumption that the perpendicular diffusion coefficient is much smaller than the parallel diffusion coefficient is not always valid. Dwyer et al. (1997) and Zhang et al. (2003) reported that the ratio of these two diffusion coefficients κ ⊥ /κ can occasionally approach or even exceed unity. Dwyer et al. (1997), for instance, found results for low rigidities usually corresponding to short parallel mean free paths. For two-component turbulence one can then use Eq. (6.45) to compute the ratio κ ⊥ /κ . According to this formula the ratio can be around one or even larger depending on what the magnetic field ratio is. This statement is true regardless whether we have a 2 = 1/3 (CLRR diffusion) or a 2 = 1 (fluid limit).

Fig. 23
The perpendicular mean free path λ ⊥ of particles propagating through the solar system. Shown are Ulysses measurements of galactic protons (Burger et al. 2000, dots), the value for the perpendicular mean free path as suggested by Palmer (1982, horizontal gray line), and observational determinations of Jovian electrons (Chenette et al. 1977, square). For the theoretical results we have shown the diffusive UNLT result represented by the black solid line and the (non-linear) FLRW limit as given by Eq. (11.1) for a 2 = 1/3 and a 2 = 1 (dashed lines). The dotted line represents the used parallel mean free path as given by Eq. (4.34). The UNLT result was obtained by combining Eq. (4.34) for the parallel mean free path with the simple analytical form given by Eq. (8.33)

Particle Acceleration at Perpendicular Interplanetary Shock Waves
Different shock waves are believed to accelerate particles up to high energies through the mechanism of diffusive shock acceleration. The latter process is a first-order Fermi acceleration process (see, e.g., Fermi 1949; Axford et al. 1977;Krymsky 1977;Bell 1978a,b;Blandford and Ostriker 1978;Drury 1983;Blandford and Eichler 1987;Malkov and Drury 2001). In the following we review the importance of the perpendicular diffusion coefficient in the problem of particle acceleration at perpendicular interplanetary shock waves. Examples of such interplanetary shocks are those driven by coronal mass ejections (CME's). The theory of particle acceleration at quasi-parallel shocks appears to be reasonably well understood, and has been applied to solar energetic particle and energetic storm particle events by  as well as Li et al. (2003) and (2005).  have developed a theory for describing particle acceleration at quasi-perpendicular shocks. This approach is discussed in the following. Figure 24 shows parallel and perpendicular CME driven shock fronts.
Diffusive shock acceleration is described via Parker's transport equation (see Eq. (1.1) of this review) or extensions thereof. For a one-dimensional scenario Eq. (1.1) reduces to where f (x, p, t) is the cosmic ray distribution as a function of position x along the shock normal, momentum p measured in the rest-frame of the streaming plasma, and time t . The term with the x-derivative on the left hand side of Eq. (11.2) describes convection whereas the terms on the right hand side correspond to diffusion, energy change, and sources/losses.   The parameter κ ⊥ corresponds to the diffusion coefficient in the shock propagation direction. In this form of the transport equation we have neglected stochastic acceleration and other effects such as radiative losses. Equation (11.2) has a steady-state solution of the form f (p) ∝ p −3σ/(σ −1) where we have used the compression ratio σ = u u /u d . For a compression ratio of σ = 4, for instance, this turns into f (p) ∝ p 4 . The time it needs to accelerate a particle from the injection momentum p inj to the maximum momentum p max is, in general, given by (see Drury 1983;Webb et al. 1995 where we have used upstream and downstream diffusion coefficients κ u and κ d as well as upstream and downstream velocities u u and u d . This relation is valid for non-relativistic as well as relativistic particles. Figure 25 shows a simple sketch of upstream and downstream regions, respectively. There are different ways of simplifying relation (11.3). Sometimes it is assumed that the downstream diffusion coefficient is much smaller than the upstream one. Using this and κ ⊥ = κ u allows us to simplify Eq. (11.3) to where we have used the shock velocity u sh = u u . Sometimes a more complicated approach is used for replacing the diffusion coefficient in Eq. (11.3). In a perpendicular magnetohydrodynamic shock one expects a jump in the magnetic field. 20 Keeping in mind that the perpendicular diffusion coefficient depends on magnetic fields, this could lead to a jump of the diffusion coefficient. One can, for instance, assume that κ d = κ u /σ (see, e.g., Kang et al. 2009). Since the compression ratio is usually assumed to be around σ ≈ 4, this means that the downstream diffusion coefficient is indeed smaller than the upstream coefficient as Fig. 25 Sketch showing the shock with its upstream and downstream regions in the rest frame of the shock assumed above. However, if we employ the relation κ d = κ u /σ , we find an increase of the time in Eq. (11.4) by a factor two (see, e.g., Ferrand et al. 2014 for more details and different relations between perpendicular diffusion coefficients and acceleration times). The aim of this subsection is to compute the maximum energy a charged particle can get due to shock acceleration. This energy will be calculated as a function of time requiring to not just specify the analytical form of κ ⊥ (p) in Eq. (11.4), but also the time-dependence of other parameters such as magnetic fields. In the following we employ the approach presented in  which is based on the ideas of . First we rewrite Eq. (11.4) by replacing the diffusion coefficient by the mean free path (κ ⊥ = vλ ⊥ /3) to find (11.5) In the last step we have used the relativistic momentum p = mvγ as well as dp = mγ 3 dv with the Lorentz-factor γ . We can evaluate the expression given by Eq. (11.5) for the velocity independent perpendicular mean free path derived in Sect. 8 (see Eq. (8.18) of this review). We obtain (11.6) The integral therein can be expressed by inverse hyperbolic tangent functions so that (11.7) After using the relation tanh(x − y) = tanh(x) − tanh(y) 1 − tanh(x) tanh(y) (11.8) we find where we used the dimensionless time (11.10) Since the injection velocity v inj is much smaller than the maximal velocity v max and the speed of light c (see, e.g., Fig. 7 of ), we can approximate Eq. (11.9) by v max c = tanh(τ ). (11.11) For the maximum particle momentum, on the other hand, we can use and for the maximum (relativistic) kinetic energy E max , one can easily derive (11.13) Equations (11.11)-(11.13) allow us to compute maximum velocity, momentum, and energy. In order to compute the maximum energy that a particle can get due to its interaction with an interplanetary shock wave, we also have to specify the shock position r(t) and the shock velocityṙ(t). After the shock has swept up sufficient mass, a self-similar Sedov-Taylor blast wave solution (see Taylor 1950;Sedov 1959) is assumed to begin. The Sedov-Taylor solution is where E blast is the total blast wave energy (typically 5 · 10 32 ergs ≈ 5 · 10 25 kg m 2 /s 2 for a CME driven blast wave) and ρ 0 is a characteristic mass density. From Eq. (11.14) we derive for the shock velocityṙ The mass density therein can be estimated by using ρ 0 = m p · n 0 where we have used the proton mass (m p ≈ 1.67 · 10 −27 kg) and the particle density n 0 ≈ 10 6 /m 3 . For the velocity we then findṙ c ≈ 4 · 10 −3 1 AU r  where we have also assumed that 80% of the magnetic energy is in the two-dimensional modes as usual. For the mean magnetic field we assume a Parker spiral field (Parker 1958) where ω 0 is the solar angular velocity, v sw is the solar wind speed, and B 0,co is the magnetic field at the heliocentric corotation radius r = r co . Typically, we have ω 0 = 2π/(25.4 days), v sw = 400 km/s, B 0,co = 1.83 · 10 −6 T, and r co = 46 · 10 −3 AU (we assumed that the corotation radius is approximately 10 times larger than the radius of the Sun) leading to v sw /ω 0 ≈ 1 AU. By employing these values Eq. (11.25) becomes B 0 (r) = B 0,co r co r 2 ≈ 1830 nT 46 · 10 −3 AU r 2 1 + r 1 AU 2 .
(11.26) By combining Eqs. (11.24) and (11.26) we find, The parameter λ ⊥,0 is given by Eq. (11.21) and for the compression ratio we use σ = 3.7. By combining Eqs. (11.22)-(11.27) with these values we can derive τ = 0.7 1 + 1 AU r 2 . (11.28) The parameter τ can easily be combined with Eq. (11.13) to compute the maximum kinetic energy. The maximum energy a proton can get by interacting with a CME driven shock is illustrated in Fig. 26 together with the curve obtained in . The results shown there were obtained by combining Eqs. (11.13) and (11.28). For small heliocentric distances we have E max ≈ 4.5 GeV and for r ≈ 1 AU we find E max ≈ 500 MeV. These values are significantly higher compared to previous results due to the updated turbulence and particle parameters such as the two-dimensional bendover scale ⊥ . Since the exact values for those parameters are not known in the heliosphere, one could easily obtain different particle energies. Then, on the other hand, the used perpendicular diffusion coefficient is in coincidence with observations (compare Eq. (11.21) with Eq. (11.1)). However, we have also plotted the result for a 2 = 1 and obtained a significantly lower energy. One could argue that for the considered energy the perpendicular diffusion coefficient is in the FLRW limit and, thus, a 2 = 1 is indeed correct. The result obtained here is consistent with that derived originally in  who also obtained GeV energies for protons during the early phase of the shock. Both here and in , the reason for the high energies is the combination of the large magnetic field strength close to the Sun and the speed and strength of the shock. The decay in the maximum particle energy accelerated at the shock is a consequence of the magnetic field strength weakening with increasing heliocentric distance, as well as the shock slowing down.
Similar work was done by Dosch and Shalchi (2010) where perpendicular diffusion coefficients based on different spectra in the energy range were used and it was explored how this impacts the maximum energy at a perpendicular interplanetary shock. It was shown that the large turbulence scales have a strong impact on the maximum energy. The largest energy was clearly found for the case where the spectrum increases in the energy range and where the corresponding spectral index agrees with the Matthaeus et al. (2007) conditions. In this case the perpendicular diffusion coefficient at high energies corresponds to the non-linear field line random walk limit where the perpendicular mean free path does not depend on particle rigidity or energy. This is exactly the case considered above. For flat or even decreasing turbulence spectra, the perpendicular mean free path depends on rigidity leading to much smaller energies due to shock acceleration (see, e.g., Figs. 2 and 3 of Dosch and Shalchi (2010)). Particle acceleration at an oblique CME-driven shock was explored in Li et al. (2012).

Cosmic Ray Propagation in the Nearby Starburst Galaxy NGC 253
In Buffie et al. (2013) the propagation of cosmic rays in an external galaxy was explored theoretically. Radio continuum observations of the nearby starburst galaxy NGC 253 can be used to measure distribution and transport of cosmic ray electrons (see Heesen et al. 2009aHeesen et al. , 2011. The current understanding is that these electrons are accelerated in the disk by shock waves induced by supernova explosions (see Reynolds et al. 2012). The particles are then transported away by either convection in a galactic wind or spatial diffusion. Radio continuum studies can be used to measure the length scale of cosmic ray transport. The time scale of cosmic ray electrons is determined by loss processes by which the electrons are losing their energy, such as synchrotron and inverse Compton radiation. The diffusion coefficient can then be estimated via (11.29) where L diff is the diffusion length scale and τ is the electron lifetime. Depending on the magnetic field structure, the observed diffusion coefficient is either along the magnetic field or perpendicular with respect to it. Ideally, one can measure the distance to the star-formation sites, where the cosmic rays are accelerated and injected into the interstellar medium, to obtain the transport length scale. This is for instance the case when observing galaxies in so-called edge-on geometry, where the observer is located in or close to the disk plane of the galaxy. Cosmic ray acceleration in supernova remnants is confined to the relatively thin disk plane, where the formation of massive stars happens. The geometry in this case is simple: cosmic rays are transported away from the star formation sites over their lifetime and their transport length scale is equal to the vertical electron scale height. Typical scale heights of galaxies are 1.8 kpc 21 at observing wavelengths of 6 cm (see Krause 2009). Typical cosmic ray electron lifetimes in the interstellar medium are between 1 and 10 Myr. This basically determines the typical parallel diffusion coefficient κ in the interstellar medium, where it is assumed that the halo magnetic field is opening up from the disk parallel to a vertical direction with increasing distance from the disk (see Beck 2012). Therefore, vertical diffusion is predominantly along the magnetic field lines hence allowing us to measure κ . Heesen et al. (2009b) confirmed this scenario for NGC 253. For the parallel diffusion coefficient one can assume κ = 1.0 × 10 29 cm 2 s −1 at a magnetic rigidity of 3 × 10 12 Volts, equivalent to an electron energy of 3 GeV (Heesen et al. 2009a).
To measure the perpendicular diffusion coefficient is more difficult because it requires to have variations of the cosmic ray distribution perpendicular to the magnetic field orientation. The observed radio continuum emission is mostly smooth in the disk-halo interface and filamentary structures are rare. Magnetohydrodynamic simulations suggest that the disk halo interface is dominated by filamentary magnetic fields (see Breitschwerdt et al. 2012), but line-of-sight confusion and limited spatial resolution hampers their detection. An exception are starburst galaxies such as NGC 253, where the spatially concentrated star formation activity results in exceptionally high radio continuum surface brightness, allowing us to employ high spatial resolution (≈ 100 pc). The spatially concentrated star formation can result into outflows of hot X-ray emitting gas in a galactic wind. The magnetic field is then concentrated and amplified by expansion of the hot gas until a pressure equilibrium is reached. This very specific geometry allowed Heesen et al. (2011) to determine the perpendicular diffusion coefficient across the magnetic field in the walls of the nuclear outflow cone in NGC 253. From the observations of particles with a magnetic rigidity of 3 × 10 12 Volts, they found a perpendicular diffusion coefficient of κ ⊥ = (2.6 ± 0.6) × 10 28 cm 2 s −1 . It is important to measure the diffusion coefficients at roughly the same electron energy. Buffie et al. (2013) have reduced the measured perpendicular diffusion coefficient of Heesen et al. (2009a) to account for a possible contribution of convection.
In order to reproduce the observed perpendicular diffusion coefficient we employ diffusive UNLT theory as given by Eq. (6.25). The theory requires the specification of the spectral tensor describing magnetic turbulence in interstellar media. In the following we employ the tensor based on the Goldreich-Sridhar model as given by Eq. (2.39) with (2.40). Using this in Eq. (6.25) and employing the integral transformations x = k and y = 1/( k ⊥ ) yields for the ratio of perpendicular and parallel mean free paths × 2x 2 y 2 + 1 x 2 y 2 + 1 y 7/3 x 2 y 4 λ /λ ⊥ + y 2 + 4λ λ ⊥ /(3 ) 2 .
(11.30) Motivated by Fig. 22 we set a 2 = 1 because this value gives good agreement between simulations and UNLT theory for Goldreich-Sridhar turbulence. 21 In the present section we measure distances in parsecs (pc). Approximately we have 1 pc ≈ 3 · 10 16 m ≈ 3.3 ly.

Fig. 27
The perpendicular diffusion coefficients of cosmic rays in the nearby starburst galaxy NGC 253. We show the theoretical perpendicular diffusion coefficient for different correlation scales . For the parallel diffusion coefficient we have used κ = 1.0 × 10 29 cm 2 /s as proposed by Heesen et al. (2009a). The theoretical values were calculated for δB/B 0 = 1 (dashed line) and δB/B 0 = 2 (solid line). The dot represents the perpendicular diffusion coefficient from the observations (see Heesen et al. 2011), where we have κ ⊥ = 2.6 × 10 28 cm 2 /s at = 50 pc Obviously, the perpendicular mean free path depends on the magnetic field ratio δB/B 0 and the length scale . The values for these parameters are discussed in the following. According to Heesen et al. (2011) the ordered magnetic field strength is B 0 = 21 µG whereas the turbulent field is δB = 41 µG. Therefore, the magnetic field ratio is δB/B 0 ≈ 2 in the case considered here. It should be noted, however, that this is an upper limit as the spatial resolution that was available for the polarization measurements may be not be enough to resolve the filamentary magnetic fields. Thus, Buffie et al. (2013) have also studied the implication of a lower value of δB/B 0 = 1 to explore how it changes the theoretical perpendicular diffusion coefficient.
More difficult to estimate is the scale . The width of the cone walls in which the magnetic fields are confined is equal or less than 40 pc (see Heesen et al. 2011). This suggests that the upper value cannot be larger than ≈ 50-100 pc. However the value for is very uncertain and, therefore, Buffie et al. (2013) have computed the perpendicular diffusion coefficient for a whole range of correlation lengths. Beck (2007), for instance, suggested that the largest scales of turbulence are in the order of 10-100 pc. Such largest scales, however, are not necessarily equal to the turbulence correlation scale. Actually they can be seen as the maximum of the scale . Sometimes it is assumed that the correlation length of interstellar turbulence is 1 pc. Therefore we compute the perpendicular diffusion coefficient for 1 pc ≤ ≤ 1000 pc to explore the values of which lead to agreement between theory and observations. Turbulence spectra in the local interstellar medium are discussed in Armstrong et al. (1995). The values found there for the scale are consistent with the values discussed above.
In Fig. 27 we show the perpendicular diffusion coefficient versus the correlation scale . The theoretical results are based on the numerical solution of Eq. (11.30). Mean free paths and diffusion coefficients are related to each other via Eq. (5.37). In those relations we replace the particle speed v by the speed of light c since we deal with relativistic particles. The observational value of the perpendicular diffusion coefficient is κ ⊥ = (2.6 ± 0.6) × 10 28 cm 2 s −1 for a correlation length of approximately = 50 pc. The theoretical perpendicular diffusion coefficient is calculated for different values of the correlation length. For a correlation length of 50 pc we find a κ ⊥ which is very close to the observations depending on the value of δB/B 0 . As shown in Fig. 27 there is very good agreement between theory and observations. Of course one has to keep in mind that in in-terstellar media the exact values of turbulence parameters and diffusion coefficients are not available. Thus, it is easy to reproduce interstellar measurements of cosmic rays compared to the interplanetary scenarios. Nevertheless, the nice agreement between theoretical and observational diffusion parameters is another confirmation of our improved understanding of perpendicular diffusion.

Cosmic Ray Acceleration at Perpendicular Shocks in Supernova Remnants
It is believed that galactic cosmic rays get their energy due to diffusive shock acceleration at supernova remnants. Diffusive shock acceleration at an interstellar shock wave is similar compared to the interplanetary scenario discussed in Sect. 11.2. It was shown in Ferrand et al. (2014) how advanced analytical results for the perpendicular diffusion coefficient are useful in the theory of diffusive shock acceleration at supernova shocks. In the current section we consider a perpendicular interstellar shock and, thus, the diffusion coefficient in Eq. (11.2) is again the perpendicular diffusion coefficient. Furthermore, the acceleration time is still given by Eq. (11.3).
In Ferrand et al. (2014) the so-called MARCOS code, developed and explained in Ferrand et al. (2008), was used in order to perform hydro-kinetic simulations of the coupled shock + particles system. The code jointly solves the three conservation equations that describe the thermal fluid and the transport equation (11.2) that describes the energetic particles. The hydro equations have been solved using an explicit Godunov (1959) scheme and the particle transport equation, written for g = p 4 f as a function of y = ln (p), is solved by using a semi-implicit Crank-Nicolson scheme (see, e.g., Press et al. 2007). Furthermore, the code includes the effect of Alfvénic heating on the flow and of Alfvénic drift on the particles (see, e.g., Kang 2013). The geometry is one-dimensional in space and one-dimensional spherically symmetric in momentum. The shock was generated by a constant velocity piston and position and properties of the shock front are diagnosed at each hydro time step. There, particles are injected continuously in time from the thermal pool. For simplicity we consider that this occurs in a fixed momentum bin p inj and at a constant rate η expressed as a fraction of the flux crossing the shock. The diffusive shock acceleration solver then spreads the distribution of particles over space and momentum. The pressure of particles is computed from their momentum distribution in each space cell and is included in the hydro equations at each time step leading to modified shocks and spectra.
In previous explorations of diffusive shock acceleration at supernova remnants isotropic Bohm diffusion has been assumed. 22 For parallel diffusion the Bohm limit was derived systematically in Shalchi (2009b) where it was shown that in this case the parallel mean free path is given by λ Bohm = R L B 0 /δB meaning that where we replaced the unperturbed Larmor radius R L by using Eq. (4.7). In Shalchi (2009b) the Bohm diffusion coefficient was obtained in a highly non-linear regime where the turbulent magnetic field is very strong. The Bohm limit, as given by Eq. (11.31), was confirmed numerically by Hussein and Shalchi (2014a). For perpendicular diffusion, on the other hand, there can be a regime where κ ⊥ ≈ κ but this requires isotropic turbulence with a very strong turbulent magnetic field (see, e.g., Shalchi and Dosch 2009;Hussein and Shalchi 2014a). However, in the case of perpendicular diffusion, one still expects to find anisotropic scattering in the high-energy limit. In this case the symmetry is broken because if the particles move very fast, they only experience the mean magnetic field even if the turbulent field is strong (see again Shalchi and Dosch 2009). Therefore, the Bohm limit is not a good approximation for particle diffusion in the general case. Below we follow Ferrand et al. (2014) and construct a more advanced and realistic model for the perpendicular diffusion coefficient based on the previous sections of this review. First we replace v = c p/(m p c) 1 + (p/(m p c)) 2 (11.32) in Eq. (11.31) to write whereB 0 is the background field in microGauss. Alternatively, some authors (see, e.g., Kang et al. 2009) used a power-law model of the form κ (p) = κ * p α (11.37) but this was mostly done in the context of parallel diffusion. However, below we shall link the perpendicular diffusion coefficient to this form of the parallel diffusion coefficient.
We have shown in previous sections of this review that the perpendicular mean free path increases with increasing particle momentum and that for high momenta the perpendicular mean free path becomes independent of momentum. This behavior of the perpendicular diffusion coefficient is universal meaning that it does not depend on turbulence parameters and can also be seen in test-particle simulations (see, for instance, Figs. 20 and 22 of this is reduced, so is the spectrum curvature. For a very low p c , the spectra are very steep. As can be seen in Fig. 28, as p c is lowered, the spectrum first sees its concavity being restricted to a smaller range of momenta while still reaching similar values at p max , leading to highly modified spectra for that particular parameter range. We note that the shape of the cosmic ray spectrum will be reflected in the shape of the spectrum of non-thermal photons radiated via the production and decay of pions for energetic protons. Recent γ -ray observations of supernova remnants do not support the existence of concave spectra of the underlying cosmic ray population, on the contrary they require rather steep spectra (see Caprioli 2011). The approach reviewed in this subsection and originally presented in Ferrand et al. (2014), naturally produces such spectra, even though particles may be efficiently accelerated and take a large part of the shock energy.

Summary and Conclusion
In the context of space science and astrophysics the diffusion of energetic particles across the mean magnetic field was first discussed in the pioneering work of Jokipii (1966). In the latter paper the perpendicular diffusion coefficient was computed by employing quasi-linear theory which has to be understood as a first-order perturbation theory. It is well-known that quasi-linear theory has its limitations. This comes due to the fact that real orbits are not unperturbed and that parallel diffusion can have a strong influence on perpendicular diffusion as described by compound sub-diffusion models (see, e.g., Webb et al. 2006 for a comprehensive analytical description of compound sub-diffusion). Furthermore, there is the large Kubo number regime where the behavior of the magnetic field line random walk itself is highly non-linear (see, e.g., Matthaeus et al. 1995). From a more modern point of view, quasi-linear theory should provide an accurate description of perpendicular transport for small Kubo number turbulence and very long parallel mean free paths corresponding to high rigidities. In general the theory does not work. During the past two decades, however, progress was achieved mainly because of two reasons: 1. Test-particle simulations were developed often using super-computing. Such simulations allow for an accurate description of particle transport but they are, of course, based on numerical approximations (e.g., finite wave number grid-size and number of time-steps) but they do not require ad-hoc assumptions necessary in analytical treatments. Of course, one still needs to specify the properties of magnetic turbulence such as spectral anisotropy and the form of the spectrum. Apart from that, simulations allow to compute particle diffusion parameters accurately and, thus, simulations can be used in order to test our understanding of the particle motion through magnetic turbulence. 2. Non-linear transport theories were developed based on Corrsin's independence hypothesis. The problem in analytical treatments of the transport is the emergence of higher order correlations involving particle positions, velocities, as well as magnetic fields. In reality all those quantities are somehow correlated. Corrsin's approximation allows for a strong simplification of such correlations by writing them as a product of correlations only containing magnetic fields and correlations containing only particle properties. Based on this approximation advanced theories such as the non-linear guiding center theory or the unified non-linear transport theory were developed. Whereas the former theory can be seen as an important step in this field, the latter theory is superior because it contains previously known results and agrees with simulations for a wider range of parameters. Right: spectrum of the particles (g = p 4 f ) as a function of momentum at different times. For the diffusion coefficient Eq. (11.41) was employed where we have set α = 1 and ν = 1. Furthermore, p c is increasing from the bottom to top such that p c = 10 0 , 10 2 , 10 4 , 10 6 , and 10 8 in m p c units. Therefore, the topmost plot corresponds to the Bohm case usually used in this type of work. The plot at the bottom shows the result for a more realistic perpendicular diffusion coefficient. Reprinted with permission from The American Astronomical Society- Ferrand et al. (2014) Non-linear theories can be tested by comparing their results with test-particle simulations. Therefore, the combination of both, simulations and non-linear tools, allow for a realistic description of particle transport. A major step forward in the analytical description of perpendicular transport was made due to the development of the non-linear guiding center (NLGC) theory by . Compared to other theories developed in the past the latter theory agreed with simulations performed for slab/2D turbulence. However, this agreement was only achieved with the help of a correction factor a 2 = 1/3. Furthermore, NLGC does not work for slab turbulence and three-dimensional turbulence with small and intermediate Kubo numbers. The reason for this is approximation (5.39) as shown in Shalchi (2005b). In order to fix the problem related to slab and small Kubo number turbulence, the unified non-linear transport (UNLT) theory was developed in Shalchi (2010). A time-dependent version of UNLT theory was presented in Shalchi (2017) as well as Lasuik and Shalchi (2017) which is, compared to previous descriptions of perpendicular transport, no longer based on the diffusion approximation. One of the most important features of diffusive and time-dependent UNLT theories is that they contain several previously derived theories and equations as special limits. Those are listed in the following: 1. In the limit of long parallel mean free paths, the unified non-linear transport equation turns into the field line random walk limit and one obtains automatically the non-linear field line theory of Matthaeus et al. (1995). Therefore, the latter theory is contained in the non-linear integral equation provided by diffusive UNLT theory. 2. For magnetostatic slab turbulence time-dependent UNLT theory provides compound subdiffusion. 3. For long parallel mean free paths and small Kubo numbers, UNLT theory is identical compared to quasi-linear theory. 4. For short parallel mean free paths and small Kubo numbers, one can derive a Rechester and Rosenbluth (1978) type of scaling from UNLT theory. Since this process does not require to incorporate collisions, this limit can be called the collisionless Rechester & Rosenbluth (CLRR) scaling. This result is interesting because it provides a ratio κ ⊥ /κ which is small and does not depend on particle energy. 5. For short parallel mean free paths and large Kubo numbers one finds κ ⊥ /κ ≈ a 2 δB 2 x /B 2 0 . For a = 1 this corresponds to the so-called fluid limit whereas for a = L U / ⊥ we obtain the CLRR limit for large Kubo numbers as explained in Sect. 9 of this review.
A strength of time-dependent UNLT theory is that is explains that normal diffusion is restored in the late-time limit as soon as there is transverse complexity of the turbulence. Furthermore, UNLT theory agrees well with most test-particle simulations as summarized in Table 9. This is in particular the case for slab turbulence and three-dimensional turbulence with small Kubo numbers. For instance, UNLT theory agrees almost perfectly with simulations performed for noisy slab and Goldreich-Sridhar turbulence. However, for two-component turbulence, often used to approximate solar wind turbulence, and three-dimensional turbulence with large Kubo numbers, the theory still requires the inclusion of the correction factor a 2 = 1/3.
In order to develop an understanding of the physics of perpendicular transport based on simple and fundamental arguments a heuristic approach was presented in Shalchi (2019a). The latter approach, which is based on three rules, states that perpendicular transport is only controlled by three effects, namely field line random walk, parallel transport, and transverse complexity of magnetic turbulence. As a consequence one is able to obtain different cases of perpendicular transport and one finds different routes to normal diffusion as summarized in Tables 6 and 7 of this review. However, the heuristic approach does not only lead to an improved understanding of perpendicular transport, it also explains the reason and value of the parameter a 2 . Previous approaches such as NLGC and UNLT theories provides the so-called fluid limit (see, e.g., Eq. (6.45) of this review) in the limit of very short parallel mean free paths. However, the heuristic approach predicts that in this limit one should find the CLRR limit. The factor a 2 is able to balance out the difference between fluid and CLRR limits which is, in the large Kubo number limit, given by Eq. (9.40) leading to a 2 = 1/3 in some cases.
As also pointed out in this review, modern theories for perpendicular transport can explain perpendicular diffusion coefficients obtained from observations. This is true for interplanetary as well as interstellar scenarios. Furthermore, if such modern results are used in studies of diffusive shock acceleration interesting new results can be obtained motivating to use these analytical forms in further applications in astrophysics and space science.

Outlook
Undoubtedly tremendous progress has been achieved over the past view years in the theory of perpendicular transport. This was partly due to advanced computer simulations but mostly due to the development of systematic non-linear transport theories. However, there are still unsolved puzzles motivating future work in this field.

Further Improvement of UNLT Theories
As shown in this review, time-dependent UNLT theory works very well in most cases. There is a nice agreement between simulations, heuristic considerations, and systematic theories. However, there is one exception and that is perpendicular transport for the case of short parallel mean free paths, corresponding to strong pitch-angle scattering, and large Kubo numbers. Strictly speaking, UNLT theories provide the fluid limit in this regime whereas the heuristic approach discussed in Sect. 9 provides a CLRR type of transport. This failure of the theory can easily balanced out via the correction parameter a 2 . Ideally, however, this should not be required. Therefore, one has to achieve a further improvement of the theory. Time-dependent UNLT theory is based on the approximation v z (t)v z (0)e ik·x ≈ v z (t)v z (0)e ik z e ik ⊥ ·x ⊥ . (12.1) In the limit of small Kubo numbers this approximation becomes exact whereas for large Kubo numbers this approximation fails. In the general case there can be a strong correlation between velocities and the perpendicular position of the particle. It was already shown in Shalchi (2005b) that the assumption of uncorrelated positions and velocities omits compound sub-diffusion. Since CLRR diffusion comes after the sub-diffusive regime, theories based on approximation (12.1) are incomplete and do not provide sub-diffusion and CLRR diffusion for large Kubo numbers. The next step in the theoretical description of perpendicular transport is, therefore, the formulation of a theory which does not rely on this approximation regardless of what the Kubo number is. Furthermore, it has to be subject of future work to determine the exact value of the parameter a 2 for a variety of turbulence models and parameter regimes. Some work has already been done by Qin and Zhang (2014) providing some more detailed insight about the value of the parameter a 2 . However, the latter work was performed for a constant spectrum in the energy range of the two-dimensional modes and was based on the extended non-linear guiding center theory. Future work should explore the value of a 2 for more modern spectra and should also be based on time-dependent UNLT theory. This will help us to ensure that we have the correct understanding of the physical meaning of this parameter.

Dropping the Corrsin Approximation in Transport Theory
Theories developed in order to describe the random walk of magnetic field lines as well as particle transport are based on Corrsin's independence hypothesis. Although it is usually believed that in the context of space and astrophysical plasmas this approximation works well, one should try to replace the Corrsin approximation by a more reliable and systematic approach. Although some work has been performed in the past such as the development and application of percolation theory (see, e.g., Gruzinov et al. 1990;Isichenko 1991) or the decorrelation trajectory method (see, e.g., Vlad et al. 1998Vlad et al. , 2004Balescu 2005;Negrea et al. 2007), it is a difficult task to develop pure analytical theories without employing Corrsin's independence hypothesis.

The Strong and Compressible Turbulence Case
Most analytical theories, but also test-particle simulations, were developed with the incompressible case δB z = 0 in mind. Even if one considers the compressible case, theories of field line random walk and perpendicular particle transport such as UNLT theories should still be valid as long as δB z < B 0 meaning that the parallel component of the turbulent field is weaker than the mean field. There could be cases, however, such as supernova shock waves, where one can find effects such as wave amplification leading to a strong turbulent field. Therefore, it should be useful to develop analytical theories for such cases. Some work has been done in the past (see, e.g., Shalchi and Dosch 2009;Plotnikov et al. 2011) with some promising results. More systematic theories for this more general case should be developed in the future.

Exploration of Exotic Regimes
For most turbulence models discussed in the literature, the tools and results presented in this review should be valid. However, there could be extreme cases where this is no longer true. In particular parallel and perpendicular transport in pure two-dimensional turbulence are not fully understood and one would expect a very different behavior in this case. It could easily be that Eq. (7.8) representing time-dependent UNLT theory is still valid for two-dimensional turbulence. However, parallel transport could be very different in pure two-dimensional turbulence as, for instance, shown analytically in  and numerically in Arendt and Shalchi (2018). Therefore, the function ξ(k , t), which enters Eq. (7.8), would be different compared to the one used in this review.

Detailed Measurements of Magnetic Turbulence
All analytical theories for particle transport require to specify the spectral tensor of the magnetic fluctuations. Such tensors contain several parameters which need to be specified regardless what the considered turbulence model is. Turbulence theorists, but also observers, often focus on the inertial range corresponding to scales smaller than the bendover scale. However, perpendicular transport of energetic particles is mainly controlled by the large scales of the energy range. This is very different compared to parallel diffusion where one finds gyro-resonant interactions.
Future work should focus on the measurement of magnetic fluctuations at large scales. In particular the Kubo number is important since this number determines together with the parallel mean free path how particles experience perpendicular diffusion and what value the perpendicular diffusion coefficient has. Therefore, at least the Kubo number should be determined with high accuracy via solar wind measurements.