Calculations of the integral invariant coordinates I and L ∗ in the magnetosphere and mapping of the regions where I is conserved , using a particle tracer ( ptr 3 D v 2 . 0 ) , LANL * , SPENVIS , and IRBEM

The integral invariant coordinate I and Roederer’s L or L are proxies for the second and third adiabatic invariants, respectively, that characterize charged particle motion in a magnetic field. Their usefulness lies in the fact that they are expressed in more instructive ways than their counterparts: I is equivalent to the path length of the particle motion between two mirror points, whereas L, although dimensionless, is equivalent to the distance from the center of the Earth to the equatorial point of a given field line, in units of Earth radii, in the simplified case of a dipole magnetic field. However, care should be taken when calculating the above invariants, as the assumption of their conservation is not valid everywhere in the Earth’s magnetosphere. This is not clearly stated in state-of-the-art models that are widely used for the calculation of these invariants. The purpose of this work is thus to investigate where in the near-Earth magnetosphere we can safely calculate I and L with tools with widespread use in the field of space physics, for various magnetospheric conditions and particle initial conditions. More particularly, in this paper we compare the values of I and L as calculated using LANL*, an artificial neural network developed at the Los Alamos National Laboratory, SPENVIS, a space environment online tool, IRBEM, a software library dedicated to radiation belt modeling, and ptr3D, a 3-D particle tracing code that was developed for this study. We then attempt to quantify the variations between the calculations of I and L of those models. The deviation between the results given by the models depends on particle initial position, pitch angle and magnetospheric conditions. Using the ptr3D v2.0 particle tracer we map the areas in the Earth’s magnetosphere where I and L can be assumed to be conserved by monitoring the constancy of I for energetic protons propagating forwards and backwards in time. These areas are found to be centered on the noon area, and their size also depends on particle initial position, pitch angle and magnetospheric conditions.


Introduction
The motion of charged particles in the geomagnetic field is complicated, even if one approximates that field with only its dipole component.It is helpful to break down the total motion of the particle into three individual components: gyration around a guiding magnetic field line, bounce along the magnetic field line between magnetic mirror points, and gradient and curvature drift across the magnetic field line in an azimuthal direction around the Earth.Because these components evolve over very different time scales, they are nearly independent of each other and can thus be summed linearly to obtain the total motion (Prölss, 2004).For time variations of the magnetic field that are slow compared to the corresponding timescale of each type of motion, an adiabatic invariant is defined.
The first invariant, µ, is associated with the cyclotron motion of the particle and expresses the constancy of the magnetic flux enclosed by the particle's gyromotion.The second invariant, J , is associated with the bouncing motion along the magnetic field between mirror points and implies that the particle will move in such a way as to preserve the total Published by Copernicus Publications on behalf of the European Geosciences Union.length of the particle trajectory.The third invariant, , is associated with the particle's azimuthal drift around the Earth, and it represents the conservation of magnetic flux encompassed by the guiding drift shell of a particle.
In calculations involving the adiabatic invariants, it is often instructive to use proxy invariant parameters.In the case of calculations concerning the second adiabatic invariant, the integral invariant coordinate I (Roederer, 1970) is defined as where the integral is between two mirror points, B(s) is the magnetic field intensity and ds the line element along the field line, and B m is the magnetic field intensity at the mirror points.I is expressed in distance units (km or R E ) and depends on the length of the particle trajectory along a field line between the two mirror points.In place of the third adiabatic invariant it is convenient to use L * or Roederer's L (Roederer, 1970).L * is defined as where k 0 is the Earth's dipole moment, R E is the radius of the Earth (6370 km) and is the third adiabatic invariant and is defined as integrated along the trajectory of the particle for the entire drift shell.L * physically approximates the distance from the center of the Earth to the equatorial point of a given field line (in R E ) if we assume a dipolar magnetic field for the Earth.L * is also an invariant, since it is inversely proportional to (Roederer, 1970).
A practical way to calculate is to find the intersection C of a series of drift-shell field lines with the Earth's surface and to numerically compute over the cap delineated by C, using the following equation: where λ e (φ) is the dipole latitude of the intersection C at a given longitude φ (Roederer, 1970).

LANL*
In general, calculation of L * is very computationally expensive because it involves an integral that is both twodimensional and global.LANL* aims to address this issue by calculating L * based on a sophisticated dynamic magnetic field model at a fraction of the time required for full drift shell integration using a neural network technique.LANL* V2.0 is an artificial neural network (ANN) for calculating the magnetic drift invariant, L * , based on the Tsyganenko TS05 magnetic field model (Tsyganenko and Sitnov, 2005).The TS05 model is an empirical best-fit representation for the magnetic field, based on data taken from a large number of satellite observations.The Tsyganenko model suite includes subroutines for the current (IGRF) and past (DGRF) internal geomagnetic field models as well as for a dipole internal field.
The artificial neural network consists of two layers.The first layer provides 19 nodes, one for each input parameter for the TS05 model plus additional nodes to help specify the drift shell especially for low Earth orbit.The hidden layer in the neural network contains 20 neurons that are connected to each input node and one output node to produce L * .The ANN was trained using the latest version of the IRBEMlib library in SpacePy, a Python-based tool library for space science, to generate the input-output database, using a constrained truncated Newton algorithm to train an ANN on the input-target data.A second neural network within LANL* V2.0 called LANLmax was created to describe the last closed drift shell (maximum possible value for L * ) under the specified solar wind conditions (Koller and Zaharia, 2011;Yu et al., 2012).

IRBEM-lib
IRBEM-lib (formerly known as ONERA-DESP-LIB) is a freely distributed library of source codes dedicated to radiation belt modeling put together by the Office National d'Etudes Aérospatiales (ONERA-DESP).The library allows the computation of magnetic coordinates and fields for any location in the Earth's environment for various magnetic field models.It is primarily written in Fortran with access to a shared library from IDL or Matlab (Bourdarie and O'Brien, 2009).IRBEM-lib calculates I by tracing the magnetic field line that crosses a given point, calculating the integrand of Eq. (1) along all line elements of the field line1 .The third invariant is evaluated in IRBEM-lib using the numerical method described by Roederer as given in the introduction, where the magnetic drift shell is defined as a set of magnetic field line segments where all the segments are characterized by the same shell parameter and the same mirror-point magnetic field intensity.Roederer's shell parameter L * is then deduced directly from the value of the third invariant using Eq.(2)2 .The latest version of IRBEM-lib can be found at http:// sourceforge.net/projects/irbem/.

SPENVIS
The European Space Agency (ESA) Space Environment Information System (SPENVIS) provides standardized access to models of the hazardous space environment through a World Wide Web interface (Heynderickx et al., 2004).SPENVIS includes magnetic field models implemented by means of the UNILIB library for magnetic coordinate evaluation, magnetic field line tracing and drift shell tracing.Among these models are TS05 external and IGRF internal models.In UNILIB, the integral invariant I is evaluated using a Runge-Kutta integration technique to evaluate Eq. (1) for a temporary magnetic field line, also traced in UNILIB (Schmitz et al., 2000).The third invariant is evaluated in UNILIB using Roederer's numerical method as for the case of IRBEM-lib3 .

ptr3D v2.0
The calculations of I and L * were also performed using a 3-D particle tracing code that was developed for this study (ptr3D v2.0).This code traces the full 3-D Lorentz motion of single charged particles by integrating the relativistic Lorentz equation in the same geomagnetic field model that was used in the above simulations, the TS05, for direct comparison between all the models considered here.In the particle tracing model the integration is performed by means of Hamming's modified predictor-corrector method in conjunction with a fourth-order Runge-Kutta method for initialization (Ralston and Wilf, 1977;Ralston, 1962).For the calculations of the magnetic field, the GEOPACK-2008 implementation of the TS05 magnetic field model was used (Kuznetsova, 2006).ptr3D v2.0 calculates I by directly evaluating Eq. ( 1) for each step of the simulation.The third invariant is calculated following the method described by Roederer (1970).Roederer's shell parameter L * is then deduced directly from the value of the third invariant using Eq. ( 2).

Calculations of I
The integral invariant I was calculated for various geocentric distances (in R E , GSM) using IRBEM and SPENVIS, for particles starting at magnetic local noon and magnetic local midnight, during quiet and disturbed magnetospheric conditions.Four initial pitch angles (15, 30, 45 and 60 • ) and five initial distances (4-8 R E in steps of 1 R E ) were used.Using the ptr3D particle tracer, I was calculated for 3 particle energies (500 keV, 1 MeV and 4 MeV), 5 initial distances (4-8 R E in steps of 1 R E ) and 12 initial particle gyrophases, also during quiet and disturbed magnetospheric conditions (in a static magnetic field), and the final I was estimated as the median of the results for all gyrophases.
In Figs. 1 to 3 the integral invariant I is shown as a function of the distance of the starting point on the GSE x axis, in R E .Four families of curves are plotted, one for each initial pitch angle.The calculations in Figs. 1 and 2 were performed for quiet magnetospheric conditions (23 February 2008, 17:55 UT was selected), for initial starting points at noon MLT and midnight MLT respectively, whereas the calculations in Fig. 3 were performed for disturbed magnetospheric conditions (8 September 2002, 01:00 UT), for initial starting points at noon MLT.
In Fig. 1, it can be seen that overall there is good agreement in the calculations of I between all three models for all four initial pitch angles.It can also be seen that I gradually increases as one moves towards greater geocentric distances and also that I is larger for greater pitch angles.Both of these cases can be explained if the connection between I and the bounce path length is considered: in the first case, the geocentric distance of the particle increases, and the magnetic field lines become longer between two given mirror points.Therefore the particle's path length increases accordingly.For the second case, particles with smaller equatorial pitch angles are mirrored further along a magnetic field line and therefore traverse greater distances along said field line than particles with greater equatorial pitch angles.A small deviation is observed in the results from ptr3D, which increased for increasing particle energy.This deviation is more pronounced for smaller pitch angles and for the 4 MeV case, whereas it is very small for the other two energy cases, and it becomes negligible for larger pitch angles.
In Fig. 2, where the calculations at midnight MLT are shown, the results from IRBEM and SPENVIS agree quite well for all cases.The results from ptr3D deviate from the results from IRBEM and SPENVIS in the following ways: there is a spread in the results that increases for increasing distances.The results for I deviate more from those of SPENVIS and IRBEM the greater the energy of the particle.This spread also becomes wider as the pitch angle decreases.
In Fig. 3 a similar trend is observed, where for small pitch angles there is good agreement between SPENVIS and IRBEM and a deviation of the results from the 3-D tracer that is proportional to the particle energy.In this case though only the results for the 4 MeV particle deviate significantly from the rest.Furthermore, I appears to be larger as calculated by ptr3D for higher particle energies, contrary to the figures above.
There is generally good agreement between the results from IRBEM and SPENVIS and those from ptr3D for 500 keV and 1 MeV particles, except for the case of 15 • pitch angle, where there is a small deviation at 5 and 6 R E .Again  it can be seen that the results for the 4 MeV particles deviate significantly and that this deviation is a function of distance and pitch angle, even though in this case we get a larger deviation for larger pitch angles and the values for I are larger than those calculated with IRBEM and SPENVIS.
I was also calculated during disturbed conditions, for particles initiating from midnight MLT.For these conditions the results from IRBEM and SPENVIS appear to agree fairly well, whereas the calculations from ptr3D deviate significantly for the various particle energies.Furthermore, a large number of particles precipitated into the atmosphere in the simulation, or otherwise failed to complete the necessary trajectory between mirror points in order for Eq. ( 1) to be calculated.The results for this case are therefore not presented here.

Calculations of L *
In the following, L * was calculated using LANL*, IRBEM, SPENVIS and ptr3D for particles initiating their trajectory at 12:00 MLT on the XGSE axis, for five initial distances from 4 to 8 R E and for three initial pitch angles of 30, 60 and 90 • .The value of L * was also calculated for the last closed drift shell, called L * max , using LANL* for the three pitch angles listed above.The results for the calculations during quiet and disturbed magnetospheric conditions are shown respectively in Figs. 4 to 9.These figures give the calculated L * as a function of distance in R E (in GSM) of the particle starting point on the magnetic equator and the point of calculation for all the other models used.L * max calculated through LANL* is also shown as a horizontal line.
Generally, for the quiet conditions case, the results from all the models tend to agree more at smaller distances (4-  6 R E ) and less further away (7-8 R E ).Also, the larger the initial pitch angle the greater the spread of the calculated L * .For example, the standard deviation becomes close to 2 for a distance of 8 R E and a pitch angle of 30 • .L * max is calculated to be around 9 R E for all initial pitch angles.
For the disturbed conditions case, similar trends are observed, albeit more accentuated.The results from LANL* and IRBEM agree relatively well, as do those from SPEN-VIS, where available.The results from ptr3D deviate significantly from those of the other models for distances greater than 4 R E .The standard deviation was found to be as high as 3 for the particles initiated at a distance of 8 R E with a 60 • pitch angle, and is significant for distances greater than 4 R E .
Similarly to the simulations above, some of the 12 particles of various initial gyrophases precipitated or otherwise failed to complete a full revolution around the Earth; these particles were not taken into consideration when averaging  the results for each initial gyrophase.If more than half of the particles failed to complete a rotation around the Earth, no L * was calculated.

Mapping regions of constant I
Next, we demonstrate at which magnetic longitude the conservation of I is broken, for different particle starting conditions.We thus map the areas where I and therefore L * cannot be safely calculated.
Using ptr3D, I was calculated for particles propagating both forwards and backwards in time, during the same two periods of quiet and disturbed solar wind conditions as above (23 February 2008, 17:55 UT and 8 September 2002, 01:00 UT respectively) starting at local noon, for 2 initial pitch angles (30 and 60 • ), 5 initial distances (4-8 R E ) and for 12 initial particle gyrophases.For each pitch angle, the values of I were plotted for each initial distance and initial gyrophase, both for forwards-and backwards-traced particles, starting at local noon, as a function of the particle's azimuth angle.In the resulting plot, for both directions of propagation, a dashed vertical line marks the approximate point where I stops being constant (see Figs. 10 and 11).Subsequently, for each case of solar wind conditions and initial pitch angles, a map was created, depicting the areas where I remained relatively constant for each case of initial particle position (see Fig. 12).The purpose of this section is therefore to demonstrate at which magnetic longitude the adiabaticity of I is broken, for different particle starting conditions as well as for different geomagnetic conditions.To this extent, in Fig. 11 we map the areas where I and therefore L * cannot be assumed to remain constant throughout a particle drift shell.With this map we demonstrate in a graphic representation the magnetic longitudes and distances from the Earth where I ceases to be adiabatic.In these areas the generalpurpose models and tools described above, such as IRBEM, LANL* and SPENVIS, cannot be safely used to calculate the values of the adiabatic invariant I and therefore L * .

Quiet conditions
In the case of the 30 • initial pitch angle, I remains constant throughout the path of the particle around the Earth for an initial particle distance of 4 and 5 R E .For other initial distances there appears to be a region in the nightside where I is no longer conserved.This region becomes larger with increasing distance.In the case of the 60 • initial pitch angle, I remains constant throughout the path of the particles around the Earth for initial particle distances of 4-6 R E .Similar to the case of particles with a 30 • initial pitch angle, there are regions where I is not constant and these regions are larger the longer the initial distance.Generally, the extent of these regions is smaller in the case of the 60 • initial pitch angle particle.

Disturbed conditions
In the case of disturbed solar wind conditions, the regions of constant I are generally smaller than in the case of quiet conditions.The symmetry of the regions around the XGSE axis is no longer there as for the quiet conditions case, since the magnetic field is also no longer symmetric along the same axis.For initial distances of 4 and 5 R E the regions are larger in the case of 60 • initial equatorial pitch angle, following the trend exhibited in the quiet conditions case, where the regions of constant I are larger for larger pitch angles.For distances of 6-8 R E this trend is reversed and the regions of constant I of the 60 • pitch angle case become smaller than those for the 30 • pitch angle case.The breaking of the adiabaticity of I observed here is due to the high curvature of the magnetic field lines in the respective areas.Therefore since the particle does not follow constant flux tubes, it cannot be assumed that it conserves I .

Conclusions
Using the ptr3D v2.0 particle tracer, LANL*, IRBEM-lib and SPENVIS, we quantified the variations in the calculations of I and L * between these models, for various particle initial starting positions in geocentric distances, for various initial pitch angles, for both quiet and disturbed magnetospheric conditions and for particles initiating their motion both in the dayside and nightside.
The results for the calculations of I in the dayside show that the models used are in good agreement for all geocentric distances of the particle starting positions, all pitch angles, and both for quiet and disturbed magnetospheric conditions.In the nightside and for quiet magnetospheric conditions, there is good agreement between models only for small geocentric distances of the particle starting positions, for all initial pitch angles.For larger distances, there is an increasing disagreement between these results, and differences are more accentuated for smaller pitch angles.
Generally, the same trends are observable for the calculations of L * between the various models.For quiet magnetospheric conditions the results from the models are in relative agreement for smaller geocentric distances of the particle starting positions and start to deviate with increasing distances and initial pitch angles.For disturbed magnetospheric conditions this deviation is more accentuated.
Using ptr3D we mapped the areas in the Earth's magnetosphere where I , and consequently also L * , can be assumed to be conserved, for two initial pitch angles, and for both quiet and disturbed magnetospheric conditions.This too was performed by monitoring the constancy of I for energetic protons propagating forwards and backwards in time.Results for quiet magnetospheric conditions show that the regions where I cannot be assumed to be conserved appear between a GSM distance of 5-7 R E in the nightside, centered at the midnight point depending on the pitch angle, and those areas expand on the nightside for larger distances.These areas are more extensive for larger particle pitch angles and appear to be symmetrical around the plane defined by the midnight-noon line and the Earth magnetic dipole axis.For disturbed magnetospheric conditions, the areas where I cannot be assumed to be conserved start to appear between a GSM distance of 3-5 R E on the nightside.
In the discussions of particle transport, energization and loss in the Earth's radiation belts, a major question concerns the relative contribution between wave-particle interactions vs. radial diffusion, which is generally best discussed in terms of phase-space density, calculated at constant adi-abatic invariants.From the discussions herein, it is evident that caution should be exercised when considering the second and third adiabatic invariants to remain constant across all Lshells and local times within the radiation belts as well as for all particle energies and all geomagnetic conditions.In particular, in regions where the results from the various models diverge from the results from the particle tracer, which most closely follows the calculations of the invariants, we can conclude that the models should be used with caution, the lack of confidence in them being analogous to the magnitude of this divergence.It has been demonstrated that under extreme curvature of the magnetospheric magnetic field, particles of high energy and low pitch angles cannot be considered to remain adiabatic in terms of their second and third invariants.
The physical mechanism that leads to breaking of the invariants in the regions illustrated does not involve temporal variations in the magnetic field of timescales shorter than the associated timescales of the second and third invariants, i.e., the bounce period and drift period, as the fields used in the simulations above are all static.Instead, the breaking of the invariants in the above is associated with deviations of the magnetic field from a dipole configuration: in the definition of the invariants, in order for the second adiabatic invariant to remain constant it is required that the magnetic field between two mirroring points does not change much in one bounce period as the particle's guiding center drifts across field lines.Similarly, in order for the third adiabatic invariant to remain constant, it is required that the magnetic flux through the guiding center orbit of a particle around the Earth should remain constant.However during active geomagnetic conditions the curvature of the field lines in the nightside of the Earth in combination with the large gyroradii of large-energy particles leads to deviations from these conditions that need to be taken into account.
The present paper by no means aims to serve as a guideline of the adiabaticity of particles at all energies, pitch angles and geomagnetic conditions; instead, it aims to raise awareness and caution in using general-purpose models and tools, such as IRBEM, LANL* and SPENVIS, to calculate the values of the adiabatic invariants in regions and cases where they are not well defined.its next version (ptr3D V3.0) will be released as a general particle tracer code that can be used for any range of particle energies, times or regions.

Figure 1 .
Figure 1.Calculations of I as a function of initial distance (in R E ) at 12:00 MLT, for pitch angles of 15, 30, 45 and 60 • and quiet solar wind conditions.

Figure 2 .
Figure 2. Calculations of I as a function of initial distance (in R E ) at 00:00 MLT, for pitch angles of 15, 30, 45 and 60 • and quiet solar wind conditions.

Figure 3 .
Figure 3. Calculations of I as a function of initial distance (in R E ) at 12:00 MLT, for pitch angles of 15, 30, 45 and 60 • and disturbed solar wind conditions.

Figure 4 .
Figure 4. Calculations of L * as a function of initial distance (in R E ) at 00:00 MLT, for an initial pitch angle of 30 • during quiet solar wind conditions.Error bars represent the standard deviation between results from the various models.

Figure 5 .
Figure 5. Calculations of L * as a function of initial distance (in R E ) at 00:00 MLT, for an initial pitch angle of 60 • during quiet solar wind conditions.Error bars represent the standard deviation between results from the various models.

Figure 6 .
Figure 6.Calculations of L * as a function of initial distance (in R E ) at 00:00 MLT, for an initial pitch angle of 90 • during quiet solar wind conditions.Error bars represent the standard deviation between results from the various models.

Figure 7 .
Figure 7. Calculations of L * as a function of initial distance (in R E ) at 00:00 MLT, for an initial pitch angle of 30 • during disturbed solar wind conditions.Error bars represent the standard deviation between results from the various models.

Figure 8 .
Figure 8. Calculations of L * as a function of initial distance (in R E ) at 00:00 MLT, for an initial pitch angle of 60 • during disturbed solar wind conditions.Error bars represent the standard deviation between results from the various models.

Figure 9 .
Figure 9. Calculations of L * as a function of initial distance (in R E ) at 00:00 MLT, for an initial pitch angle of 90 • during disturbed solar wind conditions.Error bars represent the standard deviation between results from the various models.

Figure 10 .
Figure 10.I as a function of the particle's azimuth angle for the case 1 MeV protons starting at GSM coordinates [8, 0, 0] (in R E ) with initial pitch angles of 30 • , and initial gyrophases of 0-330 • (30 • step).Particles propagating forwards in time are shown in blue, while those propagating backwards are shown in red.

Figure 11 .
Figure 11.The Lorenz trace of the forwards (blue) and backwards (red) propagating particle is plotted.The region where I is constant according to Fig. 10 is shown in magenta.The contours of constant magnetic field strength are also plotted.

Figure 12 .
Figure 12.Plots of the regions of constant I , for quiet and disturbed solar wind conditions, and 30 and 60 • initial equatorial pitch angles, for particle starting distances 4-8 R E .