The Mag.num core field model as a parent for IGRF-13, and the recent evolution of the South Atlantic Anomaly

We present the GFZ candidate field models for the 13th\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$13{\mathrm{th}}$$\end{document} Generation International Geomagnetic Reference Field (IGRF-13). These candidates were derived from the Mag.num.IGRF13 geomagnetic core field model, which is constrained by Swarm satellite and ground observatory data from November 2013 to August 2019. Data were selected from magnetically quiet periods, and the model parameters have been obtained using an iteratively reweighted inversion scheme approximating a robust modified Huber norm as a measure of misfit. The root mean square misfit of the Mag.num.IGRF13 model to Swarm and observatory data is in the order of 3–5 nT for mid and low latitudes, with a maximum of 44 nT for the satellite east component data at high latitudes. The time-varying core field is described by order 6 splines and spherical harmonic coefficients up to degree and order 20. We note that the temporal variation of the core field component of the Mag.num.IGRF13 model is strongly damped and shows a smooth secular variation that suits well for the IGRF, where secular variation is represented as constant over 5-year intervals. Further, the external field is parameterised by a slowly varying part and a more rapidly varying part controlled by magnetic activity and interplanetary magnetic field proxies. Additionally, the Euler angles of the magnetic field sensor orientation are co-estimated. A widely discussed feature of the geomagnetic field is the South Atlantic Anomaly, a zone of weak and decreasing field strength stretching from southern Africa over to South America. The IGRF and Mag.num.IGRF13 indicate that the anomaly has developed a second, less pronounced eastern minimum at Earth’s surface since 2007. We observe that while the strong western minimum continues to drift westwards, the less pronounced eastern minimum currently drifts eastward at Earth’s surface. This does not seem to be linked to any eastward motion at the core–mantle boundary, but rather to intensity changes of westward drifting flux patches contributing to the observed surface field. Also, we report a sudden change in the secular variation measured at two South Atlantic observatories around 2015.0, which occurred shortly after the well-known jerk of 2014.0.


Introduction
The observed geomagnetic field is a superposition of the dominating field generated in Earth's outer core (the main field), remanent fields originating in the lithosphere, various contributions from electric current systems in the ionosphere and magnetosphere, and secondary fields that are induced mainly in the mantle. Geomagnetic core field models serve a wide variety of purposes, from practical applications such as navigation to scientific applications such as studies of core dynamo processes deep inside the Earth (e.g., Hulot et al. 2015). A widely used model in civil and scientific applications is the International Geomagnetic Reference Field (IGRF) that is updated every 5 years by the International Association of Geomagnetism and Aeronomy (IAGA). Here, we introduce the Mag. num.IGRF13 core field model (thereafter called Mag. num) that serves as parent for GFZ's IGRF-13 candidate Open Access *Correspondence: rother@gfz-potsdam.de 1 GFZ German Research Centre for Geosciences, Telegrafenberg, 14473 Potsdam, Germany Full list of author information is available at the end of the article models, i.e. for a definitive model for epoch .0 (DGRF 2015, a model for 2020.0 (IGRF 2020), and a prediction of the secular variation from 2020 to 2025 (SV). Data selection and processing are described in the data section, followed by the description of the modelling method and how the IGRF candidates were obtained from the parent model.
The South Atlantic Anomaly (SAA) is an area of low geomagnetic field intensity, strong westerly declination, and of a complex spatial pattern of inclination. The SAA has been developing since at least the 1600s in the South African region and has moved westward over the South Atlantic (Mandea et al. 2007;Hartmann and Pacca 2009). Its minimum is currently located over South America (Terra-Nova et al. 2017). The SAA is caused by magnetic flux patches at the core-mantle boundary (CMB) with magnetic flux opposite to the flux of the dominating dipole field direction in the Southern hemisphere.
The growth of the SAA might play an important role for the ongoing dipole decay (e.g. Bloxham and Gubbins 1985;Finlay et al. 2016a) and the magnetic field's shielding is clearly reduced in that area (e.g. Heirtzler et al. 2002;Domingos et al. 2017). Some studies propose that it has an influence on climate and sea level change through an increased inflow of radiation energy in this weak field area (De Santis et al. 2012;Campuzano et al. 2018). Moreover, the SAA has been suggested as indicator for an impending field reversal or excursion (Gubbins 1987;De Santis et al. 2013;Tarduno et al. 2015). Note, however, that recent palaeomagnetic modelling results indicate that similar weak field structures in the vicinity of the South Atlantic occurred several times in the past 100,000 years and did not directly lead to excursions (Shah et al. 2016;Brown et al. 2018;Panovska et al. 2019). Here, we examine the most recent geomagnetic field development in the SAA region based on geomagnetic observatory data and the Mag.num model.

Satellite data selection
Level 1B data of the three Swarm satellites were included from the beginning of the mission on 27th November 2013, up to 3rd August 2019, using baseline 05 as provided by ESA. These data are provided at 1 Hz temporal resolution, and have been sub-sampled to intervals of at least 20 seconds. In the following, low latitude refers to regions between ± 55 solar magnetic (SM) latitude. For these regions, data are rotated into the SM coordinate system. For the complementary high-latitude regions, data are expressed in the Earthcentred Earth-fixed North, East, Centre (NEC) coordinate system. We apply the following selection criteria: • The z-component of the interplanetary magnetic field (IMF-B z ) is positive. • The large-scale magnetospheric field MMA_SHA_2F (Hamilton 2013) Swarm Level 2 product and its time derivative are lower than 40 nT and 40 nT/day, respectively. • Various quality flags indicate valid magnetic field and attitude information.
Additionally, only nighttime data are selected for low-latitude regions such that • Local time is between 23:00 and 05:00.
• The sun is below the horizon at 100 km above Earth's reference radius of 6371.2 km, corresponding to ionospheric altitudes.

Satellite difference data
Apart from the satellite vector measurements themselves, the along-track and cross-track differences between Swarm A and C can be used as a source of magnetic field information. External field contributions are notably reduced in these data, and we take advantage of this fact to construct a starting model for the iterative modelling process for Mag.num. As described in the method section below, the coefficients of the final Mag.num model are obtained by inverting solely the vector field measurements themselves, not their differences.
The velocity of the Swarm satellites is sufficiently constant to calculate along-track differences for all three satellites using pairs of data with a constant temporal separation of 20 seconds. However, the estimation of cross-track differences requires a variable temporal separation to keep the inter-satellite distance constant. For a given data point i and the chosen distance D 0 , the paired data point j is found by minimising the following functional φ: where f = 5 is a weighting factor, and D ij and dt ij are the actual distance and temporal separation between the two data points. Choosing D 0 = 160 km leads to a time separation of about 20s.
In total, we used about 5.2 million cross-track and along-track difference data. Outliers caused by data gaps have been removed, leading to a rejection of 12% of the original data set.

Observatory data selection
Additionally, we use observatory hourly mean values (HMV), spanning the period between January 2013 and July 2019, to better constrain secular variation. These data are distributed by ESA as Swarm AUX_OBS (baseline 01, version 20) (Macmillan and Olsen 2013), and are based on data stored at INTERMAGNET (https ://www. inter magne t.org/) and the World Data Centre Edinburgh (http://www.wdc.bgs.ac.uk/). Observatory data are selected to reduce external field contributions by applying the following criteria: • The Dst index (in form of the Swarm auxiliary data product AUX_Dst_2) is between ± 30 nT. • Local time is between 23:00 and 05:00.
• The sun is below the horizon at 100 km above Earth's reference radius of 6371.2 km.
Obvious individual outliers were manually removed after a few iterations, as soon as the observatory offsets and external field parameters have converged. The temporal coverage of the selected observatory HMV is shown in Fig. 1.

Parameterisation
In the source-free region, the magnetic field can be described as the negative gradient of a scalar potential V associated with sources of internal and external origin, i.e.
(2) B = −∇V i (θ , φ, r, t)−∇V e (θ, φ, r, t), where and Here, Y m l (θ , φ) are the Schmidt semi-normalised spherical harmonics (SH) of degree l and order m, θ , φ, r and a are the geocentric co-latitude, longitude, satellite altitude and model reference radius, respectively, and g m l and q m l are the Gauss coefficients. We use the convention that negative orders are associated with sin(|m|φ) terms; whereas, positive orders ( m ≥ 0 ) are associated with cos(mφ) terms.
The internal sources consist of the time-varying core magnetic field, the quasi-static lithospheric field, and induced fields. The spatial signals of the core and lithospheric fields overlap, and we assume that the core magnetic field's temporal evolution can be resolved up to SH degree and order L i = 20 (c.f. Eq. 3). This temporal evolution is parametrised using order 6 B-splines basis functions ψ 6 i (t) (de Boor 1978) which we define on the interval from December 31, 2012 to January 1st, 2020. This implies that the boundary nodes are not directly constrained by data. Overall, the time dependence of the Gauss coefficients is then given by where g m lj are a set of spline coefficients for each Gauss coefficient g m l . Here, we use a half-year knot separation, and N t = 9 spline functions are used to represent each Gauss coefficient. Further, we set the reference radius of this large-scale field to a = 3485 km , the radius of the outer core. This choice has no influence on the inversion as we do not apply a regularisation constraint on the Gauss coefficients in spectral domain. For the core and the static SH degrees 21 to 30, the reference radius is set to the mean Earth radius of a = 6371.2 km . We do not aim to model the lithospheric field, and we subtract the Swarm auxiliary lithospheric field model (AUX_LIT_2) provided by ESA for SH degrees from 25 to 80 from the satellite data to avoid aliasing from unmodelled short wavelength field. We further co-estimate the remaining constant offsets for observatory data, so called observatory biases (Langel et al. (1982); see also Mandea and Langlais (2002); Macmillan and Thomson (2003)), which account for the unmodelled static lithospheric field of smaller spatial scales (Additional file 1: Table S1). Note, however, that due to the previous subtraction of the lithospheric field model,  Table S1 links the observatory numbers to observatory names and locations these biases do not reflect the full deviations from the core field. Moreover, induced fields are parameterised with one SH coefficient (degree one and order zero), and assumed to be proportional to the induced part, MMA i , of the Swarm MMA index (Hamilton 2013). The corresponding coefficient g 0,MMA 1,j serves as a scaling factor, and is estimated in N e = 14 6-month intervals such that where H j (x) = x, t ∈ [t j , t j+1 ] and is zero otherwise. For observatory data, we use a similar parameterisation, but the I st index is used instead of the MMA i index. The I st index represents the induced part of the D st index, and is derived using a 1D mantle conductivity profile as described by Maus and Weidelt (2004).
The external field is separated into a slowly varying part which is assumed to be constant over 30-day time intervals, and a more rapidly varying part. Concerning the slowly varying part, it is separated into a component that is described in Geocentric Solar Magnetic Coordinates (GSM) and another component that is described in Solar Magnetic (SM) coordinates. Both components are assumed to be large scale, and we use a maximum SH degree l = 1 such that where R GSM and R SM are matrices rotating vectors defined in GSM and SM reference frames into the geocentric Earth fixed reference frame, respectively. Regarding the more rapidly varying part, it is separated into a component that scales with the Y-component of the Interplanetary Magnetic Field (IMF By) at a temporal resolution of 60 minutes, and another component that scales with the external part of the Swarm MMA index, MMA e , varying with the temporal resolution of one Swarm orbit, i.e. 90 minutes (Hamilton 2013). These components are completely described in SM coordinates, and given by , and q 0 MMA 2j are assumed to be constant on 30-day time intervals. The external field coefficients were estimated independently for satellite and observatory data. For the latter, we additionally impose that q 0 SM 1j = 0 to avoid colinearities with the observatory crustal offsets. Also, we use the E st index instead of the MMA e index, analogous to Eq. 6.
In addition, the rotation angles between magnetic field vectors in the sensor reference frame and the satellite coordinate system are co-estimated in bins of fixed length (Rother et al. 2013). A bin-size of 27 days gives a good trade-off between maximising goodness-of-fit and minimising temporal variability of the estimated Euler angles. At low-latitude regions (between ± 55 solar magnetic latitude, see "Satellite data selection" section), however, the number of data per bin is not sufficient, and the local time selection criterion has been loosened for the purpose of Euler angle estimation in order to obtain a sufficient number of data (Rother et al. 2013). Any remaining bins with less than 1000 data points are ignored.

Iterative inversion scheme
The relationship between magnetic field vector measurements and the Gauss coefficients is linear and can be expressed as where d is the data vector, m is the vector of Gauss coefficients, and G is the matrix linking the data to the Gauss coefficients according to Eqs. 2-8.
Due to unmodelled physical signals in the data, which will be considered as noise here, this equation cannot be solved directly. Instead, the Gauss coefficients are determined by minimising where W contains the estimated data weights w i . Here, we neglect any cross-correlations between data errors and assume that W is diagonal.
Data weights should ideally represent the noise associated with data. The least-square fit yields the maximum likelihood solution if the data noise follows a Gaussian distribution. Here, we account for non-Gaussian noise distributions with larger tails, i.e. larger data outliers, using a modified Huber norm (Morschhauser et al. 2014;Lesur et al. 2015). In practice, this norm is approximated using iteratively re-weighted least-squares (Farquharson and Oldenburg 1998), and the weight of data point d i at iteration j + 1 is given by where σ i corresponds to the assumed standard deviation (see Table 1), m j refers to the estimated Gauss coefficients at iteration j, and k i and a i are parameters that account for the expected variability among different data subsets for this norm. In particular, this norm is based on a Gaussian distribution for data residuals lower or equal than k i , and a i controls the shape of the underlying distribution's tail such that lower values will have a larger tail. These parameters are chosen such that they represent the distribution of residuals from a preliminary model using an L2 norm (Table 1).
Noise, i.e. unmodelled signals in the data, will lead to an unstable model which will show unrealistically strong spatial and temporal field variations when downward continued. Therefore, we regularise our model by damping the third time derivative of the radial field component at Earth's core radius c, i.e. by additionally minimising over the whole time interval from beginning ( t 1 ) to end ( t 2 ) of the model. Further, the temporal end points of the model are explicitly damped to avoid oscillations due to the lack of data near the end points, by minimising at times t = t 1 and t = t 2 , i.e. both the second and first derivatives of the radial field at the temporal start and end points of the model. The damping parameters are often chosen such that the data misfit corresponds to the estimated noise level. Here, we choose a rather high damping parameter to obtain a smoothed model that suits well to IGRF's constant secular variation within 5-year intervals, unlike a model that follows rapid SV variations more closely.
As described above, Mag.num is derived using a robust norm which is approximated using an iterative algorithm (Farquharson and Oldenburg 1998). Further, the estimation of Euler angles also introduces non-linearities, and requires an iterative inversion scheme (Rother et al. 2013). As starting model for the iterative process, we decided to use a model which entirely relies on Swarm satellite difference data (see data Section). The parameterisation of this Delta-model does not include any external fields as they are largely eliminated by considering differences only. As a result, the damping of the third time derivative (Eq. 12) was significantly reduced. Overall, the Delta model turned out to be a suitable and efficient starting model as it effectively reduces external field contributions and allows to reduce the number of required iterations.

IGRF-13 candidate models
The GFZ candidate model for DGRF-13 is a snapshot at epoch 2015.0 of the Mag.num parent model, truncated to degree and order 13. This epoch fully lies within the Mag.num time interval and is well covered by satellite and observatory data. Concerning the candidate for IGRF-13 at epoch 2020.0, data coverage did not reach the year 2020, still, we use a snapshot of Mag.num at 2020. This snapshot is a direct forecast which is considered relatively robust since temporal oscillations in the model have been strongly damped. Secular variation is not yet predictable by means of modelling. Therefore, we decided for a simple and conservative choice: the SV for 2020 to 2025 is approximated by the SV given by Mag. num at epoch 2019.0. This choice is justified by the good data coverage at this epoch. As an example, our secular variation prediction of g 0 1 and g 0 2 is shown in Fig. 2 (open red circle) along with the SV of the Mag.num parent

Table 1 Parameters used in data weighting
Sat and Obs stand for satellite and observatory data, LM indicates the low-and mid-latitude subset and HL the high-latitude subset (see "Satellite data selection" section). Parameters k and a refer to the applied modified Huber norm for the X, Y and Z component data (Eq. 11). The standard deviations σ (in nT) for the satellite data are given separately for Swarm satellites A, B and C

Sat
Obs For comparison, we also show SV as obtained from IGRF-12, with its predicted SV for 2015-2020 and the previous values calculated as linear 5-year SV from the main field coefficients (blue circles, see eq. 2 in Thébault et al. 2015). Mag.num is explicitly designed to describe the SV variabilities as they can be represented by IGRF. However, we note that the actual SV is more variable which is illustrated by the CHAOS model (grey line, Finlay et al. 2016b) that shows considerable non-linearities even within 5-year intervals. Our candidate models are compared to all other candidate models and the final IGRF-13 products by Alken et al. (2020a).

Geomagnetic power spectra
Geomagnetic power spectra for main field (MF), secular variation (SV) and secular acceleration (SA) are commonly used to compare different geomagnetic field models. In Fig. 3, the spectra from the Mag.num IGRF parent model for 2015.0 are shown in comparison to those of CHAOS6-x9 (Finlay et al. 2016). There is a close agreement between the models, both with respect to the MF and SV spectra. The much stronger regularisation of Mag.num, i.e. its weaker temporal variability, is reflected by a less energetic SA at low degrees (2-4) than that of CHAOS. The higher power at high SV and SA degrees (above 8) is of no concern as these degrees are not considered for the IGRF candidate models.

Data misfit and residuals
In Figs. 4 and 5, we show statistics and maps of the input satellite data residuals (after data selection and filtering) to the Mag.num IGRF parent model. The number of data, the mean of the residuals and the root mean square (rms) misfit between the data and the model are given in Table 2.
As seen in Fig. 4, the distributions of residuals are close to but not fully Gaussian, all with positive excess kurtosis between 0.95 and 2.26. All have a nearly zero mean with very small skewness, indicating that no systematic biases are present for any of the components or satellites. The mapping of residuals in Fig. 5 confirms that the data are fit similarly well over all longitudes. For low and mid latitudes, the satellite data misfit amounts to no more than 3.2 nT for all components (Table 2). Higher residuals at high geomagnetic latitudes are expected as the high-latitude input data are more strongly influenced by external field signals, which should not be fit by the core field model. The average misfit here ranges from 14.1 nT in Z to 44.1 nT in the Y component (Table 2), and there are slight biases towards negative values in the X and Z and towards positive values in the Y component. The absolute average misfit of the additional satellite data that have been used for the Euler angle determination is slightly larger, but in the same order of magnitude as the input data, also with nearly zero mean. These data have not been used for the inversion of the core field model, and hence may serve as an independent benchmark. For the observatory data, the average misfit at low and mid latitudes lies between 3.6 and 5.1 nT, reaching up to 12.9 nT at high latitudes ( Table 2). The mean and root mean square misfits for each observatory are given in Additional file 1: Table S1. SHE was installed in 2007, and data were recently re-calibrated to account for the effects of strong gradients in the crustal magnetic field on Saint Helena. SHE data are not contained in Swarm AUX_OBS, and hence were not used to construct the Mag.num model. We compare these two data sets with field predictions from Mag.num.

Comparison to independent observatory data
In particular, we calculate monthly mean values (MM) from calibrated hourly mean value (HMV) data of both observatories, and convert to geocentric coordinates. To approximate secular variation, finite differences are calculated based on these monthly means by where MM t+6 is the value 6 months after the time of interest and MM t−6 the value 6 months before. This was done to avoid seasonal effects. The resulting approximated secular variation is shown in Fig. 6 and indicated by blue dots. For low-latitude observatories such as SHE and TDC, all components, and particularly the X (geocentric North), are affected by the magnetospheric ring current. Therefore, we additionally subtract the monthly mean values of the ring current magnetic field as predicted by the CHAOS6-x9 magnetospheric field model (Finlay et al. 2016b) from the estimated SV t value. This results in much smoother estimates of secular variation (red dots). We note that a secular variation pulse (geomagnetic jerk) is observed in 2014 and is visible in the Y and Z components of both observatories. This jerk has been first reported by Torta et al. (2015), and was responsible for poor predictions of the SV in the IGRF-12 model. Additionally, a less prominent, but sudden change in secular variation is also observed in early 2015 for the (14) SV t = MM t+6 − MM t−6 ,   Figure 6 shows how Mag.num (black solid line) and CHAOS6-x9 (black dashed line) fit the data of both observatories from 2015 to 2020. The good fit between Mag.num and this independent observatory data set, not included in the construction of the model, serves as a validation of the Mag.num secular variation estimates over this interval of time. The fact that Mag.num, unlike CHAOS-x9, does not fit the more rapid variations, such as the 2014 jerk, is expected and in agreement with the model's design for low temporal variability.

At Earth's surface
The South Atlantic Anomaly (SAA) is characterised by a longitudinally elongated minimum at Earth's surface. The minimum of the SAA has been drifting westwards from the south-western coast of Africa in 1700 to its present position over Brazil (Hartmann and Pacca 2009). In recent years, the development of a second, eastern minimum within the SAA at Earth's surface has been observed at around 5 • W, 40 • S (see also Fig. 7) and it was first reported by Terra-Nova et al. (2019). Figure 7a, b illustrate the SAA and the two minima in field intensity F as predicted by Mag.num at Earth's surface for epochs 2013.0 and 2020.0, respectively. The (field intensity) isodynamic lines of the values at the saddle points (24.48 and 24.06 µ T, respectively) between the minima are marked in black. In addition, Fig. 7c shows the field for epoch 2020.0 at the reference satellite altitude of 450 km. The best-fit great circle arc connecting the time-averaged (from 2013.0 to 2020.0) locations of the two minima and of the saddle point at Earth's surface is also shown (red line).
Realistic model uncertainties cannot be derived directly from standard spherical harmonic models. To estimate the robustness of the second minimum and the reliability of the following analysis based on the Mag.num model, we compared the results from all IGRF candidate models for the locations of SAA minima and saddle point. Figure 8 shows close agreement among all models, with the different results falling within ±0.1 • of latitude and mostly ± 0.4 • of longitude from the average location for all three points. The Mag.num results lie well within this range in all cases. The sharp main minimum seems most robustly resolved, while uncertainties in longitude are somewhat larger for the flatter saddle point and second minimum, but still well separated in all models. Figure 7d shows the intensity profiles along the red line of panels a-c for the two reference altitudes and the two epochs, indicating the intensity evolution of the SAA. Due to geometric attenuation with distance from the source, i.e. the geodynamo in Earth's core, a plateau, but no clear second minimum, is observed at the typical satellite altitude of 450 km (see Fig. 7c, d). To study the evolution of the SAA over a longer time scale, we plot the longitudes of the SAA minima and of the saddle point between them as predicted by the IGRF model from 1900.0 and 2020.0 in the 5-year IGRF intervals (Fig. 9). The western, main minimum is depicted by cyan circles. The secondary minimum (light blue triangles) and the saddle point (green triangles) are not present before 2005.0. In addition, the corresponding predictions by Mag.num are shown by black, purple, and orange solid lines with an annual resolution for the epochs 2013.0 to 2020.0, and the magenta  Fig. 7a, b). The drifts of western minimum and saddle point most likely continue according to both the Mag.num and also the final IGRF-13 (Alken et al. 2020b) predictions, as shown by the red squares and blue crosses, respectively. The drift of the eastern minimum has already slowed down in recent years, and the models predict that it will not drift further eastward in the next 5 years. The field strength in the eastern minimum decreased faster, by about 0.5 µT over the 7 years from 2013 to 2020, compared to less than half of that in the western minimum (Fig. 7d). However, we note that with values around 23.5 both Mag.num and IGRF-13 predict values similar as today for the eastern minimum in 2025. The field intensity of the western minimum, on the other hand, lies at 22.25 µ T in 2020 according to Mag.num and will further decrease to about 22.1 µ T according to both Mag.num and IGRF-13.

At the core-mantle boundary
The SAA is associated with reverse flux patches (RFPs) at the core-mantle boundary (Fig. 10a), but Terra-Nova et al. (2017) showed that it is not straightforward to link RFPs and the SAA minimum. Their analysis suggests that the location of the surface intensity minimum is determined by the combined influence of reverse and normal flux patches. Here, we restrict ourselves to a qualitative analysis of the RFPs at the core-mantle boundary. Downward continued maps of the radial field at the CMB derived from Mag.num predictions show four-five RFPs (positive Z, reddish colours) in the area from 0 • to 60 • S and 80 • W to 20 • E (Fig. 10). With the exception of a nearly stationary RFP near the tip of South America (c.f. RFP 1 in Table 3, see also next section) and a similarly stationary strong normal flux patch at the West coast of Africa (c.f. NFP 1 in Table 3, see also next section) all of them have drifted westward by a few degrees over the 7 years. The eastern local minimum at Earth's surface (ca. 0 • N, 40 • S) is located above an area of weak normal flux at the CMB. Its apparent eastward drift seems to be caused mostly by the intensification of the RFP to its east (under southern Africa) compared to weaker changes of the RFP in the middle of the Atlantic and the normal patch south of that.

Linking the CMB and surface observations
Under the assumption of a source-free mantle, the surface field intensity expression of the South Atlantic Anomaly (SAA) is connected to the flux pattern at the CMB by the Laplace equation, and the mantle acts as a wavelength-dependent filter. Finlay et al. (2020) used an analytical method to identify the features at the CMB that relate to the observed SAA minima at the surface. Here, we follow a more qualitative approach by tracking the field minima of the SAA down to the CMB. For this purpose, we investigate the field morphology at different levels within the mantle as predicted by Mag.num. We note that this analysis does not provide any physical interpretation, but is rather a discussion of downward continuation of a potential field as described by the Laplace equation.
In Fig. 11, we show the field intensity (F) (panels a-c) and the vertical field component (Z) (panels d-f ) in 2020, both at levels of the Earth's reference surface radius (panels a and d), at 1600 km (b and e), and at 2400 km (c and f ) depth. Additional file 2: Figure S1 provides a more detailed view with a comparison of F in 2013 and 2020, and Z in 2020, all at levels of 0, 800, 1600, 2000, 2400 km and at the CMB depth of 2886.2 km. In maps showing F (Fig. 11a-c, Additional file 2: Figure S1), minima (blue) and maxima (orange) are labelled by numbers 1-7. In maps showing Z, patches are numbered 1-10, with normal flux patches (NFP, Z < 0 ) shown in green and blue, and reversed flux patches (RFP, Z > 0 ) shown in orange. Note that, as expected, some patches merge or disappear with increasing distance from the CMB. Table 3 lists information about the labelled extrema in F and flux patches in Z for the investigated depths, in particular including changes in their locations and strengths from 2013 to 2020.
Minimum 1 corresponds to the westward drifting main minimum of the SAA and is visible throughout the mantle. Looking at Z, it is linked to RFP 1 and its extension to the north. Minimum 2 corresponds to the eastward drifting, new eastern minimum in the SAA. It is also visible throughout the mantle, and grows to a northeast-southwesterly elongated zone of low F at 2400 km depth that also includes minimum 7. This zone is spatially linked to RFP 2, 7 and 8 at the CMB (see Additional file 2: Figure S1, the latter two have merged at 2400 km). We further note that RFP 5 and NFP 6 cancel out at 1600 km and above, leading to the rather flat area of F in the centre of our map at 1600 km depth (Fig. 11b). Interestingly, the new eastern minimum (2) is significantly larger in area and lower in F than the main minimum (1) at depths of 800-2400 km, i.e. throughout the mantle. Following table 3, minimum 2 reaches a value in F as low as 1363 nT at 2000 km depth, less than 10% of F at all other considered depth levels. Regarding the development of field strength of minimum 2 from 2013 to 2020 (c.f. Table 3), we observe a decrease in F at depths down to 2000 km (i.e. a deeper minimum in 2020 than 2013). Overall, F in minimum 2 decreased faster than in minimum 1, by at least a factor of 2 (dF in Table 3). Note, however, that the interpretation of weak F as an anomaly becomes less suitable with depth due to the increasing complexity of the field towards the CMB.
We expect that minimum 2 will probably either disappear with time, or it will grow (as it currently does) and then, controlled by the underlying westward drifting minimum, automatically turn into a westward drifting feature. The fast decrease of F in the eastern minimum 2 suggests that the latter scenario could be more likely. In this scenario, the region of anomalously low geomagnetic field strength stays in the South Atlantic region despite the general westward trend of the geomagnetic field, by a mechanism of formation of new minima of geomagnetic field strength at the eastern edge of the South Atlantic Anomaly.

Conclusions
We have presented Mag.num, the GFZ parent model for IGRF-13, derived using Swarm satellite and ground observatory data. The modelling method closely follows that of earlier developments (in particular the GRIMM models, see Lesur et al. 2008;Rother et al. 2013). A main modification is that an interim model based on Swarm along-and cross-track differences is used as a starting model for the iterative modelling process. This starting model has minimum external field influence and reduces the required number of iterations. Considering that the IGRF consists only of snapshots in 5-year intervals, we present a strongly regularised model version. This model does not reproduce rapid secular variation changes such as the jerk around 2014.0 and the sudden change in secular variation around 2015.0 that we report here for the Fig. 10 a Global map of the vertical field component downward continued to the core-mantle boundary from Mag.num for 2020.0. Local extrema of normal (red plusses) and reverse (red stars) flux are marked in the more detailed regional maps for 2013.0 (b) and 2020.0 (c) Fig. 11 Intensity (F, a-c) and vertical component (Z, d-f) for three depths, namely the reference Earth surface radius (0 km), 1600 km, and 2400 km. Several extrema in both field components are labelled and further information about them is given in Table 3. The black line in (a) is the contour line through the saddle point, the black line in (d) indicates the magnetic equator, and the red line in (a, d) indicates the great circle profile as shown in Fig. 7 South Atlantic observatories on the islands St. Helena and Tristan da Cunha. We chose the Mag.num parent model secular variation at epoch 2019.0 as estimate for the secular variation from 2020 to 2025.
Using the IGRF and Mag.num, we investigated the recent evolution of the SAA and in particular a second, less pronounced eastern minimum which has developed since 2007. A plateau is found in field intensity above the eastern minimum at satellite altitude, but due to the attenuation with increasing distance from the source no clear minimum is seen there. Contrary to the movement of the pronounced western minimum and the SAA as a whole, which have clearly and continuously drifted westward over the past decades, the eastern minimum has drifted slightly eastward at Earth's surface. However, an analysis of the magnetic flux distribution at the CMB revealed that all reverse flux patches drift slightly westward, and the apparent eastward drift at Earth's surface must be a result of intensity changes of different flux patches. Our investigation of magnetic field maps at 0, 800, 1600, 2000 2400 km depths and at the CMB for 2013 and 2020 provides a comprehensive description of the geomagnetic field variation in the SAA region between the CMB and Earth's surface, as well as in time. The secondary eastern minimum is found to be a relatively weak surface expression of a remarkable subsurface area of Table 3 Information about extrema in F and Z labelled in Fig. 11 and Additional file 2: Figure S1 Nr is number as given in the figures, and dlon and dlat are the temporal change from 2013 to 2020 in longitude and latitude, respectively (positive indicates eastward and northward, respectively). In column min/max it is indicated whether the extremum is a minimum or maximum in F, or refers to a reverse (RFP) or normal (NFP) flux patch in Z at the CMB. The values dF (given for 0-2400 km) or dZ (for the CMB) indicate the field change from 2013 to 2020, and the corresponding trend of increasing or decreasing field intensity is additionally shown by ( ↑/↓)