The CHAOS-7 geomagnetic field model and observed changes in the South Atlantic Anomaly

We present the CHAOS-7 model of the time-dependent near-Earth geomagnetic field between 1999 and 2020 based on magnetic field observations collected by the low-Earth orbit satellites Swarm, CryoSat-2, CHAMP, SAC-C and Ørsted, and on annual differences of monthly means of ground observatory measurements. The CHAOS-7 model consists of a time-dependent internal field up to spherical harmonic degree 20, a static internal field which merges to the LCS-1 lithospheric field model above degree 25, a model of the magnetospheric field and its induced counterpart, estimates of Euler angles describing the alignment of satellite vector magnetometers, and magnetometer calibration parameters for CryoSat-2. Only data from dark regions satisfying strict geomagnetic quiet-time criteria (including conditions on IMF \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$B_z$$\end{document}Bz and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$B_y$$\end{document}By at all latitudes) were used in the field estimation. Model parameters were estimated using an iteratively reweighted regularized least-squares procedure; regularization of the time-dependent internal field was relaxed at high spherical harmonic degree compared with previous versions of the CHAOS model. We use CHAOS-7 to investigate recent changes in the geomagnetic field, studying the evolution of the South Atlantic weak field anomaly and rapid field changes in the Pacific region since 2014. At Earth’s surface a secondary minimum of the South Atlantic Anomaly is now evident to the south west of Africa. Green’s functions relating the core–mantle boundary radial field to the surface intensity show this feature is connected with the movement and evolution of a reversed flux feature under South Africa. The continuing growth in size and weakening of the main anomaly is linked to the westward motion and gathering of reversed flux under South America. In the Pacific region at Earth’s surface between 2015 and 2018 a sign change has occurred in the second time derivative (acceleration) of the radial component of the field. This acceleration change took the form of a localized, east–west oriented, dipole. It was clearly recorded on ground, for example at the magnetic observatory at Honolulu, and was seen in Swarm observations over an extended region in the central and western Pacific. Downward continuing to the core–mantle boundary, we find this event originated in field acceleration changes at low latitudes beneath the central and western Pacific in 2017.


Introduction
The Earth's magnetic field is a fundamental part of our planetary environment and an integral component of many modern navigational systems, providing a natural and readily available source of orientation information. To make use of the geomagnetic field for navigation one requires a good-quality magnetometer to measure it, and a reference field model that relates the magnetic vector, at the location and time of the measurement, to the geographic directions. The International Geomagnetic Reference Field (IGRF) is a prominent example of such a reference model, and a trusted source of information on the Earth's magnetic field for the wider scientific community including space physicists, high-energy particle physicists, exploration geologists, engineers and biologists. This article reports on the parent model of candidates submitted by the Technical University of Denmark (DTU) for IGRF-13 in October 2019.
This parent model, called CHAOS-7, is the latest in a series of time-dependent geomagnetic field models developed at DTU over the past 15 years (Olsen et al. , 2014Finlay et al. 2016). CHAOS-7 spans a 21-year period from 1999 to 2020 for which both satellite and ground observatory data are available. For the past 6 years, satellite data have been delivered by the Swarm satellite trio, providing a particularly complete and homogeneous data coverage. Here, we take this opportunity to report in detail on two particularly intriguing aspects of recent geomagnetic field change. First we document changes since 2014 in the South Atlantic weak field anomaly which has important implications for the radiation dose experienced by satellites, and second we investigate patterns of rapid field change observed in the Pacific region over the past 6 years.
The South Atlantic Anomaly (SAA) was originally detected by early low-earth satellite missions in the late 1950s as a region of enhanced flux of energetic charged particles (Yoshida et al. 1960;Vernov and Chudakov 1960;Ginzburg et al. 1962). It has been well documented over the intervening years as a region of geospace where satellites systematically experience an enhanced radiation dose (e.g. Gledhill 1976;Heirtzler et al. 2002). The depth to which charged particles in the radiation belts penetrate, known as their bounce point along a field magnetic field line (Walt 2005), depends on the intensity of the geomagnetic field. Anomalously weak magnetic field in the South Atlantic (compared with the field of a centred dipole) thus gives rise to enhanced charged particle flux in this region, i.e. a radiation anomaly. By tracking the evolution of the weak magnetic field region we are therefore able to track the development radiation anomaly. Here, we document the evolution of the South Atlantic Anomaly as observed by the Swarm satellites, both in the magnetic field and in single-event electronic upsets monitored by onboard instruments. We go on to investigate changes in the field intensity mapped down to Earth's surface using CHAOS-7. Ultimately processes in the Earth's core determine the future evolution of the South Atlantic Anomaly; we use Green's functions relating the radial magnetic field at the core-mantle boundary to changes in the field intensity at satellite altitude to study the origin of these processes at the outer edge of the fluid outer core.
A second focus point in this article is rapid secular variation observed since the launch of Swarm in the Pacific region. Traditionally the Pacific has been thought of as a quiet region for secular variation (Vestine and Kahle 1966), but observatory records from Honolulu have indicated large amplitude changes in secular variation in the past decade (see the discussion in Finlay et al. (2016) and the section "Field acceleration changes in the Pacific region since 2014" below). With 6 years of data now available from Swarm, we are able to track a change in sign of the second time derivative or secular acceleration in this region, and to study its spatial signature at the Earth's surface and at the core surface.
Compared to its predecessor, CHAOS-6 (Finlay et al. 2016), CHAOS-7 uses a stricter criteria for selecting geomagnetically quiet times, it makes use of uncalibrated vector magnetic field data from the CryoSat-2 satellite between 2010 and 2014 by means of co-estimating magnetometer calibration parameters, and the temporal regularization is relaxed at higher spherical harmonic degrees. We recall that compared to other modern geomagnetic field models (e.g. Lesur et al. 2010;Maus et al. 2005;Sabaka et al. 2015), the CHAOS models are of intermediate complexity. They involve co-estimation of alignment parameters, internal, magnetospheric and lithospheric fields. A variety of strategies have been pursued by other recent field models, for example Alken et al. (2020) did not co-estimate magnetospheric fields while Sabaka et al. (2020) and Ropp et al. (2020) co-estimate estimate ionospheric, magnetospheric and induced fields. No attempt is made in CHAOS-7 to deterministically predict the future field evolution, rather our strategy is to provide regular updates, typically every 4 to 6 months, using the latest satellite and ground data.
In the "Data" section below, we present details of the ground-based and satellite geomagnetic measurements used to derive the CHAOS-7 model. In the "Field modelling" section, we describe the CHAOS-7 model parameterization and its estimation, including details of the co-estimated calibration model used for CryoSat-2 data. Special attention is given to how we treat fields resulting from currents induced in the electrically conducting Earth by the time-varying magnetospheric fields. In the section "Results", we report diagnostics for the CHAOS-7 field model including data misfit statistics, comparisons with ground observatory time series, spatial power spectra and maps of the internal field. Validation tests against data not used to determine the field model are presented. We also provide details on how the IGRF-13 candidate models were extracted. In the "Evolution of the South Atlantic Anomaly as seen by Swarm" section, we report on the recent changes in the region of weakest field intensity and their origin at the core surface. In the section "Field acceleration changes in the Pacific region since 2014", we focus on rapid secular variation during the Swarm era and place this in historical context. A summary and some final remarks are given in the "Conclusions" section.

Data
CHAOS-7 is based on magnetometer data collected onboard the satellites Ørsted, CHAMP, SAC-C, Cryo-Sat-2, and most importantly for the past 6 years, the three Swarm satellites, as well as an updated version of the revised monthly mean ground observatory secular variation data series (Olsen et al. 2014).

Satellite data
From the Ørsted mission, we used vector data between March 1999and December 2004and scalar data between March 1999 and June 2013 for quasi-dipole (QD) latitudes (see Richmond 1995, for a definition) poleward of ±55 • , or if attitude data were not available, each with 1-min sampling. In addition, along-track differences of Ørsted scalar data, separated by 15 s along track and with 1-min sampling, were used as a source of scalar gradient information. From the CHAMP mission, we used vector data between August 2000 and September 2010 (converted to scalar data at quasi-dipole latitudes poleward of ±55 • ) with 1-min sampling, as well as vector and scalar along-track differences (again based on measurements separated by 15 s along track), again with 1-min sampling. Vector and vector gradient data were used only when data from two star trackers were available. From the SAC-C mission, we used scalar data with 1-min sampling between January 2001 and December 2004. Alongtrack differences were not used from SAC-C. From the CryoSat-2 mission, we used uncalibrated vector data from August 2010 to December 2014, from magnetometer FGM1, but with corrections applied for temperature variations, magnetotorquer currents, and other spacecraft effects . Data were averaged to 1-min values using a robust linear fit in the magnetometer frame.
From the Swarm mission, we used the MAGX _ LR _ 1B 1 Hz calibrated data product, baseline 0505/0506, with an initial 1 minute sampling from the three satellites, Alpha, Bravo and Charlie from November 2013 to the end of August 2019. In addition, along-track gradient information was obtained from each satellite using differences between data 15 s apart, and east-west gradients were estimated from the lower pair Alpha and Charlie, using the 1-Hz data on Charlie with geocentric latitude closest to that of Alpha, sampled each minute, with the condition that the time difference was less than 50 s. Field differences are less affected by correlated noise than the data themselves and their use has been shown to improve the quality of both core and lithospheric field models . Finally, Swarm data were further down-sampled by a factor of three (i.e. effectively to a 3-min sampling rate) to account for the fact that there were three satellites contributing data during this time interval.
The following data selection criteria were applied to all data sets in an effort to focus on the internal field of interest for IGRF: • Kp ≤ 2 0 ( 3 0 for gradients) and RC-index (Olsen et al. 2014), changing at most by 2 nT/h (3 nT/h for gradients); • Merging electric field at the magnetopause averaged over the previous 2 h, E m ≤ 0.8 mV/m; • IMF B z at the magnetopause averaged over the previous 2 h is positive; • IMF B y at the magnetopause averaged over the previous 2 h is less than +3 nT when QD latitude is positive (northern QD hemisphere), i.e. −∞ < B y < 3 nT , while when the QD latitude is negative (southern QD hemisphere) it is greater than − 3nT, i.e. −3 nT < B y < ∞ (Friis-Christensen et al. 2017); • Only data from dark regions (sun at least 10 • below horizon), except for CryoSat-2 where calibration parameters and Euler angles are co-estimated using vector data from both dark and sunlit regions; • Vector field data and vector field gradients used only equatorward of ±55 • QD latitude. The distribution for B z peaks towards zero while the distributions for B y in the northern and southern hemispheres peak towards − 3 and + 3 nT, respectively. Figure 2 shows the locations of the selected Swarm scalar and vector data between September 2018 and September 2019. This illustrates the excellent geographical coverage available within a year during the Swarm era, despite the rather strict selection criteria applied. Figure 3 shows the number of data contributing to CHAOS-7 as a function of time, separated according to the source. Scalar data are only used at high QD latitudes for Swarm (see Fig. 2), CryoSat-2 and CHAMP. The variations with time in the number of data reflect the availability of data, especially for the first decade, and variations with solar cycle of the amount of data satisfying the selection criteria listed above.

Ground observatory data
Annual differences of revised monthly means of ground observatory data (Olsen et al. 2014) for the time interval January 1997 to July 2019 were utilized as a further source of information on the core field secular variation. Revised monthly means were derived from hourly mean values of 182 observatories, see Fig. 4 (including 11 with site changes during the considered time interval), which were checked for trends, spikes and other errors (Macmillan and Olsen 2013). Monthly means were calculated using a robust method relying on Huber weights (Huber 2004), from all local times after the removal of hourly estimates of the ionospheric (plus induced) field predicted using the CM4 model (Sabaka et al. 2004) and hourly estimates of the large-scale magnetospheric (plus induced) field, predicted by the CHAOS-6x9 model.

Model parameterization
The basic parametrization of the CHAOS-7 model is the same as that of earlier versions in the CHAOS series, with some minor extensions which we describe below. We assume measurements are collected in a region free from electric currents so, under the quasi-static approximation of electromagnetism, the vector magnetic field B can be Histograms showing number of selected satellite data (combination of all the scalar, scalar gradient, vector and vector gradient data) distributed according to IMF B z (top left), IMF B y (top right) and the merging electric field at the magnetopause E m , as estimated by coupling function 0.33v 4/3 B 2/3 t sin 8/3 (|�|/2) mV/m with v the solar wind speed in km/s, B t = B 2 y + B 2 z the Interplanetary Magnetic Field magnitude in the y − z plane in GSM coordinates in nT, and � = arctan(B y /B z ) . IMF and E m values are averages of 1-min values for 2 h prior to the time of the observation. y-axis shows the number of observations per bin Finlay et al. Earth, Planets and Space (2020) 72:156 represented by a scalar potential such that B = −∇V . The magnetic scalar potential V = V int + V ext consists of internal (primarily core and lithospheric) sources, and external (assumed here to be magnetospheric) sources and their internal Earth-induced counterparts. Both the internal and external parts are expanded in spherical harmonics (SH).

Internal potential fields
For the internal field, in a geographic Earth-centered Earth-fixed (ECEF) coordinate system we adopt a spherical harmonic expansion where a = 6371.2 km is chosen as the Earth's spherical reference radius, (r, θ, φ) are spherical polar coordinates, P m n are the Schmidt semi-normalized associated Legendre functions (Winch et al. 2005), g m n , h m n are the Gauss coefficients describing internal sources, and N int is the maximum degree and order of the internal expansion. The internal coefficients {g m n (t), h m n (t)} up to n = 20 are time-dependent; this dependence is represented using a basis of order 6 B-splines (De Boor 2001) with a 6-month knot separation and fivefold knots at the endpoints t = 1997.1 and t = 2020.1 . Internal coefficients for degrees 21 and above are static, a maximum degree of N int = 70 was used for the parent model estimated here. For the distributed versions of the CHAOS-7 model, at degree 25 we merge the static field estimated here to the (1) high resolution LCS-1 lithospheric field model (Olsen et al. 2017) which is provided out to degree 185.

External potential fields
Turning to the external part of the potential, we adopt an expansion in the solar magnetic (SM) coordinate system (up to n = 2 , with specific treatment of the n = 1 terms, see below) of the near magnetospheric sources and in the geocentric solar magnetospheric (GSM) coordinate system (also up to n = 2 , but restricted to order m = 0 ) of remote magnetospheric sources, (e.g. magnetotail and magnetopause currents): where  of the solid harmonics (spherical harmonics with the well-known radial scalings) in SM and GSM coordinate frames taking account of the induced field based on the diagonal part of the Q-matrix (Olsen 1999) for an assumed 3D Earth conductivity model-see below for more details. The degree 1 SM terms have the specific time-dependence: where the terms in brackets are designed to describe the magnetic field contribution due to the magnetospheric (3) ring-current and its Earth-induced counterpart. These are prescribed using the RC index which is derived from ground observatory hourly means (Olsen et al. 2014), RC(t) = ǫ(t) + ι(t) (see the section "Treatment of induced fields" below for further details). We estimate the static regression factors q 0 1 ,q 1 1 ,ŝ 1 1 and the time-varying "RC baseline corrections" q 0 1 , q 1 1 and s 1 1 in bins of 30 days. These allow for differences between the groundbased estimate of the degree 1 order 0 external magnetic signal (the RC index) and the degree 1 field seen by low-Earth orbit satellites.
It should be remembered that the CHAOS-7 magnetospheric field model is designed to represent the field during geomagnetically quiet times that were considered during the model construction. The magnetospheric field is known to have a more complex structure, in particular related to its local-time Swarm-AC denotes east-west field differences between Swarm Alpha and Charlie dependence, for higher levels of geomagnetic disturbance, as measured for example by the Kp index (see e.g. Lühr and Zhou 2020).

Treatment of induced fields
In this section, we describe in more detail how magnetic fields induced in the electrically conducting oceans and interior of the Earth, due to time-changing magnetospheric fields, are represented in CHAOS-7.
In CHAOS-7 induced fields are not parameterized separately in the form of additional internal sources since these are known to be difficult to separate from core field variations. Instead we consider induced fields calculated based on an assumed Earth-conductivity model, via the Q-responses which couple in the time domain internal (induced) and external Gauss coefficients (e.g. Schmucker 1985;Price 1967). Generally, for a 3-D Earth conductivity distribution σ (r, θ , φ) , each external coefficient induces infinitely many internal coefficients. In the frequency domain, the relation between external and induced internal coefficients for a given angular frequency ω reads (Olsen 1999): (4) ι l k (ω; σ ) = Here, Q lm kn are transfer functions, that can be arranged into the so-called Q-matrix, given by where Y l k (θ , φ) = P |l| k (cos θ) exp (ilφ) is a SH of degree k and order l, * denotes complex conjugation, dS is an elementary spherical surface area, r a = (a, θ, φ) is the position vector at the Earth's surface. B m n,r is radial magnetic field which is (numerically) computed for a given Earth conductivity model driven by a unit amplitude ( ε = 1) SH source, and is the corresponding external part of the radial magnetic field. B m n,r has in this study been numerically computed using a finite element code (Arndt et al. 2020;Grayver and Kolev 2015;Grayver et al. 2019). Formulae for the (5) unit amplitude SH source can be found, for example, in Guzavina et al. (2018). Substituting Eq. (6) into (5) yields The Q-matrix used in CHAOS-7 to couple induced and external fields is derived from an electrical conductivity model consisting of a mantle with 1-D (radial) conductivity distribution overlaid by a surface layer of laterally variable conductance. The latter approximates heterogeneous oceans and continents. In order to retain the computational advantages of a 1-D approach, we here take only the diagonal part of the Q-matrix as an approximation, since the off-diagonal elements of the Q-matrix are generally much smaller (Püthe and Kuvshinov 2014). We use the resulting diagonal Q-matrix in both the frequency domain and the time domain to compute the induced counterparts to different parts of our magnetospheric field model, these are discussed separately below. The conversion of Eq. (4) into the time domain yields a convolution integral (Maus and Weidelt 2004;Olsen et al. 2005a;Grayver et al. 2020) which couples induced and external coefficients via impulse responses of the corresponding transfer functions in the frequency domain (given by the Q-matrix). For example, for a diagonal Q-matrix and degree 1 zonal external field, the external and internal coefficients in time domain are related through a convolution integral as: where −∞ limit is replaced by a finite value of 1 year in practice.
The 1-D conductivity model used in CHAOS-7 and the related Q-kernels are presented in Fig. 5, along with a similar profile designed to illustrate a hypothetical global conductivity anomaly in the lowermost mantle. This anomalous case is unlikely for the entire lower mantle (Karato and Wang 2013), but even in this case the real part of the Q 00 11 is still relatively small at the frequencies overlapping with core field secular variation, less than 0.2 at periods of 1 year and longer.
The induced fields resulting from time variations in the degree 1 external field in SM coordinates, are calculated in the time domain from the RC index (Olsen et al. 2014) as follows. The input RC(t), which contains both external and induced parts RC(t) = ǫ(t) + ι(t) , is first detrended, then convolved with a time-domain IIR filter (Maus and Weidelt 2004;Olsen et al. 2005a) derived from the diagonal elements of the Q-matrix, rotated into geomagnetic coordinates as appropriate for working with RC(t), and based on a window of length 1 year.
In performing the convolution, we truncate the IIR filter to a length 1 year, so the response to signals older than a year are neglected, although all frequencies are present within the window are considered. ǫ(t) is then obtained by subtracting ι(t) from the original RC(t). This is done prior to deriving the CHAOS-7 model in which ǫ(t) and ι(t) parameterize the time-dependence of the degree 1 SM external and induced counterparts, only the parameters describing static regression (Eq. 3) and offset parameters (Eqs. 2a-2d) are solved for during the model estimation.
Induced fields are also accounted for in CHAOS-7 due to the static fields in SM and GSM coordinates which are time-dependent in the Earth-centered Earth-fixed (ECEF) frame, due to the "wobble" between the frames (depending on the solar position, season, etc.). A Fourier decomposition of the time-dependence resulting from each SM and GSM coefficient in the ECEF frame is carried out, then the diagonal Q-matrix described above, based on the same conductivity model of Grayver et al. (2017), is used in the frequency domain to determine the amplitude of the induced internal response to a unit excitation. Collecting these responses for all frequencies (here we used uniformly sampled frequencies with periods between 1 h and 4 years, there was little power at periods beyond 4 years) and inverse Fourier transforming provided a corresponding time-dependent induced field that scales linearly with the amplitude of the static coefficients.
To summarize, in CHAOS-7 we have made an effort to account in a consistent fashion for the induced response of all the parameterized magnetospheric sources, with the exception of the time-dependence of the offset parameters for SM degree 1 (Eq. 2a). In CHAOS-7, these offsets show variations of up to 2 nT on timescales close to 1 year and up to 4 nT on a timescales of around 10 years. Given the real part of the Q-response at these periods are likely to be less than 0.2 and 0.1, respectively (see Fig. 5), the corresponding induced responses are expected to have amplitude less than 0.5 nT.
Another important driver of induced field variations is ionospheric sources, for example the Sq current system. In CHAOS-7, we use satellite data only from the nightside and do not explicitly model this source or its induced response. It is therefore possible that we map the nightside-induced response into our internal field model (Olsen et al. 2005b). This signal is expected to leak predominantly into the zonal terms, particularly the coefficients g 0 1 and g 0 3 . One way we can assess the magnitude of this effect is through comparisons with other models which do seek to model this process. Of particular relevance here is the CM6 model (Sabaka et al. 2020). Comparisons of g 0 1 (t) and g 0 3 (t) from CHAOS-7 and CM6 show rms differences between 2000 and 2019 of 2.25 nT and 1.27 nT, respectively, with largest differences in the years 2002 and 2014. The SV in the two models generally agrees fairly well for the low-degree secular variation with rms differences in dg 0 1 /dt and dg 0 3 /dt of 0.71 nT/ year and 0.20 nT/year, respectively. Differences in the high degree SV between CHAOS-7 and CM6 are primarily due to differences in their core field temporal regularization schemes. Other models, such as that recently described by Ropp et al. (2020), which seek to directly coestimate the induced field arrive at rather different results regarding, for example dg 0 1 /dt (see the section "Timedependence of SV coefficients"). Further work is needed on this topic. ] and a hypothetical case with increased conductivity in the lower mantle (blue line, labelled "anomalous"). The Earth's core, at depths below 2850 km, is assumed to have a finite conductivity of 10 6 S/m, which is off the edge of the plot. The magnitude of the corresponding discrete impulse response (right) for the Q 00 11 transfer function. Bottom row: Q 00 11 transfer function in frequency domain, real part (left) and imaginary part (right) Finlay et al. Earth, Planets and Space (2020) 72:156

Magnetometer alignment and in-flight calibration parameters
In addition to the above spherical harmonic representation of internal and external potential fields, we co-estimate Euler angles describing the rotation between the vector magnetometer frame and the star tracker frame for Ørsted, CHAMP, CryoSat-2 and the three Swarm satellites. For Ørsted, for historical reasons, we employed two sets of constant Euler angles, implementing a break point at 00.00 on 25th January 2000. This takes account of an update of the onboard software of the star tracker that took place during 25th January 2000 at 04.05:26 (see also Olsen 2002, but note the date was incorrectly given there as 22nd January 2000). For CHAMP, CryoSat-2 and each Swarm satellite we solve for Euler angles in bins of 10 days.
In order to use the uncalibrated CryoSat-2 data, we coestimate 9 standard calibration parameters: 3 scale factors, 3 non-orthogonality angles and 3 offsets (see, for example, Olsen et al. 2003), in a series of bins of length 21 days. In each bin, these parameters relate the measured vector field in engineering units E to the calibrated magnetic field B in units of nT as follows: where the matrix describing the non-orthogonalities is that describing the scale factors for the magnetometer sensors in the three directions is while the vector containing the offsets is Further details can be found in Olsen et al. (2020) and a more detailed analysis of the coestimation procedure will appear in an upcoming study.

Summary of parameters defining the model setup
Details of the chosen SH truncation levels and of the temporal parameterization of the various parts of the model are summarized in Table 1. In all the model consists of 31,757 parameters that are simultaneously estimated from 4,007,404 magnetic field observations (counting each and gradient vector component separately).

Model estimation
The model parameters m = p T , q T , e T T , where p represents the spherical harmonic coefficients comprising the field model, q are the Euler angles and e is a vector of the calibration parameters, are determined by iteratively minimizing the following cost function using a Newtontype algorithm: where g(p) are model predictions based on field model coefficients, d(q, e) are the data, rotated to the geocentric frame using the model Euler angles q and calibrated (relevant only for CryoSat-2) using the model calibration parameters e . C d is the data error covariance matrix constructed as in previous versions of the CHAOS model series based on a priori data error estimates for each satellite, with the vector error estimates specified in the frame of the star tracker which allows the allocation of anisotropic pointing errors. Details of the regularization matrix are given below. Additional data weights proportional to sin θ were implemented for the satellite data in order to approximate an equal area data distribution. Huber data weights (Huber 2004;Constable 1988) were calculated after each iteration and used to re-weight the data; this enables robust estimation in the presence of long-tailed error distributions. Data error estimates for the ground observatory SV data were derived from the residuals to a previous model CHAOS-6x9, after detrending and taking account of data error covariances between the components of the vector triples.
In order to calibrate the CryoSat-2 magnetometer data, we include data from both sunlit and dark conditions, but only the dark data contribute to the determination of the spherical harmonic coefficients of the field model. A vector calibration is carried out at mid and low latitudes and a scalar calibration at polar latitudes (Olsen et al. 2003. More details about our approach to co-estimate calibration parameters will appear in a forthcoming study. Since scalar data are included, and because Euler angles and calibration parameters are co-estimated, the relation between the model parameters and the data is nonlinear. The cost function above is therefore iteratively minimized using a Newton-type descent method, with the Huber data weights updated at each step. The starting model was chosen to be a static internal field from a previous field model, CHAOS-6x8, evaluated in May 2015. The external field was initialized to zero. The Euler angles were initialized to the values determined in pre-flight tests, implemented via a pre-rotation step. The calibration offset and nonorthogonality parameters for CryoSat-2 were initialized to zero, while the scale factors were initialized to values of 1. Nine iterations from this starting model were carried out by which stage we judged that the model had converged to a satisfactory level; the maximum percentage change in a model parameter during the final iteration was 0.1425 %. There was no noticeable change in predictions of the internal field at Earth's surface (i.e. the IGRF relevant part of the model) during the final three iterations.
is a block diagonal temporal regularization matrix which was derived by adding contributing sub-matrices, each of which implements a quadratic measure of the temporal complexity of a certain aspect of the model. These are i3 , which implements a quadratic measure of the 3rd time derivative of the internal radial field integrated over the core surface and throughout the model timespan, i2e which implements a quadratic measure of the second time derivative of the internal radial field integrated over the core surface but only at the model endpoints t start = 1997.1 and t end = 2020.1 , sm which implements a quadratic measure of the time derivative (approximated by bin-to-bin differences) of the offset terms in the SM expansion of the magnetospheric field at Earth's surface integrated throughout the timespan, and cs , cu , cb implement quadratic measures of the time derivative (again using bin-to-bin differences) of the Cry-oSat-2 calibration scale factors, non-orthogonality angles and offsets, respectively. Each of these temporal regularization sub-matrices are scaled by regularization parameters, denoted by i3 , i2e , sm , cs , cu , cb .
There is a special treatment for i3 , which we allow to vary with the spherical harmonic degree and order (n, m). As was already the case in CHAOS-5 and 6, the zonal ( m = 0 ) terms are regularized more strongly than the non-zonal terms; in CHAOS-7 i3 (n, 0) = 10 i3 (n, m > 0) . A summary of the regularization parameters used in CHAOS-7 is given in Table 2. Preliminary test models showed that the large value of i3 required to ensure stability at low degree (where there is leakage from external fields) results in strongly suppressed time-dependence of the higher degree coefficients. The inability to retrieve information concerning the time-dependence of small length scales is particularly disappointing during the past 6 years when high-quality Swarm data are available. In order to relax the temporal regularization at higher degrees, a degree dependence of i3 was implemented. It takes its largest value i3 (n low , m) at low degree, n low < n tp min , then gradually reduces, eventually by a factor 5 · 10 −3 by degree n tp max . We set n tp min = 3 and n t max = 11 and implemented the reduction with degree using a Tukey cosine taper:    1, n < n tp min τ (n), n tp min ≤ n ≤ n tp max 0.005, n > n tp max , Table 2 Choice of regularization parameters in CHAOS-7

Fit to satellite data
We begin reporting results by presenting the fit of CHAOS-7 to its dominant contributing data source, the satellite data. Tables 3, 4, 5, 6 and 7 collect the Huberweighted means and rms residuals for different categories of data for each satellite. Overall, CHAOS-7 provides a satisfactory fit to the contributing satellite data, with no evidence for large biases and with residual histograms compatible with the assumed long-tailed error distributions. Compared to CHAOS-6, the Huber-weighted residuals between CHAOS-7 and its contributing data are slightly lower, by between 0.2 and 0.69 nT considering the nonpolar scalar data and the vector components, with the largest improvement in east-west (EW) components of the vector data, for example the weighted rms residual for Swarm A is 1.91 nT in CHAOS-7 compared to 2.49 nT in CHAOS-6. There was also a small decrease in the residuals to the along-track (AT) and EW differences, on the order of 0.05 nT. We attribute the slightly lower rms residuals in CHAOS-7 compared to CHAOS-6 to our stricter data selection, and to the relaxation of the temporal regularization in CHAOS-7.
Of particular interest is the fit to the new dataset provided by CryoSat-2. CHAOS-7 is able to fit this data, along with simultaneous ground observatory SV data, to a Huber-weighted rms level of 4.21 nT for non-polar scalar data, with vector components fit to between 4.0 and 5.25 nT. Mean residuals are less than 0.25 nT indicating little evidence for remaining biases. Figure 6 presents histograms of residuals between CHAOS-7 predictions and vector field data from CHAMP, CryoSat-2 and Swarm (top panel) and AT and EW differences of vector field (bottom panel), from CHAMP and Swarm. Considering the vector data, histograms for the Swarm data are most peaked, followed by CHAMP and then CryoSat-2. In each case the radial    components of the data are best fit. The southward field components have noticeable biases, at around 0.1-0.3 nT these biases are much smaller than the rms levels, and are likely due to imperfectly modelled ring current variations. Considering the AT difference residuals, the histograms for the Swarm data are again most peaked, with the rms residuals values of only 0.23 to 0.31 nT. The histograms for Swarm EW difference residuals show a slightly larger dispersion, as expected when taking differences between two different instruments. Table 8 presents the Huber-weighted mean and rms residuals to annual differences of the revised monthly mean ground observatory data between 1997 and 2020, that were used in constructing CHAOS-7, considering all latitudes. Also shown for reference are similar statistics for the CM6 model (Sabaka et al. 2020), model A of Alken et al. (2020) and model MCO_SHA_2Y, version 0101, an early version of the core field part of the MCM model described by Ropp et al. (2020). These three models cover similar time spans to CHAOS-7, but adopt different data selection and modelling strategies. Note that only CHAOS-7 was directly constrained by this dataset, the other models used fits to hourly mean or daily mean observatory data rather than to annual differences  It is also of interest to revisit the CHAOS-6 model in this context, to ascertain the extent to which the changes in the satellite data selection criteria and temporal regularization in CHAOS-7 have affected the fit to the ground observatory SV. The Huber-weighted rms misfits of CHAOS-6x9 to the ground observatory SV data were 3.78 nT/year, 3.62 nT/year and 3.33 nT/year for the radial, southward and eastward field components, respectively. The close agreement with the CHAOS-7 misfit levels, to within 0.05 nT/year, indicates the good agreement between the time-dependent large-scale internal field in CHAOS-6 and CHAOS-7, despite the differences in their construction. Figure 7 presents the fit of CHAOS-7 and the three other field models introduced above to time series of annual differences of revised monthly means at some example observatories. Fits to the three geocentric vector components are shown in the three columns; the latitude of the selected observatories moves from high northern latitudes down through the equator to high southern latitudes going down the rows. The model predictions are generally in agreement regarding the long-term trends but there are differences on timescales of 1 to 2 years. The model of Alken et al. (2020) generally shows more high frequency fluctuations. The CM6 model (Sabaka et al. 2020) and CHAOS-7 are generally in good agreement, despite the different treatment of induced fields on the nightside in CM6. Differences between the models are most prominent at high latitudes where there is also much larger scatter in the observatory data. Overall there is encouraging agreement between the models, and it is difficult to prefer one model over another based on these comparisons.

Co-estimated magnetometer calibration parameters for CryoSat-2
Figure 8 documents the co-estimated CryoSat-2 calibration parameters as a function of time. The nonorthogonalities are rather stable throughout the four years used. The offsets slowly vary by up to 4 nT. Largest variations are seen in the scale value of sensor one, S 1 , where especially in earlier years there are variations over about 9 months. Such time variations are not found when this instrument is calibrated against a fixed field model, CHAOS-6-x9 , suggesting that an increase in the regularization of the sensitivities may be required. This has been implemented in the most recent version of the CHAOS model, CHAOS-7.2 which was released in April 2020; the calibration parameters for CHAOS-7.2 are shown as dashed lines in Fig. 8.

Validation tests against independent ground and satellite data
In an attempt to test the performance of CHAOS-7, we have performed comparisons with independent data collected up to February 2020, more than 5 months after the construction of the model. These consist of (i) new secular variation data from ground observatories: annual differences of revised monthly means, derived from newly reported hourly mean values, using version 0122 of the ground observatory hourly mean database AUX_OBS prepared by BGS (Macmillan and Olsen 2013)-an earlier version 0121 was used in the construction of CHAOS-7, and (ii) new data collected from the Swarm satellites between September 2019 and February 2020. We consider only satellite data that fulfil the same selection criteria used in the construction of CHAOS-7 and for simplicity we focus here on scalar data. In making predictions for the new satellite data we use both the

Table 8 Model statistics of the misfit between CHAOS-7 and other field models and 30448 vector triples of ground observatory SV (annual differences of revised monthly mean) data, between 1997.5 and 2019.5
Mean and rms refer to Huber-weighted values in units nT/year. Similar misfits to other field models, the CM6 model (Sabaka et al. 2020), model A20 (Alken et al. 2020) and the MCO_SHA_2Y model (Ropp et al. 2020) are also reported for reference; note these did not directly use revised monthly means in their derivation. For CM6 28,986 SV vector triples between 1999.0 and 2019.0 were considered for the comparison, while for A20 and MCO_SHA_2Y  internal part of CHAOS-7, with a linear extrapolation after the formal end of the model in 2020, and the external part of CHAOS-7 with input from the an updated version of the RC index (see http://www.space cente r.dk/ files /magne tic-model s/RC/, and note that RC tapers to Dst for real-time values so is less accurate in this test for times after September 2019). Histograms of the residuals for the newly reported ground observatory secular variation data and the newly collected Swarm scalar data are presented in Fig. 9 and misfit statistics (unweighted mean and rms residuals) are given in Table 9, similar unweighted statistics for data actually used in the construction of CHAOS-7 are presented for reference. We find that the independent ground secular variation data are fit to within 4 nT/year, as good as the fit to the data actually used in the model construction. However,  Table 9 Statistics of the misfit between CHAOS-7 predictions and independent validation data, that was not used in the CHAOS-7 model construction Mean and rms residuals are calculated without Huber-weighting, in contrast to earlier tables. Similar statistics for the data used in the construction of the CHAOS-7 model are also shown for reference it should be acknowledged that the selection of new data is not random, only 313 of the 909 newly reported vector triples of SV data were recorded after the date of the last ground SV data used in CHAOS-7 (the remainder fill gaps in data coverage at earlier times), and there is certainly a bias towards observatories that regularly report new data. In addition, there is a different distribution of station latitudes in the validation dataset compared with the original dataset. Nevertheless, the histograms of residuals for the SV validation data in Fig. 9 show no evidence for systematic bias and there are few residuals greater than 10 nT/year. There is a trend towards slightly larger mean residuals for the new SV data collected after 2019.0 (mean residuals post-2019.0 are 0.74, − 1.18, − 0.70 nT/year for the radial, southward and eastward components, respectively); this is not a consequence of the impact of spurious offsets at a small number of observatories.
The rms misfit to independent Swarm scalar data, considering all latitudes is 5.76 nT, which is 1.58 nT higher than for the data used in the model construction. There is also a noticeable negative skew to the residual distribution. This is to be expected, given the limitations of using the CHAOS-7 external field model in predictive mode, due to RC merging with D st at the end of the RC series, and because time-dependent SM offset parameters (the "RC baseline correction") are not available.
Despite their limitations, we are encouraged by these tests that CHAOS-7 does a satisfactory job in predicting the observed magnetic field values on ground and at satellite altitude up to 5 months after the model construction, a length of time comparable to that between consecutive CHAOS model updates.

Extraction of IGRF candidate models
DTU's candidate models for IGRF-13 were obtained from the parent model CHAOS-7 parent as follows: We did not provide uncertainty estimates along with our candidate models since we are not able to calculate such estimates in a rigorous fashion; the formal error will be unrealistically small as data error correlations are ignored and the model is incomplete. The largest errors will be due to biases caused by sources not modelled (e.g. due to the polar electrojet).

Time-dependence of SV coefficients
Having established the ability of CHAOS-7 to represent adequately magnetic measurements made onboard satellites and at ground observatories over the past 21 years, we now proceed to present it's predictions concerning the structure of the global geomagnetic field and its evolution during this period. We begin by presenting in Fig. 10 time series of the first time derivatives of the Gauss coefficients dg m n /dt and dh m n /dt between 1999 and 2020. For reference purposes we also present similar series for model A of Alken et al. (2020), model MCO_SHA_2Y version 0101-an early version of the model described by Ropp et al. (2020), and the CM6 model of Sabaka et al. (2020). The top two rows show a selection of zonal harmonics, the bottom two rows examples of sectoral harmonics and the middle two rows selected tesseral harmonics, in each case increasing in degree within the two rows.
Due to the difficulties in accurately modelling the rapidly varying magnetospheric field, and in separating externally driven signals induced in the electrically conducting Earth from core field variation on timescales of months to years, the estimated SV of the zonal terms show interesting differences between the four models. Considering the axial dipole, they show differences of 1-2 nT/year. Compared to CHAOS-7, the MCO_SHA_2Y model (Ropp et al. 2020), which seeks to co-estimate induced signals, shows less time variation in its presented core field part between 2000 and 2012; during this time its axial dipole SV is almost constant. On the other hand, CM6 (Sabaka et al. 2020), which seeks to estimate induction using an a priori conductivity model, shows prominent oscillations on periods around 1 year between 2003 and 2008. Model A of Alken et al. (2020) shows larger oscillations than CHAOS-7, with periods of 1-2 years between 2010 and 2019. Here, it should be remembered that CHAOS-7 was constructed using enhanced temporal regularization of the zonal terms; this is not the case for the other models. The MCO_SHA_2Y model displays an oscillation in dg 0 7 /dt with a period 3 to 5 years which is not present in the other models. After 2005, there is good agreement for dg 0 10 /dt between the four models.
For the tesseral components there is typically rather good agreement across the four models at least up to degree 10, with the MCO_SHA_2Y model and model A of Alken et al. (2020) showing relatively more variability and CM6 being relatively smoother as the spherical harmonic degree increases. Turning to the sectoral harmonics, there is impressive agreement between the models for large-scale harmonics such as dg 3 3 /dt . The MCO_SHA_2Y model and model A of Alken et al. (2020) show sharper variations in dh 7 7 /dt , all models show similar trends in dh 15 15 /dt over the past 20 years. There are larger differences in the time-dependence of the SV across the models on going to high degree; the time changes seen in CHAOS-7 are generally slightly stronger than those seen in CM6 and weaker than those seen in the MCO_SHA_2Y model and in model A of Alken et al. (2020). These differences reflect differences in the selected data, the level to which the modellers seek to fit the data and the prior assumptions or regularization applied in order to control the field's time-dependence. Having a variety of different modelling approaches to compare thanks to the CHAOS-7, Earth's surface impetus of IGRF-13 is very informative, providing insight into which aspects are robust across the various approaches and which aspects are more challenging.

Spherical harmonic spectra of the field and its time derivatives
Spherical harmonic spectra for the internal field, its first time derivative (the SV) and its second time derivative (the SA) are displayed, colour coded by epoch, at the Earth's surface and at the core-mantle boundary (CMB) in Figs. 11 and 12, respectively. CHAOS-7 and the final CHAOS-6 model, CHAOS-6-x9, constructed used data up to May 2019, are shown for comparison. The spectra for the main field are very similar, but CHAOS-7 shows more time-dependence in its high degree SV, particularly at degree 14 and above. The SA spectra decrease more gradually at the Earth's surface in CHAOS-7, and at the core-mantle boundary they continue to trend upwards in contrast to the sharp decrease in the SA spectrum above degree 14 seen in CHAOS-6-x9. The reason for the different behaviour of the high degree SA in CHAOS-7 is the decrease in the strength of its high degree temporal regularization compared to CHAOS-6x9. The relaxation of the temporal regularization at high degree was made in order to enable the study of high degree SA during the Swarm era; there was a concern the strong regularization imposed at high degree in CHAOS-6 was preventing such changes from being recovered. This change had the desired result: the high degree SA since 2014, when Swarm data begin to provide constraints, and which was the focus for CHAOS-7 in order to provide candidate models for IGRF-13, is stable and reasonable. However, prior to 2005 the SA at degrees above 10 shows evidence for instability. In a more recent updates of the CHAOS model, attempts have been made to remedy this by increasing the upper limit to the tapering of the temporal regularization from degree 11 up to 15. Figure 13 presents maps of the radial component of the main field (MF), its first time derivative (secular variation, SV) and its second time derivative (secular acceleration, SA) from CHAOS-7 at the Earth's surface, up to SH degree 20, in 2019.0. Figure 14 presents similar maps, but downward-continued to the CMB and truncated at degrees 13, 17 and 15, respectively, for the MF, SV and SA.

Maps of the field, secular variation and secular acceleration
The surface radial magnetic field component has the well known features of strong high latitude patches and a weaker radial magnetic field in the South Atlantic. The radial field SV shows the region of largest radial field increase lies in the north-east corner of South America. On the other hand, there is a band of decreasing radial field extending south west from Southern Africa. The radial SA at the surface shows an intense localized dipolar structure in the Pacific (seen here at opposite sides of the map), with negative acceleration in 2019 in the central Pacific, including near Hawaii, and a positive acceleration in 2019 in the western Pacific north-east of Australia.
Descending to the CMB involves making the assumption that induction in the mantle plays a minor role on the timescales of several years and longer that are captured by the time-dependent internal field model. The truncation degrees (of, respectively, degree 13 (MF), 17 (SV) and 15 (SA) were chosen to ensure the maps at the CMB in 2019 were stable and well behaved. The strongest CHAOS-7, Core-mantle boundary CMB SV features are found in the equatorial region between Africa and South America; they result from the westward motion of intense radial field features this region (Olsen et al. 2014;Finlay et al. 2016). Enhanced CMB SA is seen at eastern longitudes in a band 60° to 90° east. There are also large-scale accelerations in the Pacific region, corresponding to the surface features. We return to this topic in the section "Field acceleration changes in the Pacific region since 2014".

Evolution of the South Atlantic Anomaly as seen by Swarm
In this section, we focus on changes in the South Atlantic Anomaly since the launch of the Swarm satellites in 2014 and discuss their origin in changes of the CMB radial field. Figure 15 (top panel) presents contours of the field intensity at 450 km altitude in August 2017 The mean value has been removed from each field intensity series and the maximum field change is 492 nT Finlay et al. Earth, Planets and Space (2020) 72:156 (approximately the middle epoch of the available Swarm data), as given by the internal part of the CHAOS-7 field model, with so-called single-event upsets (SEUs), recorded onboard the Swarm Alpha, Bravo and Charlie satellites between November 2013 and August 2019, superimposed. These SEUs are routinely recorded as part of bit checking procedures in the onboard electronics and they indicate when the satellite instrumentation has been affected by collisions with high-energy charged particles. The occurrence of SEUs generally increases with latitude towards the polar regions where high-energy charge particles are guided along magnetic field lines coupling the magnetosphere and ionosphere. Nevertheless, the highest concentration of SEUs is clearly observed at mid and low latitudes in the South Atlantic Anomaly weak field region. This provides a vivid illustration of the impact of the SAA on low-Earth orbit space infrastructure. Considering a sequence of such maps, the weakest field region at satellite altitude, shown by contours of blue colours in Fig. 15, has over the past 6 years slowly extended eastwards from South America towards South Africa, at latitudes between 30 and 45° south (see also Rother et al. 2020). SEU instances also show a tongue extending eastwards from the centre of the anomaly in this region. The bottom panel of Fig. 15 shows the field intensity from CHAOS-7, again at 450-km altitude in August 2017, but now overlaid with time series showing the change of the observed field intensity at Swarm altitude between 2014 to 2019 at a network of 300 geomagnetic virtual observatories (GVOs) Olsen and Mandea 2007;Barrois et al. 2018). The data presented in these time series are derived from observations within a radius of 700 km of the red target points by fitting a local potential every 4 months. Each time series has had its mean value removed, and the maximum recorded intensity change was 492 nT. The GVO series show that field intensity at the altitude of the Swarm satellites has generally decreased since 2014 over the Americas, with the most rapid decline occurring over North America and over the Pacific to the west of South America. On the other hand, the field intensity has increased over the Indian ocean and Asia. Of particular interest is what has happened in the South Atlantic Anomaly region. Here the Swarm virtual observatory series show that the field intensity has decreased on the western edges of the anomaly, leading to its westwards expansion. There have been more modest decreases in intensity within the anomaly, and intensity increases at its northeastern edge in the central Atlantic towards northern Africa. A striking fall in the field intensity was also seen on Anomaly's south-eastern edge, in the region towards southern Africa around 45°S on the Greenwich meridian, leading to an eastward expansion of the anomaly in this direction. This demonstrates that the development of the South Atlantic Anomaly is more complex than a simple westward motion and expansion of a single anomaly. Figure 16 presents the change in the field intensity at Earth's surface according to up to spherical harmonic degree 20, over the past 6 years when there is excellent data coverage provided by the three Swarm satellites. The top panel shows the field intensity in 2014, the middle panel the field intensity in 2020 and the bottom Bottom row shows the difference in intensity between 2020 and 2014 in nT, blue colours show a decrease in intensity, red colours an increase, contour lines here again show the field intensity in 2020 in steps of 500 nT from 22,500 to 27,000 nT as black contours Finlay et al. Earth, Planets and Space (2020) 72:156 Fig. 17 Core-mantle boundary origin of the South Atlantic Anomaly and its recent changes analysed using Green's functions for the Laplace equation under Neumann boundary conditions. Top row: black areas show regions of lowest intensity (under 24,000 nT) selected for analysis, the main anomaly (left column) and the new secondary minimum (right column). 2nd row: combined sensitivities G F (see Eq. (24)) showing the sensitivity of the field intensity, relative to the dipole, to the core-mantle boundary radial field, for each defined region. 3rd row: taking these sensitivities as weights and multiplying by the radial field at the core-mantle boundary. This shows the parts of the radial field at the core-mantle boundary that in combination are responsible for the main South Atlantic anomaly (left) and the new secondary minimum (right), respectively, in 2020. Bottom row: change in the sensitivity-weighted radial field at the core-mantle boundary between 2014 and 2020, showing the origin of field intensity changes for the selected regions at Earth's surface panel the difference or accumulated change over the 6 years. Contours in steps of 500 nT between 22,500 and 27000 nT are used to highlight detailed changes in the structure of the South Atlantic Anomaly. The region of weakest magnetic field over central South America has expanded during these years, and a distinct secondary minimum has developed around 40° South on the Greenwich meridian (Terra-Nova et al. 2019;Rother et al. 2020). This secondary minimum is seen even if the field is truncated at spherical harmonic degree 9. The change of intensity displayed in the bottom panel of Fig. 16 illustrates that field intensity decrease is not uniformly distributed across the South Atlantic Anomaly at Earth's surface, but is presently happening fastest at its western and southern edges, as well as in the region south west of South Africa where the secondary minimum has developed. The appearance of this distinct secondary minimum of field intensity at Earth's surface indicates that the South Atlantic anomaly must result from the combined action of a number of underlying non-dipolar flux features, and cannot be due to a single flux feature at the CMB. Due to the attenuation of small-scale features with altitude, the secondary minimum is not yet directly visible at satellite altitude although the changes responsible for its appearance at the Earth's surface are observed in the form of the eastward extension of the weak field region from the central Atlantic toward Africa. Figure 17 presents an investigation into the origin of these recent changes in the South Atlantic Anomaly. Our analysis makes use of the Green's function for the Laplace equation under Neumann boundary conditions, which formally links the core-mantle boundary radial field to the observed field (Gubbins and Roberts 1983). The top row shows the two regions of weakest intensity at the Earth's surface in 2020 (with intensities below 24000 nT) considered for the analysis, the second row shows for each region a map of the sensitivity (see below for details) of the average intensity in each region to the underlying radial field at the CMB, multiplying the CMB radial field in 2020 by these weights gives the maps in the third row. The bottom row shows changes in the CMB radial field with these weights applied between 2014 and 2020.
The sensitivity maps shown in the second row of Fig. 17 were derived as follows. Considering only sources in the Earth's core, the magnetic field B measured at a position r on or above Earth's surface may be written as a weighted integral of the radial magnetic field at the CMB B r (ŝ) (Gubbins and Roberts 1983): where dS = sin θ s dθ s dφ s is a surface element at the core-mantle boundary and G = G r , G θ , G φ are the (17) B(r) = � G(r,ŝ)B r (ŝ) dS, following Green's functions or sensitivities that link B r (ŝ) to the spherical polar vector components of the field B r (r), B θ (r), B φ (r) : where the derivative with respect to µ = cos θ is and h = c/r , where c is the radius of the core-mantle boundary.
Here, we are interested not in the vector field components themselves, but in the field intensity F = B 2 r + B 2 θ + B 2 φ , at the Earth's surface. This is a nonlinear function of B r (ŝ) . However, following Johnson and Constable (1997) and Terra-Nova et al. (2017), we may apply the chain rule and differentiate F (r) with respect to B r (ŝ) we obtain We can make use of this expression by considering changes of the magnetic field about some defined (known) background reference field B(r, t 0 ) ; Equation (22) then defines an appropriate linearized Green's function. Here, we take the background reference field to be the dipole part of the magnetic field in 2014.0 according to the CHAOS-7 model, B(r, t 0 ) = B dip (r, 2014) , and we consider the sensitivity of departures from this to the CMB radial field. In particular, we consider with �F (r, t) = F (r, t) − F dip (r, 2014) and G F (r,ŝ) calculated using Eq. (22) and the reference field B(r) = B dip (r, 2014).
Rather than considering sensitivities due to the intensity at a single point (e.g. the position of minimum intensity, as in Terra-Nova et al. (2017)), we instead consider combined sensitivities for extended low-intensity regions, by integrating over regions with intensity below 24,000 nT. The integral is performed numerically using a dense approximately equal area grid on the Earth's surface, by summing the sensitivities obtained for all gridpoints within the region of interest. To aid comparisons we normalized the summed sensitivities by the number of grid-points considered where i are the indices of grid-points within the chosen region of interest. G F defines the combined sensitivities shown in the second row of Fig. 17. When considering the region corresponding to the main minimum, there is averaging of spatially varying sensitivities from a larger number of locations, hence the resulting combined sensitivity is smoothed and of lower peak amplitude, compared to the more focused combined sensitivities obtained for the smaller secondary minimum region.
The third row and fourth rows in Fig. 17 show how the departure in intensity from the reference dipole field depends on the radial magnetic field at the CMB in 2014 and 2020, respectively. The integrated values from these maps give, for the region considered, the average departure of the intensity from the intensity of the 2014 dipole. Weak intensity (compared to the dipole) results from the integral being dominated by negative rather than positive contributions. The most important features for producing the weak intensity in the region of the secondary minimum are reversed flux features (i) under South Africa and (ii) below the Southern Atlantic between Africa and Antarctica. Turning to the main minimum, the weak intensity originates from a large reversed flux region that extends beneath the western side of South America, combined with reversed flux regions underneath the central Atlantic east of Brazil and under the Southern Atlantic.
The bottom panel in Fig. 17 shows the weighted changes in the CMB radial field responsible for the changes in the average field intensity in the two black regions. Changes in the CMB radial field at the southwestern corner of Africa, associated with the westward movement and development of the reversed flux patch under South Africa are seen to be the dominant influence on the field intensity where the secondary minimum has developed. The change of intensity in the main minimum is dominated by a negative change in the radial field under the region in the Pacific to the west of South America, and also a region of negative field change under the eastern edge of Brazil. Overall, the change to the main minimum appears to be a result of the reversed flux patch east of South America gathering towards the larger reverse flux limb that extends down western South America and is moving slowly westwards. Based on this analysis using the Green's functions, our interpretation is that the westward motion and deepening of the main minimum of the South Atlantic anomaly is a result of the westward motion of reversed flux patches under South America and the mid-Atlantic, and their gathering under South America due to their different speeds of westward drift. Large intensity variations at Earth's surface caused by the migration of flux patches have also been observed in numerical geodynamo simulations (Davies and Constable 2018). The growing secondary minimum observed to the south-west of Africa appears to primarily be due to the westward movement of an intense reversed flux feature below South Africa, which is converging towards reversed flux patches under the Southern Atlantic between Africa and Antarctica. of lower acceleration amplitudes were seen at Honolulu in the late 1970s and around 2004 at M'Bour. The recent changes seen in Guam, which is on the edge of the region of where the rapid acceleration changes have occurred, do not seem extraordinary. It therefore seems that the recent acceleration changes observed in the Pacific should not be viewed as surprising events, rather they are an integral part of the expected spectrum of rapidly changing SV behaviour that takes place at low latitudes.
By downward-continuing the field to the core-mantle boundary we have assumed a core origin for the changes in the Pacific. But it is worth pausing to consider whether externally driven induced currents might instead be responsible. In order not to be seen at all longitudes this would require exotic local conductivity anomalies below each region were such rapid SV events occur. As shown in Fig. 5, even if anomalous end-member conductivity profiles are considered, the real part of the Q-response for periods of 5 or 6 years is smaller than 0.2, thus to change the radial field acceleration by 10 nT/year 2 over 3 years, changes of order 50 nT/year 2 would be needed in the southward component. There is no evidence for such large changes in the southward component at HON in the past 6 years, or in SJG in around 1970. Further study is needed of the relationship between the rapid SV changes occasionally observed at low latitude stations in the radial field and the classical jerk signals most often reported at mid-latitudes in the Eastward component for European observatories.

Conclusions
We have presented the CHAOS-7 geomagnetic field model, the basis for DTU's IGRF-13 candidate field models, and used it to investigate geomagnetic field evolution, focusing on changes occurring over the past 6 years when excellent data coverage has been available from the Swarm trio of satellites. We find CHAOS-7 adequately represents data from the Ørsted, CHAMP, SAC-C, Cryo-Sat-2 and Swarm satellites over the past 20 years, and is able to follow the trends in secular variation observed at ground observatories. The DTU candidate model for DGRF 2015 is based directly on values of the internal field from CHAOS-7 in 2015.0, when data from the Swarm satellites were available both before and after the target epoch. The DTU candidate model for IGRF 2020 is based on the internal field from CHAOS-7 in 2019.75 (at the time of the last contributing data) propagated to 2020.0 using the secular variation from CHAOS-7 in 2019.0, when the final constraints from annual differences of ground data were available. The secular variation in 2019.0 also provides the DTU predictive SV model for the period 2020.0 to 2025.0; we adopted this simple approach since we know of no reliable way to forecast future SV changes for the upcoming 5 years.
We find that at low-Earth orbit the South Atlantic weak field anomaly continues to expand and deepen. This is directly seen in Swarm magnetic measurements and in single-event electronic upsets recorded on board the satellites. Mapping the field intensity down to the Earth's surface we find evidence from CHAOS-7 for the development a secondary intensity minimum near 40° South on the Greenwich meridian. This seems to be linked to the westward movement and evolution of a reversed flux feature at the core-mantle boundary under South Africa.
We find localized changes in radial field acceleration averaged over consecutive three year periods (2014 to 2017 and 2017 to 2020), of amplitude up to 12 nT/year 3 in the central and western Pacific. The pattern of field change resembles a localized, east-west oriented, dipole spanning 120° in longitude and confined to within 30° of the equator. Descending to the core-mantle boundary structures more confined in longitude are seen with highest amplitudes close to the equator. An important task is now to track the development of such features in detail, for example to ascertain whether energy propagates from high latitudes towards lower latitudes where it focuses at the equator.
The CHAOS-7 model, and its updates, which are typically released every 4 to 6 months, depending on Swarm data calibration and reprocessing activities, as well associated software, are available from: http://www.space cente r.dk/files /magne tic-model s/CHAOS -7/index .html.
Users interested in studying the CMB magnetic field are advised to truncate the field itself at SH degree 13, the secular variation, which is less polluted by the lithospheric field, at degree 14 (degree 17 after 2014) and the secular acceleration, which is more challenging to extract, at degree 9 (degree 15 after 2014). We have recently updated the original CHAOS-7(.1) model to CHAOS-7.2 and further updates will be released in due course. Henceforth we shall follow the convention for Swarm data products that the second version of a given model should be referred to as version 2. In CHAOS-7.2 we increased the temporal regularization of CryoSat-2 magnetometer sensitivities and increased the maximum degree of the temporal regularization tapering function in Eq. (14) from 11 to 15. The latter change was designed to mitigate instability observed in the original CHAOS-7 release in the high degree secular variation and acceleration prior to 2005; though unimportant for IGRF-13 this was undesirable. We recommend that users always use the most recent version of the model available on the