Integral Equations for Problems on Wave Propagation in Near-Earth Plasma

: We consider the solutions of two integrodifferential equations in this work. These equations describe the ultra-low frequency waves in the dipol-like model of the magnetosphere in the gyrokinetic framework. The ﬁrst one is reduced to the homogeneous, second kind Fredholm equation. This equation describes the structure of the parallel component of the magnetic ﬁeld of drift-compression waves along the Earth’s magnetic ﬁeld. The second equation is reduced to the inhomogeneous, second kind Fredholm equation. This equation describes the ﬁeld-aligned structure of the parallel electric ﬁeld potential of Alfvén waves. Both integral equations are solved numerically.


Introduction
The logic of our (material) world is described by differential equations. Indeed, the movement of the single-particle is governed by Newton's second law: where r is the particle's radius-vector, v = d r/dt is the particle's velocity, and F is the force. When particles are numerous, this equation as itself becomes useless, but the medium can be described by the kinetic approach with the help of the Boltzmann equation. In particular, in the collisionless plasma (such as the Earth's magnetosphere plasma) this equation takes the following form: where f is the distribution function (showing how many particles fit in the six-dimensional cube with the volume d rd v in the phase space). Again, this is a differential equation as well.
Knowing the distribution function, one can calculate the charge and current density, as follows: ( Here, the summation is applied to all sort of the charged particles in plasma, such as electrons, protons, etc. Then, one can substitute n and j into yet another set of the differential equations-Maxwell's equations-and obtain the spatio-temporal behavior of the electric and magnetic fields E, B, that is, obtain all possible knowledge of the plasma system. Thus, the world is governed by differential equations! However, this program works only in theory. In practice, to solve the set of Newton, Boltzmann, and Maxwell's equations for the non-trivial system (such as the Earth's magnetosphere) one faces insurmountable difficulties, and one has to apply various approximations. One of them works in the solution of Equation (3) in plasma physics. It is called the integration along the particle trajectories. Let us remind that in the dipole-like magnetic field (such as the Earth's magnetic field), the motion of a charged particle can be divided into three types: moving around a magnetic field line with the cyclotron frequency or gyrofrequency ω c ; moving between "mirror points" or bounce-motion with bounce frequency ω b ; and moving around the Earth with drift frequency ω d [1]. As a result, the Boltzmann equation is reduced to an integral equation, which is much easier to solve than the original Equation (3). Therefore, the differential equations give a correct and complete description of the world but in some situations they are useless, as they are almost impossible to solve, and one has to appeal for the help of the integral equation, which is not so exact but truly works. It is the situation that those who study plasma have to deal with.
Then, many phenomena in the magnetosphere flow rather slow, compared with the movement around a field line. An example is the so-called ultra-low-frequency (ULF) waves. These are the plasma perturbations with the wavelength of the same order as the field line itself (several of Earth's radii). These waves are widely observed in the magnetosphere and are responsible for the auroral intensifications in the polar regions [2][3][4] and transport and acceleration of radiation belt electrons [5][6][7], and thus represent an important factor of space weather. For studying such slow phenomena, one has to average the particle motion over the cyclotron frequency, which reduces the wave equation to an even simpler form-but they still remain the integral equations! Such formalism is called gyrokinetics and was developed in papers [8][9][10]. This paper aims to review some results in the studies of the ULF waves in the Earth's magnetosphere obtained with the integral equations [11,12].

Coordinate System and Governing Equations
To construct an adequate model of wave propagation in the inner magnetosphere, we use a dipole-like model of the magnetosphere. This model takes into account the curvature of the magnetic field lines and the inhomogeneity of the magnetic field. In addition, in such a model, we can take into account the gradients of the background plasma across the magnetic shells. We use an orthogonal curvilinear coordinate system, where the coordinates x 1 and x 2 represent the radial and azimuthal coordinates, and the x 3 coordinate marks a point on a particular field line [13]. The increment of physical length along a field line is dl = √ g 3 dx 3 , where g 3 is the metric tensor component. Similarly, for the transverse direction, we have dl 1 = √ g 1 dx 1 , and dl 2 = √ g 2 dx 2 .
The equilibrium condition of plasma can be written as follows: where B and P are the equilibrium magnetic field and plasma pressure, respectively, R is the field line curvature radius, and β = 8πP/B 2 is the plasma to magnetic pressure ratio. "Let us consider the low (but finite) pressure plasma, 0 < β 1. The plasma is assumed to be composed of core cold particles and an admixture of hot protons and electrons. The hot protons make the main contribution to plasma pressure, β p β e . We will use the Maxwellian distribution for the hot particles: where n is the number density, ε = v 2 /2 is the particle energy per unit mass, v is particle velocity, and ε 0 j is the temperature. Hereinafter, the index j indicates the belonging to the protons (p) or the electrons (e), respectively " [12].
The modes with high azimuthal wave numbers (m 1) are considered, which allows description in the transverse WKB approximation [13]. In this case, all perturbed values can be written as the following: where ω is the wave's frequency, k 1 and k 2 are the wave-vector's radial and azimuthal components, respectively. If the coordinate x 2 is represented by the azimuthal angle, then k 2 equals the azimuthal wave number m.
The ULF waves are characterized by the frequencies ω much lower than the gyrofrequency ω c . Such modes are conventionally considered in the gyrokinetics framework, where the wave's electromagnetic field is described by three variables: "φ is the electrostatic potential, b is the parallel magnetic field, and ψ is the potential related to the parallel vector potential, as A = −(ic/ω)∂ψ/∂l " [12]. In this case, the transverse and parallel electric fields of the wave are the following, respectively: The gyrokinetics equations system describing electromagnetic waves in a plasma can be written, using the following three variables [10]: Here, (9) represents the perpendicular Ampere's law: (10)-the parallel Ampere's law or vorticity equation; and (11)-the quasi-neutrality condition. In these equations,L i are integral or integro-differential operators defined in [10] (see their Equations (12)-(14), respectively). Next, we will consider a few simplified cases because this system is too complex to be solved without simplifications.

The Alfvén Mode Equation and Compressional Component of the Alfvén Wave
In the first case, we consider the system (9)-(11) neglecting the parallel electric field of the wave: E = 0. This is a plausible situation since the admixture of cold electrons in the plasma should short-out the parallel electric fields [14]. "In this case, we can consider the perpendicular Ampere's law Equation (9) and vorticity Equation (10) independent of the quasi-neutrality condition (11)" [12]. To simplify, we neglect the contribution of hot electrons since β p β e . In this case, the system (9) and (10) is written in the following simplified form:L whereL A andL M are the Alfvén and compressional mode operators, andL AM andL MA are the coupling operators defined below. "The Alfvén mode operator is as follows: where is the integral over velocity space, is an average for a bounce period τ b , where ±l 0 are reflection points of a magnetically trapped particles with energy ε and magnetic moment µ = v 2 ⊥ /(2B), q is the particle charge, m p is the proton mass, c is the velocity of light, v and v ⊥ is the parallel and transverse components of particles velocity relative to the Earth's magnetic field, respectively: is the drift frequency,Q where ω c is the gyrofrequency, g ⊥ = g 1 g 2 " [12]. The operatorsL T define the toroidal operator as follows:L andL P defines the poloidal operator for cold plasma as follows: where v A = B/ 4πρ is the Alfvén speed, and ρ is the equilibrium mass density. The compressional mode operator is defined by its action on the b function as follows: The operatorsL AM andL MA define the mode coupling as follows: Let us assume that both wave and drift frequencies of protons are much larger than the bounce frequency. Then, Equation (12) is reduced to the following form [15]: This means that in the dipolar geometry, the Alfvén mode (the ψ variable) is always accompanied by the parallel magnetic field (the b variable). Substituting b into Equation (13), we yield the single Alfvén mode equation: Here,L * P is a modified poloidal Alfvénic operator determined in the limit ω ω d p,e , ω b p as the following:L * The solution of Equation (26) gives a profile for ψ. We need the obtained results further.

Drift-Compressional Modes
Some of the observed magnetospheric ULF waves are coined as the compressional waves characterized by the variation of the parallel component of the wave magnetic field [16][17][18]. In the magnetohydrodynamics framework, these waves can be identified with fast or slow MHD modes [19]. To describe such waves in the kinetics, we consider the situation of compression resonance, as in [11]. Consider the second term of Equation (13): "Changing the order of integration, we obtain the following: where ΛK(l , l ) is the function that is obtained as a result of changes in the order of integration, while Λ does not depend on l and l . Equation (13) can formally be considered an inhomogenous Fredholm integral equation of the second kind: where K(l , l ) is the kernel of the integral equation, Λ is a parameter, l I is the distance from the equator to the ionosphere along the field line and b = (ω/c)b " [11]. The solution of Equation (30) is the following: where Λ and b N are the eigenvalues and eigenfunctions of the homogeneous equation as follows:L N is the wavenumber of the harmonic. The eigenfunctions are orthogonal to each other. The normalization condition is the following: where δ NN is the Kronecker delta.
In case Λ → Λ N , b is proportional to the corresponding eigenfunctions of the homogeneous equation b N as follows: Thus, when Λ matches to one of the eigenvalues Λ N , the field-aligned structure of b is defined by the eigenfunction b N and its absolute value goes to infinity, b → ∞. We call this situation the compressional resonance. In this case, to find the parallel structure of the wave magnetic field, we will find the solution of homogeneous Equation (32). To do this, we can change the integration variables µ to λ and reduce Equation (32) to the canonical form, as in [11]: where λ = sin 2 α = µB 0 /ε, α is the pitch angle, B 0 is the magnetic field value in the equatorial plane. Here and below, we use l instead of l . We introduce the following notations: which is the particle path length over bounce period. As shown in [1] for the dipole magnetic field the bounce averaged magnetic drift frequency and L b weakly depend on λ: Therefore we can assume that Λ(λ, ω) is independent of λ. Changing the order of integration in (35), we obtain the following integral equation: with the kernel .

(40)
For convenience, Figure 1 shows the integration area. Since K decreases rapidly away from the equator, we expand the limit of integration to infinity. Then, we reduce the kernel of Equation (39) to a symmetric formK(l, l ) = K(l, l ) B(l )/B(l). As a result, we obtain the second kind Fredholm integral equation with a real symmetric nondegenerate positive-definite kernel:b whereb(l) = b B 0 /B(l) and the following holds: "As can be seen from Equation (41), the scale of solution b along the field line is determined by the scale of kernel K. Since K rapidly decreases away from the equator, the eigenfunctions b N are localized in the equator area. In further calculations, we will use the quadratic approximation for the magnetic field" [11]: where Equation (41) is solved numerically to obtain the first three eigenfunctions b N = b N √ B/B 0 . The corresponding eigenvalues are Λ 1 = 0.5/r 0 , Λ 2 = 1.5/r 0 and Λ 3 = 2.5/r 0 . The eigenfunctions b N are strongly localized on the geomagnetic equator ( Figure 2). It agrees with the observations of the compressional ULF waves in the magnetosphere [16,20]. Other wave parameters (eigenfrequency, increments) can be found by substituting the observed parameters in Λ.

The Parallel Electric Field of the Alfvén Wave
Another problem when solving, by which an integral equation is obtained, is the determination of the Alfvén wave parallel electric field, E . This field is believed to energize the electrons responsible for the aurora [3]. "Since in the standard MHD approach, the Alfvén wave parallel electric field is assumed to be zero, E = 0, the kinetic effects must be taken into account" [12]. Most of the previous works consider hybrid models. In such models, the MHD approximation is used to calculate the frequency and eigenmode structure of the Alfvén waves in a given magnetic field. Then, the kinetic equation for electrons is solved for a given parallel electric current profile. As a result, the parallel electric field of the wave is calculated, which is necessary to maintain a given current [2,21]. The effects of plasma pressure and the wave's parallel magnetic field are usually neglected. In contrast to previous works, we obtain all equations describing the parallel electric field exclusively from the kinetic approach [12].
To obtain an equation governing the parallel electric field of the Alfvén mode, we will make some simplifications: (1) the anisotropy of the particles thermal velocities is ignored; (2) the Larmor radius is small k ⊥ ρ 1; and (3) the waves' eigenfrequency is much higher than the frequency of the diamagnetic drift and the protons' bounce-frequency (ω ω d p,e , ω b p ). In this case, the Equation (11) is the following: where q is the particle charge, m is the particle mass, c is the speed of light and δk s is distribution function perturbation due to the wave-particle interaction. The designation δk s is defined by Equation (3) in [10]. This term is different for electrons and protons. Due to their lightness, electrons are assumed to reflect many times between the mirror points during the wave period. On the contrary, protons, as relatively heavy particles, can be considered sedentary along the magnetic field lines during the wave period [15]. Therefore, we will average the motion by the bounce period only for electrons: (see Equation (16) in [10]). δk s for protons is as follows: Equations (9) and (10) are solved under the assumption of the smallness of the parallel electric potential φ ψ, b , as described above. Taking into account Equation (25), Equations (45)-(47) are reduced after some algebra to the following form: where and are the diamagnetic frequencies due to the radial gradients of density and temperature. As we did above, we change the integration variable µ → λ, where λ = µB 0 /ε. Additionally, we take into account the symmetry of the magnetic line relative to the geomagnetic equator and the ideal conductivity of the ionosphere. Ignoring the terms ω * n j /ω and ω * ε j /ω, we obtain an inhomogeneous second kind Fredholm integral equation (for more details see [12]): with kernels "Here, ς 3 (l) = 2 √ g 1 /R, l I is the reflection point at the ionosphere, B I is the magnetic field value at the ionosphere boundary. For particles with energies of the order of keV, we can consider L b ∼ 4L, where L is the distance from the center of the Earth to the middle of the magnetic field line on the selected L-shell [1]. If we assume that the value β is small (β 1) and take into account the plasma equilibrium condition (4), then we can represent the following" [12]: and After integrating the kernels (54) and (55) and changing the integration variable l → θ, in Equation (53), where θ is the magnetic latitude, we have the following: where f l (θ ) = cos θ √ 1 + 3 sin θ , is the value of electron magnetic drift frequency in the equatorial plane. The new coordinate system is schematically represented in Figure 3A, where θ I = arccos √ R I /L, R I is the ionosphere boundary and R E is the Earth's radius.
To produce a numerical solution of the system (26) and (58) in the poloidal wave case (k 1 k 2 ), we use the parameters: "v A = 1550 km/s is the Alfvén speed in the equator, E a = 10 mV/m is the azimuthal component of the electric field, m = 50 is the azimuthal wave number, n p = n e 3 sm −3 is the particle concentration, β = 0.4, B 0 = 0.34 G, R E = 6371 km, L = 6.6 R E , R I = R E + 1500 km. In this case, the average energy of the hot protons is T p = 5 keV, and the average energy of the hot electrons is T e = 1 keV" [12].
On the first step, we numerically solve Equation (26), obtaining the parallel structure of the Alfvén mode's transverse electric field ψ ( Figure 3B) and the corresponding eigenfrequencies, the principal of which has a frequency ω = 19.02 mHz.
On the next step, we perform a numerical solution of Equation (58) to obtain the parallel potential φ . The result is shown in Figure 3C. The obtained results are in good agreement with the observations [4].

Conclusions and Discussion
In the paper, we considered two cases where the solution of the kinetic equations set (9)-(11) was obtained, using integral equations. Both solutions are special cases that show only a qualitative view of the physical processes in the magnetosphere. Further research will use fewer different approximations to obtain a more accurate physical processes description. One of them is adopting a three-dimensional model of the magnetosphere, which allows one to take into account the day-night asymmetry of the circumterrestrial plasma [22][23][24]. Using the kinetic approach will lead to more complex integral equations.
The kinetic approach is also applicable to the solar atmosphere physics. In addition, a similar theory is used to solve problems of controlled thermonuclear fusion. In this approach, the magnetosphere and outer space act as a natural laboratory for high-energy plasma physics. Therefore, further theoretical research is very important.