Eccentric Dipole Evolution during the Last Reversal, Last Excursions, and Holocene Anomalies. Interpretation Using a 360-Dipole Ring Model

: The eccentric dipole (ED) is the next approach of the geomagnetic ﬁeld after the generally used geocentric dipole. Here, we analyzed the evolution of the ED during extreme events, such as the Matuyama-Brunhes polarity transition (~780 ka), the Laschamp (~41 ka) and Mono Lake (~34 ka) excursions, and during the time of two anomalous features of the geomagnetic ﬁeld observed during the Holocene: the Levantine Iron Age Anomaly (LIAA, ~1000 BC) and the South Atlantic Anomaly (SAA, analyzed from ~700 AD to present day). The analysis was carried out using the paleoreconstructions that cover the time of the mentioned events (IMMAB4, IMOLEe, LSMOD.2, SHAWQ-Iron Age, and SHAWQ2k). We found that the ED moves around the meridian plane of 0–180 ◦ during the reversal and the excursions; it moves towards the region of the LIAA; and it moves away from the SAA. To investigate what information can be extracted from its evolution, we designed a simple model based on 360-point dipoles evenly distributed in a ring close to the inner core boundary that can be reversed and their magnitude changed. We tried to reproduce with our simple model the observed evolution of the ED, and the total ﬁeld energy at the Earth’s surface. We observed that the modeled ED moves away from the region where we set the dipoles to reverse. If we consider that the ring dipoles could be related to convective columns in the outer core of the Earth, our simple model would indicate the potential of the displacement of the ED to give information about the regions in the outer core where changes start for polarity transitions and for the generation of important anomalies of the geomagnetic ﬁeld. According to our simple model, the regions in which the most important events of the Holocene occur, or in which the last polarity reversal or excursion begin, are related to the regions of the Core Mantle Boundary (CMB), where the heat ﬂux is low.


Introduction
The Earth's magnetic field is generated and sustained in the outer core by the convection that drives the geodynamo. Probably the most remarkable feature of the geomagnetic field is the occurrence of polarity reversals, which have taken place many times in the past at irregular intervals (e.g., [1]). Reversals have also been observed in geodynamo simulations (e.g., [2]), but their cause has yet eluded a convincing explanation (e.g., [3]). In addition to reversals, there are anomalous events, named excursions, where the geomagnetic direction departs greatly from the Geocentric Axial Dipole (GAD), sometimes achieving the reversed direction for a short time. Excursions have been interpreted as large secular variation events, a pair of full reversals, or as aborted reversals (e.g., [4][5][6][7][8]). 2 of 16 During reversals and excursions, the dipolar field energy decreases and becomes comparable to the non-dipole field energy, but when analyzing the behavior of the geomagnetic field during these events it is common to represent the Virtual Geomagnetic Pole (VGP) paths, that is, the apparent trajectories of the position of the poles of the dipole field as seen from one point of the Earth's surface. When using that representation, a remarkable effect is observed: the VGP paths during the last reversals seem to be constrained to the American and opposite longitudes [9] that correspond to "cold" zones of the Core-Mantle Boundary (CMB) where the velocity of the seismic waves is higher than the global average [10][11][12]. It is also the case of the VGPs observed during the Blake excursion [13], but not during the Laschamp excursion [14]. However, these results have been questioned referring to the confidence of the sedimentary data used to obtain the VGP [15,16]. Other studies relate the VGP with the GAD, evaluating the tilted dipole during the Holocene (e.g., [17]). In this context, we think that the use of VGP paths does not help to visualize what is happening within the core during these extreme events. It would seem more appropriate to represent the location and evolution of the eccentric dipole (ED) rather than the VGP. In addition, the ED provides a useful tool for detecting the hemisphere within which the Earth's magnetic field displays the strongest magnitude, which is the hemisphere towards which the ED is offset [18].
On the other hand, simple models, like the "domino" model [19]-an Ising-Heisenberg model of interacting magnetic spins-, seem to reproduce some of the statistical properties of the occurrence of reversals and excursions. In detail, the "domino" model consists of a system of interacting macro-spins with freedom to tilt aligned along a ring. To some extent, the spins of the "domino" model might be associated with the outer core convection columns. In this context, the upwellings rising from the Inner-Core Boundary (ICB) could be considered as tilted or reversed spins that progress interacting with neighboring spins.
In this paper, we analyzed the evolution of the ED path during the Matuyama-Brunhes transition (~780 ka) and the Laschamp (~41 ka) and Mono Lake (~34 ka) excursions by using current geomagnetic field paleoreconstructions: the IMMAB4 [20], the IMOLEe [14], and the LSMOD.2 [21]. In addition, we also analyzed the evolution of the ED during two major geomagnetic anomalies that occurred during the last 3 kyr: the Levantine Iron Age Anomaly (~1000 BC), covered by the SHAWQ-Iron Age [22], and the South Atlantic Anomaly (analyzed from~700 AD to present), included in the SHAWQ2k [23]. Inspired by the "domino" model, we used a simple model based on 360-dipole ring (360-DRM) to compare the previous results. Although the 360-DRM is very simple, it helps to visualize the hemisphere in which these extreme events originate and evolve.

Paleoreconstructions
For analyzing the last reversal, we used the IMMAB4 geomagnetic field paleoreconstruction [20] and for the Laschamp excursion, the IMOLEe paleoreconstruction proposed by Leonhardt et al. [14] and the LSMOD.2 paleoreconstruction published by Korte et al. [21]. This last model is also able to resolve the Mono Lake excursion.
The IMMAB4 and IMOLEe use a Bayesian iterative approach. This method tries to find iteratively the most probable evolution of the field without imposing previous physical or geometric assumptions on the data [14,20]. In each iteration, the output results are compared with independent data not used in the model. For the reversal, the IMMAB4 covers a time window between 795 and 765 ka. The input data are mainly sedimentary records. The limited volcanic data available were used to calibrate the relative paleointensity registered by the sediments. For the development of the IMOLEe to reconstruct the Laschamp excursion, one sequence of volcanic data and five high-quality sediment cores that registered the excursion were used in the iterative procedure, spanning from 45 to 35 ka. On the other hand, the LSMOD.2 spans from 50 ka to 30 ka and recovers the geomagnetic field morphology at the times of Laschamp and Mono Lake excursions. This model was constructed by using the whole sedimentary and volcanic database and For analyzing the most recent Holocene geomagnetic anomalies, we selected the SHAWQ paleoreconstructions family: the SHAWQ-Iron Age [22] which covers between 1300 BC and 0 for the Levantine Iron Age Anomaly and the SHAWQ2k [23], from 100 BC to 1900 AD, for the South Atlantic Anomaly. These paleoreconstructions are based exclusively on thermoremanent data (that is, archaeomagnetic and volcanic data) to avoid the smoothing effect of sedimentary data. In addition to statistical criteria, both paleoreconstructions include quality criteria based on the laboratory protocol (for paleointensities) and the number of specimens used [24]. A weighting scheme, which gives 10 times more weight to quality data than to the rest of the dataset, was applied. They also follow the classical SHA technique in space and cubic B-splines in time.

Eccentric Dipole (ED)
To study the evolution of the ED, we calculated its location following the expressions derived by Schmidt [25] and Bartels [26]. From these equations, the resulting dipole is parallel to the geocentric dipole with its center shifted from the Earth's center. Some authors (e.g., [27,28]) pointed out that this approach might not represent the best fitting ED for the geomagnetic field. However, the methodology they proposed to obtain the ED was more complicated and their properties only differ by 10% compared with the Schmidt and Bartels method [28]. For that reason, most of the ED studies use the definition of Schmidt and Bartels (e.g., [29][30][31][32]), and it is the one that we used in this work.
The ED location relative to the Earth's center is specified completely by the first eight Gauss coefficients (i.e., dipolar and quadrupolar harmonic degrees) given by a paleoreconstruction of the geomagnetic field: g 0 The equations used to compute the deviation from the center of the Earth in the three dimensions were [31]: where R E is the Earth's mean radius equal to 6371.2 km. The Gauss coefficients provided by the paleoreconstructions described in Section 2.1 were substituted in Equations (1) and (2) to study the position of the ED and its time evolution.

The 360-Dipole Ring Model (360-DRM)
Inspired by the "domino" model of Mori et al. [19], we designed a simplified model which consists of 360 identical and static (non-interacting) dipoles parallel to the Earth axis and aligned along an equatorial ring ( Figure 1). Dipoles (one dipole per degree in longitude) are located at a distance of 1300 km from the Earth's center (that is, close to the tangent cylinder and the ICB). The magnetic field produced by the distribution of dipoles was calculated on a 1500points mesh evenly distributed over the Earth's surface. From this synthetic dataset, a Spherical Harmonic expansion was obtained (up to the maximum harmonic degree n = 10). Thus, the geomagnetic field accepts matrix representation as where D is the geomagnetic components data vector in the 1500-points regular mesh, M the equations of condition matrix (the matrix that relates the Gauss coefficients to geomagnetic components and depends on the geographical coordinates), and G the Gauss coefficients vector. A least-squares inversion was implemented for each temporal step, whose solution is the one that minimizes the sum of the squares of the residuals [33]. The least-squares solution for the Gauss coefficients is given by Then, the location of the ED of the 360-DRM can be calculated from the Gauss coefficients using Equations (1) and (2).
Each temporal step in this model was defined by a change in the polarity of one or more dipoles (or a change in their magnetic moment). The evolution of the model was set trying to adjust simultaneously: the ED path and the evolution of the magnetic energy at the Earth's surface. This adjustment also allowed us to assign real-time to each temporal step. For simplicity, once a dipole was inverted or modified, in the next step the two neighboring pairs of dipoles on either side were modified.
Please note that it is a simple and conceptual model that aims to help visualize where the changes during the extreme geomagnetic events occur in the outer core. The 360-DRM also helps to explore the potential of the ED to address them, but not to explain in detail the dynamics phenomena that happen during these events.

The Matuyama-Brunhes Transition
The evolution of the Earth's magnetic field during the last reversal, the Matuyama-Brunhes 780 kyr ago, was analyzed from the IMMAB4 (up to its maximum degree 4). The location of the ED during the transition is shown in Figure 2a along with the normalized The magnetic field produced by the distribution of dipoles was calculated on a 1500-points mesh evenly distributed over the Earth's surface. From this synthetic dataset, a Spherical Harmonic expansion was obtained (up to the maximum harmonic degree n = 10). Thus, the geomagnetic field accepts matrix representation as where D is the geomagnetic components data vector in the 1500-points regular mesh, M the equations of condition matrix (the matrix that relates the Gauss coefficients to geomagnetic components and depends on the geographical coordinates), and G the Gauss coefficients vector. A least-squares inversion was implemented for each temporal step, whose solution is the one that minimizes the sum of the squares of the residuals [33]. The least-squares solution for the Gauss coefficients is given by Then, the location of the ED of the 360-DRM can be calculated from the Gauss coefficients using Equations (1) and (2).
Each temporal step in this model was defined by a change in the polarity of one or more dipoles (or a change in their magnetic moment). The evolution of the model was set trying to adjust simultaneously: the ED path and the evolution of the magnetic energy at the Earth's surface. This adjustment also allowed us to assign real-time to each temporal step. For simplicity, once a dipole was inverted or modified, in the next step the two neighboring pairs of dipoles on either side were modified.
Please note that it is a simple and conceptual model that aims to help visualize where the changes during the extreme geomagnetic events occur in the outer core. The 360-DRM also helps to explore the potential of the ED to address them, but not to explain in detail the dynamics phenomena that happen during these events.

The Matuyama-Brunhes Transition
The evolution of the Earth's magnetic field during the last reversal, the Matuyama-Brunhes 780 kyr ago, was analyzed from the IMMAB4 (up to its maximum degree 4). The location of the ED during the transition is shown in Figure 2a along with the normalized (dividing by the maximum value) total energy at the Earth's surface calculated summing up the decomposition of the spatial power spectra at each epoch for each degree (following [34]): The departure of the ED from the Earth's center starts at around 780 ka. The location of its center seems to follow a meridian plane (around 20 • E) towards the Atlantic sector. The minimum of energy is reached at around 775 ka. Instabilities produced in the calculation process by higher harmonic degrees make the use of the ED approximation useless between 775.5-773.5 ka. The recovery of the field starts around 772 ka and it is displayed in Figure 2a as a return of the ED to the Earth's center following the meridian plane of approximately 200 • E. Regarding the Z-axis, when projecting the ED in the 20 • E-200 • E meridian, it seems to follow a linear path inclined around 23 • , moving southwards at the beginning of the transition and coming back from the North Hemisphere in the recovering phase (side view, Figure 2a). (dividing by the maximum value) total energy at the Earth's surface calculated summing up the decomposition of the spatial power spectra at each epoch for each degree (following [34]): The departure of the ED from the Earth's center starts at around 780 ka. The location of its center seems to follow a meridian plane (around 20° E) towards the Atlantic sector. The minimum of energy is reached at around 775 ka. Instabilities produced in the calculation process by higher harmonic degrees make the use of the ED approximation useless between 775.5-773.5 ka. The recovery of the field starts around 772 ka and it is displayed in Figure 2a as a return of the ED to the Earth's center following the meridian plane of approximately 200° E. Regarding the Z-axis, when projecting the ED in the 20° E-200° E meridian, it seems to follow a linear path inclined around 23°, moving southwards at the beginning of the transition and coming back from the North Hemisphere in the recovering phase (side view, Figure 2a). . Blue/red circumference indicates the first dipoles to change their polarity (black triangle) and the propagation directions (black arrows). Note that in (b) the evolution of the ED position is dimensionless and it is presented for a shorter time interval (the color scale is different) to remark the direction of the displacement. Also note that the ED is virtually coincident with the origin before and after the selected shorter period (see Figure S1). Coastlines are plotted only for longitude reference (note they are projected on a sphere of arbitrary radius).
The apparent path followed by the ED is surprisingly very consistent: it is confined between 10° E-20° E meridians at the beginning of the transition and by the 190° E-200° E meridians during the recovering phase, which are sectors nearly symmetric with respect to the rotation axis. This simple behavior suggests us to model this transition with our 360-dipole ring model (360-DRM). Please note that, since we wanted to keep our model simple, we located the 360 dipoles at the equatorial plane, and thus we could not reproduce the ED movement in the Z-axis.
Displacement of the ED is in the opposite direction from the sector where the dipoles start to reverse (Figure 2b). That is, the displacement of the ED towards the Atlantic sector suggests that the reversal started in the Pacific Hemisphere. In this geometric context, we modeled the reversal by a simple ring model in which the dipole that reverses, at the first step, is the one located at 190° E. In the second step, the two neighboring dipoles also . Blue/red circumference indicates the first dipoles to change their polarity (black triangle) and the propagation directions (black arrows). Note that in (b) the evolution of the ED position is dimensionless and it is presented for a shorter time interval (the color scale is different) to remark the direction of the displacement. Also note that the ED is virtually coincident with the origin before and after the selected shorter period (see Figure S1). Coastlines are plotted only for longitude reference (note they are projected on a sphere of arbitrary radius).
The apparent path followed by the ED is surprisingly very consistent: it is confined between 10 • E-20 • E meridians at the beginning of the transition and by the 190 • E-200 • E meridians during the recovering phase, which are sectors nearly symmetric with respect to the rotation axis. This simple behavior suggests us to model this transition with our 360-dipole ring model (360-DRM). Please note that, since we wanted to keep our model simple, we located the 360 dipoles at the equatorial plane, and thus we could not reproduce the ED movement in the Z-axis.
Displacement of the ED is in the opposite direction from the sector where the dipoles start to reverse (Figure 2b). That is, the displacement of the ED towards the Atlantic sector suggests that the reversal started in the Pacific Hemisphere. In this geometric context, we modeled the reversal by a simple ring model in which the dipole that reverses, at the first step, is the one located at 190 • E. In the second step, the two neighboring dipoles also reverse, and so on. The procedure continues until the last dipole located at 10 • E is reversed. To estimate the time duration of each step, i.e., the time between two consecutive reversed dipoles, we compared the normalized energy of the geomagnetic field at the Earth's surface of both IMMAB4 paleoreconstruction and the 360-DRM. This provided a mean rate of 9 reversed dipoles per kyr for the beginning of the transition and a faster recovery of 18 reversed dipoles per kyr. The resulting evolution of the ED from the 360-DRM and the normalized energy at Earth's surface are also shown in Figure 2b. Please note that the normalized total energy was calculated up to maximum degree 4 for both the IMMAB4 paleoreconstruction and the 360-DRM. It is worth noting that although the 360-DRM was generated up to degree 10, just the first 24 Gauss coefficients were significant, i.e., up to degree 4 ( Figure S2).

The Laschamp Excursion
Two paleoreconstructions were used to study the Laschamp excursion: the IMOLEe and the LSMOD.2. The LSMOD.2 originally spans from 50 to 30 ka, but to homogenize we considered the same time interval that the IMOLEe, between 45 and 35 ka. Figure 3a,b show the evolution of the ED for both paleoreconstructions (following the same structure of Figure 2). A displacement of the ED towards 180 • E near to the energy minimum (~41 ka) in both paleoreconstructions can be observed, although the departure from the Earth's center is more limited in the case of the IMOLEe paleoreconstruction. The displacement of the ED is confined in the equatorial plane except for the period related to the energy minimum in the LSMOD.2 (between 41.5 ka and 40 ka).
Geosciences 2021, 11, x FOR PEER REVIEW 6 of 16 reverse, and so on. The procedure continues until the last dipole located at 10° E is reversed. To estimate the time duration of each step, i.e., the time between two consecutive reversed dipoles, we compared the normalized energy of the geomagnetic field at the Earth's surface of both IMMAB4 paleoreconstruction and the 360-DRM. This provided a mean rate of 9 reversed dipoles per kyr for the beginning of the transition and a faster recovery of 18 reversed dipoles per kyr. The resulting evolution of the ED from the 360-DRM and the normalized energy at Earth's surface are also shown in Figure 2b. Please note that the normalized total energy was calculated up to maximum degree 4 for both the IMMAB4 paleoreconstruction and the 360-DRM. It is worth noting that although the 360-DRM was generated up to degree 10, just the first 24 Gauss coefficients were significant, i.e., up to degree 4 ( Figure S2).

The Laschamp Excursion
Two paleoreconstructions were used to study the Laschamp excursion: the IMOLEe and the LSMOD.2. The LSMOD.2 originally spans from 50 to 30 ka, but to homogenize we considered the same time interval that the IMOLEe, between 45 and 35 ka. Figure 3a,b show the evolution of the ED for both paleoreconstructions (following the same structure of Figure 2). A displacement of the ED towards 180° E near to the energy minimum (~41 ka) in both paleoreconstructions can be observed, although the departure from the Earth's center is more limited in the case of the IMOLEe paleoreconstruction. The displacement of the ED is confined in the equatorial plane except for the period related to the energy minimum in the LSMOD.2 (between 41.5 ka and 40 ka).  Following the same approach as for the Matuyama-Brunhes transition, we adjusted the 360-DRM to recreate the evolution of the ED and the energy at the Earth's surface (Figure 3c). The normalized total energy was calculated up to the maximum degree provided by the IMOLEe paleoreconstruction, which is n = 5, in the three cases (IMOLEe, LSMOD.2, and the 360-DRM).
Since excursions could be considered as aborted reversals (e.g., [4,8]), we started with the simplest model in which the reversed region is located at 0 • E and spans up to 270 • E and 90 • E. That is, 180 dipoles are reversed (with a rate of 22 reversed dipoles per kyr) and then the reversal is aborted and the dipoles return to the original configuration (with 15 reversed dipoles per kyr). The 360-DRM reproduces the path of the ED along the 180 • E meridian plane observed in the two paleoreconstructions and the general evolution of the energy at the Earth's surface.
In terms of the total field energy, our model ( Figure 3c) returned to the initial values at 35 ka, while the field energy given by the paleoreconstructions was not completely recovered at that time. This is probably related to the existence of the Mono Lake excursion (~34 ka). Since the LSMOD.2 spans from 50 to 30 ka, we also evaluated the path of the ED during this excursion (see Supplementary Material), identifying also a predominance of the ED location towards 180 • E ( Figure S3).

The Levantine Iron Age Anomaly (LIAA)
During the Iron Age, two intensity peaks of the geomagnetic field were recorded in the Levant region (between~1050 BC and~700 BC) known as the Levantine Iron Age Anomaly, LIAA [35,36]. This extreme anomalous event has been considered as a geomagnetic spike [35] and has attracted the attention of the geomagnetic community. A great number of studies have been published in the last decade on the LIAA extension and evolution (e.g., [37][38][39][40][41][42]).
Using the last global paleoreconstruction that spans this time interval, the SHAWQ-Iron Age, we estimated the location of the ED from 1300 BC to 0 (Figure 4a). It can be seen it moves from the Earth's center towards the 20 • E, reaching the greatest eccentricity around 1025 BC, and comes back to the origin following a meridian plane. For this time window, the vertical component of the ED is positive. That is, the ED departs toward the Northern Hemisphere, being its maximum distance to the Earth's center at about 570 km. The normalized total energy, calculated up to maximum degree 10 for the SHAWQ-Iron Age, shows the two peaks in the intensity of the geomagnetic field at~950 BC and~775 BC, reflecting the characteristics of the LIAA. Please note that our ring model did not have enough resolution to reproduce in detail the two intensity peaks of the LIAA. Therefore, we considered only one intensity peak around 975 BC. In addition, considering the temporal limitation of the SHAWQ-Iron Age model, we do not believe that there is enough data to constrain the evolution of the LIAA in its earlier stages, so we focused on its final phase. Therefore, the 360-DRM was designed to catch the return to the center along the 20 • E meridian, that is, to model the ED displacement at the end of the LIAA.
Since the effect of this anomaly is an increase in the intensity of the field in a localized region, in the 360-DRM we selected a sector of dipoles from 10 • E to 30 • E that increase their magnitude with 2 dipoles per 100 years increasing between 1300 BC and 950 BC and a rate of 1 changed dipole per 100 years in the recovery phase (between 950 BC and 0). In this case, the 360-DRM starts with several dipoles increased at 1300 BC since at that time the SHAWQ-Iron Age presents a normal flux patch (NFP) at the CMB according to Osete et al. [23]. The ED path from the 360-DRM seems to reproduce the recovery of the Earth's magnetic field from the LIAA event. Please note that as a simple model with equatorial antisymmetry, our model was unable to reconstruct the details of the more complex ED displacement suggested by the SHAWQ-Iron Age, such as the displacement towards the Northern Hemisphere. Furthermore, the intensity of anomalous dipoles should only be considered qualitatively. Since the effect of this anomaly is an increase in the intensity of the field in a localized region, in the 360-DRM we selected a sector of dipoles from 10° E to 30° E that increase their magnitude with 2 dipoles per 100 years increasing between 1300 BC and 950 BC and a rate of 1 changed dipole per 100 years in the recovery phase (between 950 BC and 0). In this case, the 360-DRM starts with several dipoles increased at 1300 BC since at that time the SHAWQ-Iron Age presents a normal flux patch (NFP) at the CMB according to Osete et al. [23]. The ED path from the 360-DRM seems to reproduce the recovery of the Earth's magnetic field from the LIAA event. Please note that as a simple model with equatorial antisymmetry, our model was unable to reconstruct the details of the more complex ED displacement suggested by the SHAWQ-Iron Age, such as the displacement towards the Northern Hemisphere. Furthermore, the intensity of anomalous dipoles should only be considered qualitatively.

The South Atlantic Anomaly (SAA)
Probably the most important anomaly of the present-day geomagnetic field is the South Atlantic Anomaly (SAA), related to a low geomagnetic field intensity at the Earth's surface spanning over the South Atlantic between the South American and African continents. Its emergence is linked to the presence of a reverse flux patch (RFP) at the CMB below Africa (e.g., [43,44]) around 1400 AD (although there are studies, such as [23,45], that suggest its presence at least since 950 AD) that migrated to the southwest during the following 500 years. Nowadays, the original RFP seems to be divided into two different patches: one below South America, which is vanishing, and another one in Africa, which is reinforcing [46]. The presence of this anomaly has been discussed as one of the main causes of the observed decay in the geomagnetic field intensity at a global scale (e.g., [47,48]), and even has been studied as a possible precursor of a geomagnetic polarity transition [49,50].

The South Atlantic Anomaly (SAA)
Probably the most important anomaly of the present-day geomagnetic field is the South Atlantic Anomaly (SAA), related to a low geomagnetic field intensity at the Earth's surface spanning over the South Atlantic between the South American and African continents. Its emergence is linked to the presence of a reverse flux patch (RFP) at the CMB below Africa (e.g., [43,44]) around 1400 AD (although there are studies, such as [23,45], that suggest its presence at least since 950 AD) that migrated to the southwest during the following 500 years. Nowadays, the original RFP seems to be divided into two different patches: one below South America, which is vanishing, and another one in Africa, which is reinforcing [46]. The presence of this anomaly has been discussed as one of the main causes of the observed decay in the geomagnetic field intensity at a global scale (e.g., [47,48]), and even has been studied as a possible precursor of a geomagnetic polarity transition [49,50].
We analyzed the path followed by the ED from 700 AD to the present day (Figure 5a) by using the SHAWQ2k paleoreconstruction [23] until 1900 AD and then the IGRF-13 model [51]. At 700 AD the ED starts to move away from the Earth's center westward and at 1100 AD turns to the Pacific sector. When the RFP is below southern Africa (~1200-1400 AD, [23]), the ED is found between 210 • E and 230 • E; later, when the RFP migrates westward throughout the South Atlantic, the ED moves eastward along the Pacific Hemisphere and it is located around 140 • E at the present day (in 2020 AD, according to IGRF-13, a black star in Figure 5a). This observed evolution of the ED against the displacement of the SAA was also reported by Domingos et al. [52,53].
at 1100 AD turns to the Pacific sector. When the RFP is below southern Africa (~1200-1400 AD, [23]), the ED is found between 210° E and 230° E; later, when the RFP migrates westward throughout the South Atlantic, the ED moves eastward along the Pacific Hemisphere and it is located around 140° E at the present day (in 2020 AD, according to IGRF-13, a black star in Figure 5a). This observed evolution of the ED against the displacement of the SAA was also reported by Domingos et al. [52,53]. Figure 5. Evolution of the ED and the normalized total energy for the South Atlantic Anomaly period according to (a) the SHAWQ2k paleoreconstruction (position in km) and (b) the 360-DRM proposed in this paper, reversing and multiplying of the dipoles between 300° E and 0° E. Gray and black stars correspond to the ED location at 2000 and 2020 AD, respectively, calculated with the IGRF-13. Orange/green circumference indicates the dipoles that change their polarity (the first one is marked with a black triangle), the propagation directions (black arrows), and the limits of the reversing regions (black lines). Note that in (b) the evolution of the ED is dimensionless. Coastlines are plotted only for longitude reference (note they are projected on a sphere of arbitrary radius). As our simple model is not capable of modeling in detail the two meridian planes followed by the ED associated with the evolution of the SAA, we chose to model the average migration of the ED. We consider that in each step of the 360-DRM two dipoles reverse (one for the first step) in the sector of 300° E-0° E, starting from the position 330° E. The rate of reversed dipoles is 5 per 100 years. The SHAQW2k and the 360-DRM presented consistent results in the decreasing trend of the total field energy for the last millennium ( Figure 5, last column), although the decay in the ring model was slightly overestimated. The normalized total energy was calculated up to a maximum degree of 10 in both cases.

Discussion
Two types of anomalous events were studied in this work: extreme events related to polarity transitions (last reversal and excursions) and Holocene geomagnetic field anomalies (LIAA and SAA). The path followed by the ED during these events ( Figure 6) seems to be constrained to a band around the 0°-180° E meridian (±30°). This apparently simple behavior is what inspired us to use de 360-DRM to model qualitatively the transitional paths and the Holocene anomalies. Figure 5. Evolution of the ED and the normalized total energy for the South Atlantic Anomaly period according to (a) the SHAWQ2k paleoreconstruction (position in km) and (b) the 360-DRM proposed in this paper, reversing and multiplying of the dipoles between 300 • E and 0 • E. Gray and black stars correspond to the ED location at 2000 and 2020 AD, respectively, calculated with the IGRF-13. Orange/green circumference indicates the dipoles that change their polarity (the first one is marked with a black triangle), the propagation directions (black arrows), and the limits of the reversing regions (black lines). Note that in (b) the evolution of the ED is dimensionless. Coastlines are plotted only for longitude reference (note they are projected on a sphere of arbitrary radius). As our simple model is not capable of modeling in detail the two meridian planes followed by the ED associated with the evolution of the SAA, we chose to model the average migration of the ED. We consider that in each step of the 360-DRM two dipoles reverse (one for the first step) in the sector of 300 • E-0 • E, starting from the position 330 • E. The rate of reversed dipoles is 5 per 100 years. The SHAQW2k and the 360-DRM presented consistent results in the decreasing trend of the total field energy for the last millennium ( Figure 5, last column), although the decay in the ring model was slightly overestimated. The normalized total energy was calculated up to a maximum degree of 10 in both cases.

Discussion
Two types of anomalous events were studied in this work: extreme events related to polarity transitions (last reversal and excursions) and Holocene geomagnetic field anomalies (LIAA and SAA). The path followed by the ED during these events ( Figure 6) seems to be constrained to a band around the 0 • -180 • E meridian (±30 • ). This apparently simple behavior is what inspired us to use de 360-DRM to model qualitatively the transitional paths and the Holocene anomalies.

Evolution of the Eccentric Dipole during the Last Millennia. Average Position and Characteristic Times
Olson and Deguen [54] studied the geomagnetic dipole eccentricity from geodynamo models trying to find out which of the boundary conditions control the eccentricity of the geomagnetic field: at the CMB or at the ICB. According to their simulation results, they concluded that lopsided inner core growth produces persistent eccentricity of the geomagnetic field, with the best-fitting dipole axis offset in the direction of fastest inner core solidification. Seismic anisotropy studies estimated the best fitting growth direction of the inner core from 120 • E to 300 • E, placing the regions of fastest and slowest growth under the Banda Sea (Indonesia) and Brazil, respectively [55].
When analyzing the historical archaeomagnetic and paleomagnetic models, Olson and Deguen [54] concluded that the eccentric dipole has moved into Earth's Eastern Hemisphere over the past two centuries. In contrast, the average location of the ED for the last 10 kyr had an offset to the west according to CALS3k.4b [56] and CALS10k.1b [57] paleoreconstructions. Previously (last 5 Myr), it was located in the Eastern Hemisphere. From these results, they concluded that faster growth may have occurred in the Western Hemisphere of the inner core during the past ten millennia and that a reorientation of the location of the fastest inner core growth has occurred recently. A phenomenon that seems to have also occurred over the past 5 million years. models trying to find out which of the boundary conditions control the eccentricity of the geomagnetic field: at the CMB or at the ICB. According to their simulation results, they concluded that lopsided inner core growth produces persistent eccentricity of the geomagnetic field, with the best-fitting dipole axis offset in the direction of fastest inner core solidification. Seismic anisotropy studies estimated the best fitting growth direction of the inner core from 120° E to 300° E, placing the regions of fastest and slowest growth under the Banda Sea (Indonesia) and Brazil, respectively [55].
When analyzing the historical archaeomagnetic and paleomagnetic models, Olson and Deguen [54] concluded that the eccentric dipole has moved into Earth's Eastern Hemisphere over the past two centuries. In contrast, the average location of the ED for the last 10 kyr had an offset to the west according to CALS3k.4b [56] and CALS10k.1b [57] paleoreconstructions. Previously (last 5 Myr), it was located in the Eastern Hemisphere. From these results, they concluded that faster growth may have occurred in the Western Hemisphere of the inner core during the past ten millennia and that a reorientation of the location of the fastest inner core growth has occurred recently. A phenomenon that seems to have also occurred over the past 5 million years. Figure 6. Thin lines represent the geomagnetic eccentric dipole (ED) locations in the equatorial plane relative to Earth's center (a) during the last reversal and excursions and for the last 100 kyr according to the IMMAB4 (blue), IMOLEe (orange), and LSMOD.2 (green), and GGF100k (purple) paleoreconstructions, respectively, and (b) during the Holocene according to SHADIF14k (red), CALS10k2 (green), SHAWQ2k (yellow), and SHAWQ-Iron Age (purple) paleoreconstructions. Triangles are the time-average location of the ED for each paleoreconstruction calculated as the mean position in each axis. X and Y correspond to the geocentric cartesian axes contained in the equatorial plane and pointing to 0° and 90° E longitude, respectively. Note the different scales between both figures. Coastlines are plotted only for longitude reference (note they are projected on a sphere of arbitrary radius).
In Figure 6a, the ED paths and their time-average position are represented for the last reversal and last excursions studied in the Results section. The GGF100k paleoreconstruction that covers the last 100 kyr [58] is also plotted. It seems that transitional paths are constrained along the 0°-20° E/180°-200° E meridian. Figure 6b shows the ED evolution in the equatorial plane during the Holocene according to different paleomagnetic reconstructions as well as its time-average location for the Holocene (SHA.DIF.14k [24] and CALS10k.2 [59]); for the Iron Age, the first millennium BC (SHAWQ-Iron Age [22]); and for the last 2 kyr (SHAWQ2k [23]). As expected, the mean position of the ED for a time window of 10 kyr was closer to the rotation axis than when studied for a shorter time window with the SHAWQ2k and the SHAWQ-Iron Age. In Figure 6a, the ED paths and their time-average position are represented for the last reversal and last excursions studied in the Results section. The GGF100k paleoreconstruction that covers the last 100 kyr [58] is also plotted. It seems that transitional paths are constrained along the 0 • -20 • E/180 • -200 • E meridian. Figure 6b shows the ED evolution in the equatorial plane during the Holocene according to different paleomagnetic reconstructions as well as its time-average location for the Holocene (SHA.DIF.14k [24] and CALS10k.2 [59]); for the Iron Age, the first millennium BC (SHAWQ-Iron Age [22]); and for the last 2 kyr (SHAWQ2k [23]). As expected, the mean position of the ED for a time window of 10 kyr was closer to the rotation axis than when studied for a shorter time window with the SHAWQ2k and the SHAWQ-Iron Age.
It must be also pointed out that although the average positions of the ED for the Holocene, calculated from CALS10k.2 and SHA.DIF.14k, did not differ by more than 170 km and were statistically indistinguishable, and the east-west hemisphere in which they are located was not the same. Both ED averages were in the North Hemisphere. Although it must be taken into account that the two paleoreconstructions are not based on the same type of data (thermoremanent in the case of the SHA.DIF.14k and including sedimentary data in CALS10k.2) and this could be the origin of the small discrepancy in the location of the ED mean position, it is perhaps more prudent to conclude that the paleomagnetic information is not yet accurate enough for the calculated ED averages to be significant (please note the error bars include in both cases the Earth's center).
However, when using smaller time averages the eccentricity of the geomagnetic field is more evident. During the last 3 kyr the mean location of the ED moved from the European Hemisphere towards the Pacific Hemisphere biased by the presence of two major geomagnetic anomalies that occurred during this period: the LIAA and the SAA. The ED was located around 50 • E and 300 km away from the Earth's center during the first millennium BC according to the SHAWQ-Iron Age. The most eccentric position for the ED average was obtained during the last 2 millennia by SHAWQ2k paleoreconstruction: it was located at 330 km from the Earth's center and deviated towards the Pacific Hemisphere, around 200 • E meridian. It is important to note that SHAWQ2k and SHAWQ-Iron Age have a shorter time span than the SHA.DIF.14k and CALS10k.2, which affects the mean position of the ED. In both cases the mean position is biased by the presence of the anomalies.
González-López et al. [60] analyzed the secular variation at the millennial and centennial time scale searching for characteristic periods during the last 10 kyr by using different paleoreconstructions. They found characteristic times in the equatorial dipole and quadrupole terms (coefficients that control the equatorial location of the ED approximation) of around 1000-1400 years according to the SHA.DIF.14k, CALS10k.2, and SHAWQ2k paleoreconstructions. A wavelet analysis applied to the evolution of the ED as its distance from the Earth's center (Rc), the distance from its equatorial projection to the Earth's center (Sc), and the vertical displacement (Z) for the SHA.DIF.14, CALS10k.2, SHAWQ2k, and SHAWQ-Iron Age paleoreconstructions can be found in the Supplementary Material ( Figure S4).
If the position of the ED is related to the anisotropy in the ICB as suggested by Olson and Deguen [54], these analyses would suggest that there are variations in the ICB conditions with characteristic periods of 500, 1000-1500, and 3000 years. Obviously, further investigations are needed to verify this fact.

The ED as an Indicator of Anomalous Regions at the CMB. Interpretation with the 360-DRM
Olson and Deguen [54] also modeled a lopsided ICB dynamo that spanned a dipole collapse event that resulted in a polarity change. Before the reversal, the dipole intensity was high, the dipole offset was relatively small, and the dipole axis was confined to the cold Hemisphere of the core. When the dipole collapsed, the transitional path followed by the ED in the simulation was not so clear as the ones we observed from the paleoreconstructions (Figures 2a, 3a,b, 4a, 5a and S3), but it seems that the ED drifted far into the opposite hemisphere during the transition and before the field assumed a more complex multipolar structure. Following the polarity change, the dipole intensity was recovered, the new polarity was stabilized and the dipole axis returned to the cold hemisphere.
One of the most striking results of this study is the coherence of the paths followed by the ED during the polarity transition and the excursions studied (Figures 2a, 3a,b and S3). These are located very close to the meridian 0 • E-180 • E. During the Matuyama-Brunhes transition, according to the 360-DRM, the dipole that started the reversal (changing from reversed to normal polarity) was located at 190 • E meridian (Pacific Hemisphere), and during the last excursions the first dipole that reversed the polarity (from normal to reversed polarity) appeared in the opposite meridian (the 0 • meridian, Atlantic Hemisphere).
It is important to note that reversed dipoles will originate anomalous flux patches at the CMB centered in the location where the first dipole reversed. Following the simple 360-DRM, these anomalous flux patches would grow more or less in situ until they cover the entire CMB and the reversal is stabilized. In the case of the excursion the region of abnormal flux grows to a certain point, stops growing, and subsequently vanishes.
Considering the evolution of the ED as reflecting the production and gathering of flux patches at the CMB within preferential hemispheres was already proposed by Gallet et al. [18]. They also related the archaeomagnetic jerks to episodes of maximum geomagnetic field hemispheric asymmetry (the most eccentric events). Part of our work also analyzed the most eccentric episodes that occurred in the last 3.5 kyr, which seem to be related to the evolution of the most important anomalies of the geomagnetic field: the LIAA and the SAA. The major difference between the two anomalies is the sign associated with each of them. While the SAA has been related with a set of reversing dipoles located in the Atlantic Hemisphere (Figure 5b), the LIAA has been associated with a set of anomalous high-moment dipoles (Figure 4b). It is important to note that the ED tends to move towards the region where the intensity of the Earth's magnetic field is higher. That is, toward the anomalous region if the anomaly is in the sense of the present dipole field, as the LIAA, and towards the opposite hemisphere if the event is associated to reversing dipoles, as it is the case of the SAA. We cannot reproduce in detail the ED migration associated to these anomalies with our simple ring model. However, it is important to underline that the anomalous dipoles are restricted in both cases to one quadrant: between 200 • E and 30 • E, that is, the African-Atlantic region.
The search for observables that relate deviation from the GAD of the Earth's magnetic field and its relation to the lateral variations in heat flux at the CMB has been an active topic of research (e.g., [61][62][63][64][65], among others). In this context, a simple approximation to the geomagnetic field, as the ED, can help to enhance the Earth's magnetic field lateral asymmetries. Our result showed that ED paths during anomalous Holocene events or during polarity reversal or excursions correlate well to the regions where the heat flux is low (e.g., [66]), as evidenced in Figure 7.
abnormal flux grows to a certain point, stops growing, and subsequently vanishes.
Considering the evolution of the ED as reflecting the production and gathering of flux patches at the CMB within preferential hemispheres was already proposed by Gallet et al. [18]. They also related the archaeomagnetic jerks to episodes of maximum geomagnetic field hemispheric asymmetry (the most eccentric events). Part of our work also analyzed the most eccentric episodes that occurred in the last 3.5 kyr, which seem to be related to the evolution of the most important anomalies of the geomagnetic field: the LIAA and the SAA. The major difference between the two anomalies is the sign associated with each of them. While the SAA has been related with a set of reversing dipoles located in the Atlantic Hemisphere (Figure 5b), the LIAA has been associated with a set of anomalous high-moment dipoles (Figure 4b). It is important to note that the ED tends to move towards the region where the intensity of the Earth's magnetic field is higher. That is, toward the anomalous region if the anomaly is in the sense of the present dipole field, as the LIAA, and towards the opposite hemisphere if the event is associated to reversing dipoles, as it is the case of the SAA. We cannot reproduce in detail the ED migration associated to these anomalies with our simple ring model. However, it is important to underline that the anomalous dipoles are restricted in both cases to one quadrant: between 200° E and 30° E, that is, the African-Atlantic region.
The search for observables that relate deviation from the GAD of the Earth's magnetic field and its relation to the lateral variations in heat flux at the CMB has been an active topic of research (e.g., [61][62][63][64][65], among others). In this context, a simple approximation to the geomagnetic field, as the ED, can help to enhance the Earth's magnetic field lateral asymmetries. Our result showed that ED paths during anomalous Holocene events or during polarity reversal or excursions correlate well to the regions where the heat flux is low (e.g., [66]), as evidenced in Figure 7. If the equivalence between ring dipoles and outer core convective columns is considered [19], these results support that extreme geomagnetic events are controlled by the CMB thermal flux heterogeneities (e.g., [67,68], among others).
Large heat flux is expected to concentrate the geomagnetic field (e.g., [69]), and this could be envisaged as magnetic dipoles that are anchored in regions where heat flux is high. In contrast, the anchoring of the dipoles is less in the regions of low heat flux. During reversals or excursions, the magnetic field is destabilized in the regions of low heat flux. If the equivalence between ring dipoles and outer core convective columns is considered [19], these results support that extreme geomagnetic events are controlled by the CMB thermal flux heterogeneities (e.g., [67,68], among others).
Large heat flux is expected to concentrate the geomagnetic field (e.g., [69]), and this could be envisaged as magnetic dipoles that are anchored in regions where heat flux is high. In contrast, the anchoring of the dipoles is less in the regions of low heat flux. During reversals or excursions, the magnetic field is destabilized in the regions of low heat flux. This causes, as the simple 360-DRM indicates, the ED location to move in the opposite direction to the region in which the inversion began.
Other more complex explanations have been proposed to interpret the extreme events of the Earth's magnetic field ( [68,70], among others). For example, Terra-Nova et al. [71] considered that low CMB heat flux leads to core fluid upwelling [72] which disperses the magnetic field and causes events like SAA. However, this type of mechanism would not explain the existence of spikes (as the LIAA).
In any case, the ED evolution provides valuable information about the regions where changes start in polarity transitions and the location of the sources of important anomalies of the geomagnetic field. Although it is important to remember that the 360-DRM does not try to reproduce the quantitative characteristics of the geomagnetic field, but the qualitative ones (in terms of field geometry).

Conclusions
The evolution of the eccentric dipole was analyzed during different extreme events of the geomagnetic field using paleoreconstructions: the last reversal (~780 ka), last excursions (~41 ka and~34 ka), and two important Holocene anomalies: the Levantine Iron Age Anomaly (~1000 BC) and the South Atlantic Anomaly (analyzed from~700 AD to present day).
We observed two different behaviors of the ED for the reversal and the excursion events. During the Matuyama-Brunhes transition, the evolution can be separated into two parts: before the change of polarity the ED moves towards the Atlantic sector between 0 • E and 30 • E; after the transition, it comes back to the Earth's center at opposite longitudes, 180 • E-200 • E. However, during the Laschamp and Mono Lake excursions, the ED has its movement confined to a longitude frame around 180 • E, before and after the energy minimum. Despite this difference, both reversal and excursions presented a coincidence in the meridian plane followed by the ED path (0 • -180 • E).
The two Holocene anomalies studied have different nature and origin, and these are reflected in the position of the ED. While the LIAA is related to a local enhancement of the geomagnetic field (reproduced by our simple 360-DRM as a local increase of the dipoles) and the ED moves towards the anomalous region, during the SAA, which is modeled with reversing dipoles, the ED moves opposite to the anomaly.
Despite its simplicity, the 360-DRM is able to reproduce the ED path based on the paleoreconstructions information. Following the approximation that the dipoles could be related to convective columns in the outer core of the Earth, our simple model suggests that last reversal, last excursions, and Holocene anomalies were originated in the CMB regions where the heat flux is low.
Supplementary Materials: The following are available online at https://www.mdpi.com/article/ 10.3390/geosciences11110438/s1, Figure S1: evolution of the eccentric dipole according to the 360dipole ring model proposed in this work for studying the Matuyama-Brunhes transition (whole time interval); Figure S2: Power spectra of the Gauss coefficients obtained with the 360-dipole ring model in all cases of study; Figure S3: Evolution of the eccentric dipole and the normalized total energy for the Mono Lake excursion; Figure S4: Wavelet analysis for the time evolution of the eccentric dipole using paleoreconstructions that cover the last 10 kyr.  Innovación (CGL2017-92285-EXP, CGL2017-87015-P, PID2020-117105RB-I00, PID2020-113316GB-I00). We are very grateful for the contribution from two anonymous reviewers and the academic editor of this journal that have helped us to improve this work.