text:
Dynamic responses of the Earth's outer core to assimilation of observed geomagnetic secular variation
Weijia Kuang1 &
Andrew Tangborn2
Assimilation of surface geomagnetic observations and geodynamo models has advanced very quickly in recent years. However, compared to advanced data assimilation systems in meteorology, geomagnetic data assimilation (GDAS) is still in an early stage. Among many challenges ranging from data to models is the disparity between the short observation records and the long time scales of the core dynamics. To better utilize available observational information, we have made an effort in this study to directly assimilate the Gauss coefficients of both the core field and its secular variation (SV) obtained via global geomagnetic field modeling, aiming at understanding the dynamical responses of the core fluid to these additional observational constraints. Our studies show that the SV assimilation helps significantly to shorten the dynamo model spin-up process. The flow beneath the core-mantle boundary (CMB) responds significantly to the observed field and its SV. The strongest responses occur in the relatively small scale flow (of the degrees L≈30 in spherical harmonic expansions). This part of the flow includes the axisymmetric toroidal flow (of order m=0) and non-axisymmetric poloidal flow with m≥5. These responses can be used to better understand the core flow and, in particular, to improve accuracies of predicting geomagnetic variability in the future.
Geomagnetic field observed at the Earth's surface varies significantly in time: its temporal scales range from minutes to geological time scales. Though it was first noticed by mankind over 5000 years ago (Roberts 1992), and its origin was sought as early as 800 years ago (Dibner Library 1980), the modern theory that the geomagnetic field is generated and maintained by convective flow in the Earth's outer core (geodynamo) was originated from the seminal work of Larmor (1919). Successful numerical simulation of the geodynamo was first carried out by Glatzmaier and Roberts (1995), and then followed by Kageyama and Sato (1997) and by Kuang and Bloxham et al. (1997). Christensen et al. (2010) provided a comprehensive summary of numerical geodynamo solutions and their relevances to geomagnetic observations.
Assimilation of geomagnetic observations with numerical geodynamo models started less than a decade ago. Sun et al. (2007) and Fournier et al. (2007) used simplified magnetohydrodynamic (MHD) systems, and synthetic data tested the applicability of assimilation of sparse magnetic data. Liu et al. (2007) first used observation system simulation experiments (OSSEs) with a full dynamo model and demonstrated clearly that one could use assimilation of magnetic field at the surface to estimate the dynamo state deep in the fluid core. Kuang et al. (2008) published the first working geomagnetic data assimilation system MoSST −DAS in which the Gauss coefficients of various geomagnetic and paleomagnetic field models are assimilated with their MoSST geodynamo model (Kuang and Chao 2003; Jiang and Kuang 2008) for estimation of the core state and prediction of geomagnetic field variation. Kuang et al. (2009) then used this assimilation system and 100 years of the Gauss coefficients from gufm1 (Jackson et al. 2000) and CM4 (Sabaka et al. 2004) to understand the responses of the core state to surface geomagnetic observations, and their implications to core state estimation and secular variation (SV) prediction. We refer the reader to Fournier et al. (2010) for a comprehensive review of the data assimilation algorithms for geomagnetic data assimilation (GDAS) and some of the early results.
Rapid advances have occurred in multiple facets of GDAS. Several independent assimilation systems have been developed to understand better the core dynamical state. For example, Aubert and Fournier (2011) and Fournier et al. (2011, 2013) carried out OSSEs with synthetic observations and numerical dynamo models to examine possibilities of core state determination. Fournier et al. (2011, 2013) also tested their approach with a geomagnetic field model. Aubert (2013, 2014) investigated possibilities of inverting core state properties using the observed field and SV. In addition to the sequential data assimilation systems mentioned above, there are also efforts in developing GDAS systems based on variational data assimilation techniques. For example, Li et al. (2011, 2014) have been continuing their effort on a new combined forward and adjoint system towards a full geodynamo model. Encompassed application is the contributions of assimilation results to international geomagnetic reference field (IGRF) (Kuang et al. 2010) and efforts to determine field model error statistics (Gillet et al. 2013).
Despite these advances, GDAS is still in an early stage similar to that of early numerical weather prediction (NWP) (for a more comprehensive review, see, e.g., Kalnay 2003). Many important questions are still to be fully answered, such as comprehensive assessment of numerical dynamo model biases, observation and core state covariances and error statistics, and the dynamic responses of dynamo state to the observed geomagnetic field. The latter is of particular importance to the spin-up processes of the numerical models which is characterized by the difference between the observation and the forecast, often called (\({\mathcal {O}}\)-\({\mathcal {F}}\)). These, in turn, determine how fast and how close the numerical solutions can be pulled to the true state of the core.
Concerns on the spin-up of the numerical models can be examined from the time scales of the observed field and of the numerical models. Global field model results from the past 400 years of geomagnetic data (e.g., Jackson et al. 2000; Sabaka et al. 2004, 2015; Olsen et al. 2006, 2014) show that the typical time scales τ l of the degree l components (Stacey 1992; Hulot and Le Mouël 1994; Olsen et al. 2006)
$$ \tau_{l} = \left[\frac{\sum_{m} \left({g_{l}^{m}}\right)^{2}+\left({h_{l}^{m}}\right)^{2}} {\sum_{m} \left({\dot g}_{l}^{m}\right)^{2}+\left({\dot h}_{l}^{m}\right)^{2}}\right]^{1/2} $$
((1))
varies from over 1000 years for the dipole (l=1) to less than 100 years for higher degrees. In (1), \(\left ({g_{l}^{m}}, {h_{l}^{m}}\right)\) are the Gauss coefficients of the field, and \(\left ({\dot g}_{l}^{m}, {\dot h}_{l}^{m}\right)\) are their first order time derivatives, i.e., the Gauss coefficients of the SV. Currently, the longest record for low-degree (l≤5) field coefficients is from the paleo/archeo magnetic data (e.g., Korte et al. 2011; Nilsson et al. 2014). The high-quality coefficients for up to degree l≤8 could be obtained from historical and observatory data (Jackson et al. 2000). Very high quality coefficients for degrees l≤13 are obtained in the past 50 years with satellite magnetic data (Sabaka et al. 2004, 2015; Olsen et al. 2006, 2014). In summary, the data record is no more than 10 times of the typical time scales of the geomagnetic field. This brings the very concern on whether the observational record is sufficient to spin up numerical dynamo models. OSSEs results also suggest long spin-up time in geomagnetic data assimilation. For example, Liu et al. (2007) showed, via their OSSEs with synthetic magnetic data and a fully nonlinear dynamo model, that the difference between the forecast and the truth reach the minimum in approximately 40 % of the magnetic free decay time, or 8000 years, if only the poloidal field (of the first eight spherical harmonic degree coefficients) is assimilated (see Figure four in Liu et al. 2007). Therefore, the model spin-up also has direct consequence on estimation of the core state.
How could we improve geomagnetic data assimilation systems within the observational limit? There are several areas for improvements. For example, improvements in global geomagnetic field modeling are needed since the Gauss coefficients from various field models have been used in most of the previous GDAS studies. Currently, there are many field models covering different epochs (e.g., Jackson et al. 2000; Korte et al. 2011; Gillet et al. 2013; Olsen et al. 2014; Sabaka et al. 2015). A unified field model covering the longest possible period could certainly reconcile differences in these models and thus help greatly GDAS systems. There is an ongoing effort on constructing a unified global field model for the past millennium (private communication with Korte). The field model error statistical information of such unified field models, such as those in Gillet et al. (2013), is also necessary for GDAS.
Improvement in the assimilation algorithms could also help data utilization. Some efforts were made by Kuang et al. (2010) in which a subset of the Gauss coefficients (of lower degrees) with much longer records are assimilated first to spin up the model, followed by assimilating those of higher degrees for the past 100 years. Tangborn and Kuang (2015) showed, via a set of experiments, that such assimilation methodology can have positive impact on core state and improve accuracies of predicting the subset of the Gauss coefficients not assimilated. Another example is employment of ensemble Kalman filtering (EnKF) approach (Evensen 1994). Fournier et al. (2011, 2013) used OSSEs to show the potential to speed up the transfer of information from geomagnetic data to the core state. But such speedy transfer depends on model errors (that are in general very large due to limitations of numerical dynamo models) not considered in their studies. It should also point out that GDAS is computationally very expensive. Such expense needs to be considered in the algorithm improvement.
Another improvement is on exploiting and utilizing further geodynamic information embedded in surface geomagnetic measurements. An immediate candidate for such exploitation is the geomagnetic secular variation (SV), described by the first-order time derivative \(\left ({\dot g}_{l}^{m}, {\dot h}_{l}^{m}\right)\) of the Gauss coefficients since, as we will describe in the next section, they provide additional constraints on the core flow beneath the core-mantle boundary (CMB) and on the radial variation of the magnetic field. The former is not new, as there is a long history of, started from Roberts and Scott (1965), core flow inversion from observed SV at the Earth's surface via the "frozen-flux" approximation (in which the Ohmic dissipation beneath the CMB is ignored). However, this approximation comes with the price: the core flow cannot be uniquely inverted (e.g., Roberts and Scott 1965; Backus 1968). Thus, additional constraints on core flow properties are necessary in such core flow inversion studies (for more complete reviews, please read, e.g., Holme 2007; Kuang and Tangborn 2011). If the Ohmic dissipation is retained (no "frozen flux" approximation), then the observed SV imposes the constraints on the radial variation of the field in the core, as the latter is part of the magnetic induction. Since both field advection and Ohmic dissipation are included in geodynamo modeling, both kinds of constraints can be examined in MoSST −DAS or any other GDAS system without mathematical difficulties.
Therefore, a natural expansion of data utilization in GDAS is to assimilate both the field and its SV, so that the embedded geodynamic constraints can be used to make more optimal analysis, thus speeding up the transport of information from the surface geomagnetic observations to the dynamical state in the outer core. Since the SV is not included in the state vector of numerical geodynamo models, it will be connected through a non-linear observation operator,, which transforms the model state space to the observations space. Obviously, will depend on, among others, fundamental physical properties of the magnetic field.
It should be pointed out here that assimilating the rate of change of geodynamic observables has been routinely used in NWP. For example, precipitation rate measured from a variety of satellite instruments is assimilated, despite not being a state variable in a GCM (Hou et al. 2000).
There is also a long history of attempting to invert core state from surface observations. In addition to core flow inversion initiated by Roberts and Scott (1965), Zatman and Bloxham et al. (1997) extended further the inversion of the poloidal magnetic field in the outer core. More recently, Aubert (2013, 2014) attempted a new approach to invert core dynamical state with the observed SV and the dynamo models. These inversions, however, are not data assimilation. However, the inversion results could benefit geomagnetic data assimilation. For example, the inversion by Aubert (2013, 2014) could be utilized for making "analysis" in a geomagnetic assimilation system.
In this paper, we describe in detail the results from our recent effort on assimilation of both the field and its SV. These results, from a series of experiments, will demonstrate the improvement in prediction and knowledge on core flow responses to the SV assimilation. The results also provide valuable information for further development in this direction.
This paper is organized as follows: the numerical model details and the mathematical formulation for SV assimilation will be given in the next section. Followed are the experimental results we have with this assimilation approach. Discussions and plans for further improvements are presented in the last section.
The mathematical formulation for SV assimilation depends on the numerical geodynamo models and the assimilation algorithms, in addition to the physics controlling the time variation of the magnetic field. In this section, we provide the mathematical methodologies used in MoSST −DAS employed in this study (Kuang et al. 2008; Sun and Kuang 2015). But, with some modifications, they can be applied to other GDAS systems.
Dynamo state vector and geomagnetic observation
MoSST −DAS utilizes the MoSST core dynamics model for time integration of the magnetic field (Kuang et al. 2008; Sun and Kuang 2015). In this system, the state vector x
$$ \mathbf{x} = \left(\mathbf{v}, \mathbf{B}, \delta\varrho\right)^{T} $$
includes the velocity field v and the density anomaly δ ϱ in the outer core r i ≤r≤r c (r i and r c are the mean radii of the inner core boundary (ICB) and CMB, respectively); and the magnetic field B in the outer core, the electrically conducting inner core r≤r i , and the D ″-layer r c ≤r≤r d (r d is the mean radius at the top of the layer). The superscript "T" in (2) implies the transpose. The solid mantle above the D ″-layer r d ≤r≤r s (r s is the mean radius of the Earth's surface) is electrically insulating. The whole system is defined in the reference frame fixed with the solid mantle.
The velocity field v and the magnetic field B are decomposed into the poloidal and toroidal components, with the scalars described via spherical expansions
$${} \begin{aligned} \left(\mathbf{v}, \mathbf{B}\right)^{T} &= \nabla\times \left[\left(T_{v}, T_{b}\right)^{T}\mathbf{1}_{r}\right]+ \nabla\times\nabla\times \left[\left(P_{v}, P_{b}\right)^{T}\mathbf{1}_{r}\right], \end{aligned} $$
$${} \begin{aligned} \left(P_{v}, T_{v}, P_{b}, T_{b}, \delta\varrho\right)^{T} &= \sum_{0\le m\le l}^{L_{M}}\left({{v_{l}^{m}}},{{\omega_{l}^{m}}},{{b_{l}^{m}}}, {j_{l}^{m}}, {\vartheta_{l}^{m}}\right)^{T}\\ &\quad\times{Y_{l}^{m}}(\theta,\phi) + C.C., \end{aligned} $$
where 1 r is the unit radial vector, θ is the co-latitude, ϕ is the longitude, \({Y_{l}^{m}}\) are the fully normalized spherical harmonic functions of degree l and order m, L M is the truncation order, and C.C. implies the complex conjugate part. P and T in (3) are called the poloidal and toroidal scalars. It is therefore convenient to write
$$ \mathbf{x} \,\ = \left(\mathbf{x}_{v},\mathbf{x}_{\omega}, \mathbf{x}_{b},\mathbf{x}_{j},\mathbf{x}_{\rho}\right)^{T}, $$
where the subsets are defined with the relevant spectral coefficients in (4), e.g.,
$$ \mathbf{x}_{b} = \left\{{{b_{l}^{m}}}(r_{k}) \mid 0 \le r_{k} \le r_{d};\ 0\le m\le l\le L_{M} \right\}^{T} $$
for the poloidal magnetic field. (5) and (6) can be different for other dynamo models.
In geomagnetic field modeling, geomagnetic measurements are used to obtain the magnetic field B o originated from the core (simply called the geomagnetic field hereafter) that is described as
$$\begin{array}{@{}rcl@{}} \mathbf{B}^{o} & =& -\nabla\Psi, \end{array} $$
$$\begin{array}{@{}rcl@{}} \Psi \, \ & =& r_{s}\! \sum_{0\le m\le l}^{L_{o}} \left(\frac{r_{s}}{r}\right)^{l+1}\!\left({{g_{l}^{m}}}\cos m\phi + {{h_{l}^{m}}} \sin m\phi\right) {P_{l}^{m}}(\theta) \end{array} $$
where \({P_{l}^{m}}\) is the Schmidt normalized associate Legendre polynomial of degree l and order m, \(\left ({{g_{l}^{m}}}, {{h_{l}^{m}}}\right)\) are the Gauss coefficients (slightly different from the standard notation), and L o is the maximum degree (L o ≤13 in general). Since these Gauss coefficients \(\left ({{g_{l}^{m}}}, {{h_{l}^{m}}}\right)\) are provided by different field models over the past 10,000 years (e.g., Jackson et al. 2000; Korte et al. 2005, 2011; Gillet et al. 2013; Olsen et al. 2014; Sabaka et al. 2015), they are used as the "observations" in our study.
By (3), (4), (7), and (8), we can obtain the relationship between \(\left ({{g_{l}^{m}}}, {{h_{l}^{m}}}\right)\) in (8) and the observed \({b_{l}^{m(o)}}\) in the form of (4) via the radial component B r of the magnetic field B
$$ \begin{aligned} B_{r}^{o} & = - \frac{\partial\Psi}{\partial r} = \sum_{0\le m\le l}^{L_{o}} (l+1)\left(\frac{r_{s}}{r}\right)^{l+2}\\&\quad\times\left({{g_{l}^{m}}}\cos m\phi + {{h_{l}^{m}}}\sin m\phi\right) {P_{l}^{m}} (\theta)\\ & = -\frac{{\hat L}}{r^{2}}P_{b}^{{o}} = \sum_{0\le m\le l}^{L_{o}}\frac{l(l+1)}{r^{2}}{b_{l}^{m(o)}} {Y_{l}^{m}} + C.C. \end{aligned} $$
where \(\hat L\) is the angular momentum operator. With the definitions of \({Y_{l}^{m}}\) and \({P_{l}^{m}}\), (9) requires that
$$ \begin{aligned} {b_{l}^{m(o)}}(r) &= \frac{{r_{s}^{2}}}{l}\left(\frac{r_{s}}{r}\right)^{l}G_{m} \left({{g_{l}^{m}}}-i{{h_{l}^{m}}}\right),\\ G_{m} &= \left[\frac{2\pi(1+\delta_{m0})}{2l+1}\right]^{1/2} \end{aligned} $$
for r d ≤r≤r s . The spectral coefficients of the SV are the time derivatives of (10):
$${} {\dot b}_{l}^{m(o)}(r) = \frac{{r_{s}^{2}}}{l}\left(\frac{r_{s}}{r}\right)^{l}G_{m}\left({\dot g}_{l}^{m}-i{\dot h}_{l}^{m}\right) \quad\text{for}\quad r_{d}\le r\le r_{s}, $$
where \(\left (\dot \ \right)\) means the time derivative.
SV and core state
Geomagnetic observations only provide the time series of \(\left ({{g_{l}^{m}}}, {{h_{l}^{m}}}\right)\). The SV coefficients \(\left ({\dot g}_{l}^{m}, {\dot h}_{l}^{m}\right)\) are actually derived. Assimilation of the SV thus raises two major concerns: could the SV be approximated as "instantaneously" measured and whether it is redundant to the assimilation of the field?
Answers to the first concern depend on the significance of numerical errors in SV calculation. Consider, for example, a central difference scheme is used:
$${\dot g}_{l}^{m}(t) \approx \frac{{{g_{l}^{m}}}(t+\delta t)-{{g_{l}^{m}}}(t-\delta t)}{2\delta t}. $$
Then, the relative numerical error can be estimated as
$$\frac{{{g_{l}^{m}}}(t+\delta t)-{{g_{l}^{m}}}(t-\delta t)}{2\delta t} = {\dot g}_{l}^{m}(t) + {\mathcal{O}}\left[\frac{\partial^{3} {{g_{l}^{m}}}}{\partial t^{3}}\left(\delta t\right)^{2}\right] $$
Since \(\partial ^{3} {{g_{l}^{m}}}/\partial t^{3} \approx {\dot g}_{l}^{m}/{\tau _{l}^{2}}\) by (1) and δ t≈t o (the typical time intervals of data series), we have
$$\begin{aligned} \frac{{{g_{l}^{m}}}(t+\delta t)-{{g_{l}^{m}}}(t-\delta t)}{2\delta t} &= {\dot g}_{l}^{m}(t)\left[1+ {\mathcal{O}}\left(\tau_{o}/\tau_{l}\right)^{2}\right]\\ &\equiv {\dot g}_{l}^{m}(t)\left(1+\epsilon_{n}\right) \end{aligned} $$
In general, τ o ≤1 month in the field models using modern observatory and satellite data (e.g., Sabaka et al. 2004, 2015; Olsen et al. 2006, 2014), while τ l can be as short as 10 years (Christensen et al. 2012) for high-degree coefficients. Thus, ε n ≈10−6, which leads to an order 10−4 nT/year error in SV. On the other hand, the external field is several tens of nanotesla at the Earth's surface (Sabaka et al. 2015) on the solar cycle (∼11 years) and shorter time scales. Thus, ε n is negligible compared to those arising from, e.g., separation of the external and the internal magnetic signals. One could then argue that both the field and its SV are "concurrently" measured.
The redundancy is not an issue because the observed SV brings different knowledge of the core state x compared to the observed field. To see this, let us consider the magnetic induction of the poloidal magnetic field beneath the impenetrable and "free-slip" CMB (\(r=r_{c}^{-}\))
$${} {\dot b}_{l}^{m} = -\frac{r^{2}}{l(l+1)}\left[\nabla_{h}\cdot\left(\mathbf{v}_{h} B_{r}\right)\right]_{l}^{m} + \eta \left[\frac{\partial^{2}}{\partial r^{2}}-\frac{l(l+1)}{r^{2}}\right] {{b_{l}^{m}}}, $$
and in the D ″-layer
$$ {\dot b}_{l}^{m} = \eta_{d} \left[\frac{\partial^{2}}{\partial r^{2}}-\frac{l(l+1)}{r^{2}}\right] {{b_{l}^{m}}}. $$
In (12), the subscript "h" implies the horizontal components of the velocity field v and η is the magnetic diffusivity of the outer core fluid; η d in (13) is the magnetic diffusivity of the D ″-layer (η≤η d in general). These two equations show clearly that the observed \({\dot b}_{l}^{m(o)}\) will impose the constraint on v and on the non-potential part of the poloidal field.
The latter, i.e., (13), implies that, at the top of the D ″-layer (r=r d ), a potential poloidal field \(b_{l}^{m(p)}\) can fully recover the observed field \({b_{l}^{m(o)}}\). However, it cannot recover the observed SV \({\dot b}_{l}^{m(o)}\) since
$$\frac{\partial^{2}b_{l}^{m(p)}}{\partial r^{2}}-\frac{l(l+1)}{r^{2}} b_{l}^{m(p)} = 0. $$
In other words, the observed SV provides the information on the non-potential part of the field that is missing in the observed field \({b_{l}^{m(o)}}\). Therefore, SV assimilation is not redundant to the field assimilation.
Indeed, our assimilation results (in Fig. 2) demonstrate clearly that assimilation of \({b_{l}^{m(o)}}\) could not reduce the differences between the forecast SV \({\dot b}_{l}^{m(f)}\) and the observed SV \({\dot b}_{l}^{m(o)}\), called (\({\mathcal {O}}\)-\({\mathcal {F}}\)) of the SV, although that of the field is reduced very rapidly in the first few analysis cycles, a strong indication for the need of SV assimilation.
The rms \(({\mathcal {O}}\text {-}{\mathcal {F}})_{B}\) of the magnetic field in case II (dashed line) and case III (solid line). In both cases, \(({\mathcal {O}}\text {-}{\mathcal {F}})_{B} \ll 1\) and decays monotonically after the first three analysis cycles and then levels off in the last 20 years. This shows the continuing improvement in the forecast accuracies. In addition, the \(({\mathcal {O}}\text {-}{\mathcal {F}})\) results in case III (with the assimilation of \({b_{l}^{m(o)}}\) and \({\dot b}_{l}^{m(o)}\)) are in general more than 20 % smaller than in case II (with only the assimilation of \({b_{l}^{m(o)}}\)), showing a clear improvement in forecast accuracies
Similar to Fig. 1, but for \(({\mathcal {O}}\text {-}{\mathcal {F}})_{\dot B}\) of the SV. In case II (dashed line), \(({\mathcal {O}}\text {-}{\mathcal {F}})_{\dot B} = {\mathcal {O}}(1)\) for much of the assimilation period before decays gradually in the last 20 years, implying that there is no similarity between the forecasted SV \({\dot b}_{l}^{m(f)}\) and the observed SV \({\dot b}_{l}^{m(o)}\). But its magnitude is much smaller in case III (solid line), and it decays monotonically in time, indicating that the SV assimilation accelerates the spin-up process
New assimilation approach
We have been using the sequential assimilation approach in MoSST −DAS (e.g., Kuang et al. 2008; Sun and Kuang 2015). It can be summarized as follows: at the analysis time t a when the observation y is made, a new initial condition x a (called the "analysis") is made from the forecast x f and the observation y, future forecast for t>t a can then be made with the following initial value system:
$$ \frac{\partial \mathbf{x}^{f}}{\partial t} = \mathbf{M}\left(\mathbf{x}^{f}\right), \qquad \mathbf{x}^{f}(t_{a}) = \mathbf{x}^{a}. $$
If there is a linear observation operator H that projects x to the observation space (where y is defined), then the analysis x a is of the form
$$\begin{array}{@{}rcl@{}} && \mathbf{x}^{a} = \mathbf{x}^{f}+\mathbf{K}\left(\mathbf{y}-\mathbf{H}\mathbf{x}^{f}\right) \end{array} $$
$$\begin{array}{@{}rcl@{}} && \mathbf{K} \, = \mathbf{P}^{f}\mathbf{H}^{T}\left(\mathbf{H}\mathbf{P}^{f}\mathbf{H}^{T}+\mathbf{R}\right)^{-1} \end{array} $$
where K is called the gain matrix, P f and R are the error covariance matrices of the forecast x f and of the observation y, respectively. (15) is obtained to minimize the error |H·(x t−x a)|2 between the analysis x a and the truth x t.
In our assimilation system MoSST −DAS, the error covariance P f is calculated via three different approaches: an ensemble-based covariance analysis (Sun et al. 2007; Sun and Kuang 2015), an empirical covariance based on forecast solution properties (Tangborn and Kuang 2015), and an optimal interpolation (OI) scheme with a predefined time-invariant covariance (Kuang et al. 2009). The first approach is computationally very expensive, in which an ensemble N initial states are created x i =x f+ε i (ε i are random white noise perturbations). Their free running model solutions at a later time (often a fraction of the magnetic free decay time) are then used to calculate the covariance \(\mathbf {P}^{f} = \langle \left (\mathbf {x}_{i}-{\mathbf {\overline {x}}}\right) \left (\mathbf {x}_{i}-{\mathbf {\overline {x}}}\right)^{T}\rangle \) (\(\mathbf {\overline {x}} = \sum \mathbf {x}_{i}/N\) is the mean). The empirical covariance is assumed diagonal (no cross covariance between different degrees l and orders m) and is determined by the forecast error standard deviations and an exponentially decaying spatial correlation function. This approach can be updated at analysis time without much computational effort. The OI-type error covariance is similar to the empirical one, except that the forecast error standard deviations are assumed time-invariant (i.e., constant throughout the assimilation) and that spatial correlation does not decay exponentially with the distance. We simply use the OI-type error covariance in this study.
If only the observed field is assimilated, then
$$ \mathbf{y} = \left\{{b_{l}^{m(o)}}(r_{d})\Big\vert 0\le m\le l \le L_{o}\right\}^{T} \equiv \mathbf{y}_{b}. $$
By (2) and (5), H is linear and very simple
$$ \mathbf{H} = \left({\mathbf 0}, {\mathbf 0}, \mathbf{H}_{b}, {\mathbf 0}, {\mathbf 0}\right)^{T}, $$
where H b corresponds to the subset x b and has only non-zero entries for \({{b_{l}^{m}}}(r_{d})\) with 0≤m≤l≤L o . If the observed SV is also assimilated, then
$$\begin{array}{@{}rcl@{}} && \mathbf{y} \ = \left(\mathbf{y}_{b},\mathbf{y}_{\dot b}\right)^{T} \end{array} $$
$$\begin{array}{@{}rcl@{}} && \mathbf{y}_{\dot b} \equiv \left\{{\dot b}_{l}^{m(o)}(r_{d})\mid 0\le m\le l \le L_{o}\right\}^{T}, \end{array} $$
However, by (12) and (13), transformation between \(\mathbf {y}_{\dot b}\) and x f is a differential-functional projection and is denoted as (x f). One could of course construct an independent projection system which evaluates (x f) directly (see, e.g., Remarks 5.3.1 in Kalnay 2003). Alternatively, a linearization approximation (x f)≈H·x f could be made so that (15) can still be used.
There are different means to linearize (x f). In our current study, we create an effective observed field \({{\widetilde {b}}_{l}^{m(o)}}\) defined in the D ″-layer that matches both \({b_{l}^{m(o)}}\) and \({\dot b}_{l}^{m(o)}\). In this approach, \({{\widetilde {b}}_{l}^{m(o)}}\) comprises of a potential field that accounts for \({b_{l}^{m(o)}}\) and a non-potential field for \({\dot b}_{l}^{m(o)}\):
$$ \begin{aligned} {{\widetilde{b}}_{l}^{m(o)}}(r) &= \left(\frac{r_{d}}{r}\right)^{l}{b_{l}^{m(o)}}(r_{d})+ \frac{1}{2\eta_{d}}\\ &\quad\times\left(r-r_{d}\right)^{2}{\dot b}_{l}^{m(o)}(r_{d}) \quad \text{for}\quad r_{c}\le r \le r_{d}. \end{aligned} $$
Obviously, at the top of the D ″-layer r=r d ,
$$\begin{aligned} {{\widetilde{b}}_{l}^{m(o)}} &= {b_{l}^{m(o)}},\quad \frac{\partial{{\widetilde{b}}_{l}^{m(o)}}}{\partial r} = \frac{\partial{b_{l}^{m(o)}}}{\partial r}\\ &= - \frac{l}{r_{d}}{b_{l}^{m(o)}},\quad{{\dot {\widetilde{b}}}_{l}^{m(o)}} = {\dot b}_{l}^{m(o)}. \end{aligned} $$
The coefficients \({{b_{l}^{m}}}(r)\) in the D ″-layer can be approximated via Taylor series expansion
$$\begin{aligned} {{b_{l}^{m}}}(r) &= {{b_{l}^{m}}}(r_{d}) + {{b_{l}^{m}}}^{\prime}(r_{d}) \left(r-r_{d}\right) + \frac{1}{2}{{b_{l}^{m}}}^{\prime\prime}\\ &\quad\times(r_{d}) \left(r-r_{d}\right)^{2} + {\mathcal{O}}\left[{{b_{l}^{m}}}^{\prime\prime\prime}\left(r-r_{d}\right)^{3}\right], \end{aligned} $$
where ′ implies the radial derivative. Since the approximation (21) satisfies the first three terms in the above expansion, the errors are given by the last term. Since \({{b_{l}^{m}}}^{\prime \prime \prime }\sim {{b_{l}^{m}}}/r^{3}\) in D ″-layer, the relative errors are of the order [(r d −r c )/r c ]3. For a 20 km layer thickness, it is smaller than 10−6. (21) allows us to extend the surface observations to the CMB. The observation vector y is now of the form
$$ \mathbf{y} = \left\{{{\widetilde{b}}_{l}^{m(o)}}(r) \mid 0\le m\le l \le L_{o};\ r_{c}\le r \le r_{d}\right\}^{T}\equiv \mathbf{y}_{{\widetilde b}}. $$
The observation projection is again linear:
$$ {\boldsymbol{\cal{H}}}(\mathbf{x}^{f})=\mathbf{H}\cdot\mathbf{x}^{f} $$
with H defined in (18). However, H b now includes non-zero entries on all grid points in the D ″-layer r c ≤r≤r d .
We can use this approach to further construct an effective observed velocity field \({\mathbf {\widetilde {v}}}^{o}\) beneath the CMB \(\left (r = r_{c}^{-}\right)\). Since \({\dot b}_{l}^{m}\) is continuous across the CMB, by (12), (13) and (21), we have
$${} \begin{aligned} -\frac{{r_{c}^{2}}}{l(l+1)}&\left[\nabla_{h}\cdot\left({\mathbf{\widetilde{v}}}_{h}^{o} {\widetilde B}_{r}^{o}\right)\right]_{l}^{m} +\eta \left[\frac{\partial^{2}}{\partial r^{2}}-\frac{l(l+1)}{{r_{c}^{2}}}\right] {{\widetilde{b}}_{l}^{m(o)}} \\&= {{\dot {\widetilde{b}}}_{l}^{m(o)}}(r_{c}) = {\dot b}_{l}^{m(o)}(r_{d})\!\!\left[\!1-\frac{l(l+1)}{2{r_{c}^{2}}}\!\left(r_{c}-r_{d}\right)^{2}\!\right] \end{aligned} $$
For those geodynamo models without an electrically conducting D ″-layer, one can simply replace \({{\dot {\widetilde {b}}}_{l}^{m(o)}}(r_{c})\) in the above equation with \({\dot b}_{l}^{m(o)}(r_{c})\) continued downward from the surface observation.
Obviously, (24) is an underdetermined system, since both \({{\widetilde {b}}_{l}^{m(o)}}\) and \({\mathbf {\widetilde {v}}}^{o}\) are unknown at \(r_{c}^{-}\). But one can find the "best-fit" \({\mathbf {\widetilde {v}}}^{o}\) and \({{\widetilde {b}}_{l}^{m(o)}}\) via minimizing the following difference
$$ \begin{aligned} &\min_{{\mathbf{\widetilde{v}}}^{o}, {{\widetilde{b}}_{l}^{m(o)}}}\left|{{\dot {\widetilde{b}}}_{l}^{m(o)}}(r_{c}) + \frac{{r_{c}^{2}}}{l(l+1)}\left[\nabla_{h}\cdot\left(\mathbf{v}_{h} B_{r}\right)\right]_{l}^{m}\right.\\ &\left.\qquad\qquad\qquad\,\,\,\,- \eta \left[\frac{\partial^{2}}{\partial r^{2}}-\frac{l(l+1)}{{r_{c}^{2}}}\right] {{b_{l}^{m}}}\right|^{2}. \end{aligned} $$
We provide here only a sketch of the minimization procedure: denoting by x b and x v the vectors of the spectral coefficients of B r and of v h in (25), respectively, and then (25) is equivalent to minimize the (inner) product
$$ \left[\dot{\mathbf{x}}_{b} + {\mathbf{A}}\cdot\mathbf{x}_{b}+ {\mathbf{B}}(\mathbf{x}_{b})\cdot\mathbf{x}_{v}\right]^{T}\cdot \left[\dot{\mathbf{x}}_{b} + {\mathbf{A}}\cdot\mathbf{x}_{b}+ {\mathbf{B}}(\mathbf{x}_{b})\cdot\mathbf{x}_{v}\right] $$
with respect to x v . Note that A·x b and B(x b )·x v in (26) describe the diffusion and the advection terms in (25), respectively. Minimizing the product (26) leads to the solution
$$ \left({\mathbf{B}}^{T}\cdot{\mathbf{B}}\right) \cdot {\tilde {\mathbf{x}}}_{v} = -{\mathbf{B}}^{T}\cdot\left[\dot{\mathbf{x}}_{b} + {\mathbf{A}}\cdot\mathbf{x}_{b}\right] $$
which provides the spectral coefficients of the effective observed velocity field \({\mathbf {\widetilde {v}}}^{o}(r_{c}^{-})\) beneath the CMB. If \({\mathbf {\widetilde {v}}}^{o}(r_{c}^{-})\) is included, then the observation vector y is
$$ \mathbf{y} = \left(\mathbf{y}_{{\widetilde v}},\, \mathbf{y}_{{\widetilde \omega}},\, \mathbf{y}_{{\widetilde b}}\right)^{T} $$
where \(\mathbf {y}_{{\widetilde \omega }}\) includes, as shown in (3) and (4), the spectral coefficients \({\widetilde \omega }_{l}^{m(o)}\) of \({\mathbf {\widetilde {v}}}^{o}\) at \(r_{c}^{-}\). Again, the linearized observation projection (23) is achieved. However, H includes additional subsets:
$$ \mathbf{H} = \left(\mathbf{H}_{v}, \mathbf{H}_{\omega}, \mathbf{H}_{b}, {\mathbf 0}, {\mathbf 0}\right)^{T}, $$
where H v and H ω include only non-zero entries for \({\widetilde v}_{l}^{m(o)}\) and \({\widetilde \omega }_{l}^{m(o)}\) at \(r = r_{c}^{-}\), respectively.
Effective observation error covariance
Since the gain matrix K in (16) depends on the observation error covariance R, we need to determine the effective error covariance \({\mathbf {\widetilde {R}}}\) for \({{\widetilde {b}}_{l}^{m(o)}}\) which can be calculated from those of the Gauss coefficients \({{g_{l}^{m}}}\) and \({{h_{l}^{m}}}\). In this section, we only describe a formal procedure without going into the details.
In geomagnetic field modeling (Jackson et al. 2000; Sabaka et al. 2004; Korte et al. 2005; Olsen et al. 2006; Gillet et al. 2013), the Gauss coefficients, e.g. \({{g_{l}^{m}}}\), can be described in general as
$$ {{g_{l}^{m}}} = {\mathbf S}^{T}(t)\cdot{\boldsymbol\alpha}_{lm}, $$
where S is the vector describing deterministic, model-specific base functions in the time domain, e.g., B-spline functions, and α is the coefficient vector which includes the observation error statistics.
For illustrative purpose, we use the simplest error statistics for our derivation. Assume that geomagnetic observations (and thus α) are unbiased and with known error covariances:
$$ {\boldsymbol{\alpha}}_{lm} = {\boldsymbol{\alpha}}_{lm}^{t} + {\boldsymbol{\epsilon}}_{\alpha},\quad \langle {\boldsymbol{\epsilon}}_{\alpha}\rangle = {\mathbf 0},\quad \left\langle\boldsymbol{\epsilon}_{\alpha}{\boldsymbol{\epsilon}_{\alpha}^{T}}\right\rangle = {\mathbf{C}}_{\alpha}, $$
where \({\boldsymbol {\alpha }}_{\textit {lm}}^{t}\equiv \left \langle {\boldsymbol \alpha }_{\textit {lm}}\right \rangle \) is the truth (expectation) and C α is the observation error covariance matrix of α lm . Thus, by (30),
$$ \begin{array}{l} {{g_{l}^{m}}} = g_{l}^{m(t)} + \epsilon_{g},\quad g_{l}^{m(t)} \ = {\mathbf S}^{T}\cdot{\boldsymbol{\alpha}}_{lm}^{t},\quad \epsilon_{g} = {\mathbf S}^{T}\cdot{\boldsymbol{\epsilon}}_{\alpha}, \\ \left\langle{\epsilon_{g}^{2}}\right\rangle= {\mathbf S}^{T}\cdot{\mathbf C}_{\alpha}\cdot{\mathbf S} \equiv R_{g}^{lm}. \end{array} $$
Similar formulation applies to \({{h_{l}^{m}}}\), but with the error ε h . By (10) and (32), we have
$$\begin{array}{@{}rcl@{}} {b_{l}^{m(o)}}(r) & = & b_{l}^{m(t)}(r) + \epsilon_{b}(r), \end{array} $$
$$\begin{array}{@{}rcl@{}} b_{l}^{m(t)}(r) & =& \frac{{r_{s}^{2}}}{l}\left(\frac{r_{s}}{r}\right)^{l}G_{m}\left(g_{l}^{m(t)}-ih_{l}^{m(t)}\right), \end{array} $$
$$\begin{array}{@{}rcl@{}} \epsilon_{b}(r)\ \ \ \ & =& \frac{{r_{s}^{2}}}{l}\left(\frac{r_{s}}{r}\right)^{l}G_{m} \left(\epsilon_{g}-i\epsilon_{h}\right) \end{array} $$
This leads to
$$ \left\langle\epsilon_{b}\epsilon_{b}^{*}\right\rangle = \left(\frac{{r_{s}^{2}}}{l}\right)^{2}\left(\frac{r_{s}}{r}\right)^{2l}{G_{m}^{2}} \left[\left(R_{g}^{lm}\right)^{2}+\left(R_{h}^{lm}\right)^{2}\right] $$
One can use this equation to evaluate the covariance at any location in the mantle, including r=r d the top of the D ″-layer. If S in (32) is replaced by \(\dot {\mathbf S}\), then we can obtain the covariance \(R_{\dot g}^{lm}\) of the SV
$$R_{\dot g}^{lm} = {\dot {\mathbf S}}^{T}\cdot{\mathbf C}_{\alpha}\cdot{\dot {\mathbf S}}, $$
and therefore the variance of \({\dot b}_{l}^{m(o)}\)
$$\begin{array}{@{}rcl@{}} {\dot b}_{l}^{m(o)}(r) & =& {\dot b}_{l}^{m(t)}(r) + \epsilon_{\dot b}(r), \end{array} $$
$$\begin{array}{@{}rcl@{}} \left\langle\epsilon_{\dot b}\epsilon_{\dot b}^{*}\right\rangle(r) & = & \left(\frac{{r_{s}^{2}}}{l}\right)^{2}\left(\frac{r_{s}}{r}\right)^{2l}{G_{m}^{2}} \left[\left(R_{\dot g}^{lm}\right)^{2}+\left(R_{\dot h}^{lm}\right)^{2}\right]. \end{array} $$
The full error covariance of \({{\widetilde {b}}_{l}^{m(o)}}(r)\) can then be determined from (21), (36), and (38).
In this study, we focus only on (22), i.e., assimilation of the effective observed field \({{\widetilde {b}}_{l}^{m(o)}}\) which matches both the observed field \({b_{l}^{m(o)}}\) and the observed SV \({\dot b}_{l}^{m(o)}\) at the top of the D ″-layer, mainly for two goals: to explore improvements of the assimilation system with the observed SV, such as the model spin-up process and rms of the observed minus forecast (\({\mathcal {O}}\)-\({\mathcal {F}}\)) of the magnetic field; and to understand responses of the core state x to the observed SV, in particular changes of the velocity field v beneath the CMB. Both are critical for determination of the effective velocity field \(\mathbf {\widetilde {v}}\) in (25) and thus for implementation of the more comprehensive observation (28).
The baseline geodynamo model is the MoSST core dynamics model (Kuang and Chao 2003; Jiang and Kuang 2008) for the thermal convection of a Boussinesq electrically conducting fluid in the (rapidly rotating) outer core, confined between the electrically conducting inner core and the D ″-layer. The non-dimensional parameters include the Ekman number E (for viscosity), the magnetic Rossby number R o (for magnetic diffusivity), the modified Prandtl number q κ (for thermal conductivity), and the modified Rayleigh number R th (for buoyancy force). In our assimilation, the ICB and the CMB are assumed impenetrable, stress-free, and fixed heat fluxes. We also select the following parameter values:
$$ E = R_{o} = 1.25\times 10^{-6},\quad R_{o} = 1.0\, \quad R_{th} = 15 R_{th}^{c}, $$
where \(R_{\text {th}}^{c}\) is the critical Rayleigh number for purely thermal convection. In our assimilation, the numerical truncation order is L M =96. There are 20 grid points in D ″-layer (which is 20 km thick with η d =20 η), 80 grid points in the outer core, and 40 grid points in the inner core. With this set of the parameter values, the mean time scale of the dipole field is approximately 0.7 % of the magnetic diffusive time τ η (which is used for the time scaling of the dynamo model) and those for non-dipole components are more than an order of magnitude shorter. These are very similar to those derived from satellite magnetic data. Therefore, in our assimilation, we choose τ η =200,000 years to convert non-dimensional dynamo time to years (with this conversion, the mean time scale of the dipole field is approximately 1400 years).
We consider only the observations for the time period 1900−2000 simply because modern observatory and satellite data provide very high quality \(\left ({{g_{l}^{m}}}, {{h_{l}^{m}}}\right)\) and \(\left ({\dot g}_{l}^{m}, {\dot h}_{l}^{m}\right)\). These coefficients are from gufm1 (Jackson et al. 2000) for 1900–1962 and CM4 (Sabaka et al. 2004) for 1962–2000. We also set L o =8, lower than the highest degrees of the two models. For our research purposes, we carry out three distinct experiments:
Case I: Free-running model (no assimilation)
Case II: Assimilation of \({b_{l}^{m(o)}}\) with (17)
Case III: Assimilation of \({{\widetilde {b}}_{l}^{m(o)}}\) with (22)
Except the differences in the data y in analysis, everything else is identical in the experiments, including the original initial state at 1900. The analysis cycle is Δ t=5 years. By this design, we can identify exactly the causes of changes in the dynamo state x: the differences between the solutions of case I and case II are due to assimilation of the observed field \({b_{l}^{m(o)}}\) and the differences between the solutions of case II and case III are due to the assimilation of the observed SV \({\dot b}_{l}^{m(o)}\). These allow us to understand clearly the responses of the core state to surface observations and their dynamical consequences.
We use a modeled observation error covariance, since the actual error covariances of the field models are not yet available. The model error covariance R is assumed diagonal, with the diagonal elements defined as
$${} R_{lm} = \left\vert\epsilon_{R}(l) {{b_{l}^{m}}}\right\vert^{2},\ \ \epsilon_{R}(l) = \epsilon_{0}(t) + \left[\epsilon_{1}(t)-\epsilon_{0}(t)\right]\frac{l-1}{L_{o}-1}, $$
where ε 0 and ε 1 decreases linearly in time: ε 0 decreases from 0.01 in 1900 to 0.001 in 2000 and ε 1 decreases from 0.3 in 1900 to 0.1 in 2000. These imply that the relative errors in (40) decreases in time, but increases with the degree l.
We would like to point out here that Gillet et al. (2013) provided a global field model which includes a full error covariance of the Gauss coefficients. This model and any future model with specified error statistic knowledge are more appropriate for GDAS. However, we conjecture that (40) is sufficient for our current objectives.
Responses of the magnetic field to SV assimilation
The quantities used to understand the responses of the magnetic field are the (\({\mathcal {O}}\)-\({\mathcal {F}}\)) of the radial magnetic field B r and its SV \({\dot B}_{r}\). Instead of using traditional (\({\mathcal {O}}\)-\({\mathcal {F}}\)), we prefer the following modified definition
$$ \begin{aligned} ({\mathcal{O}}-{\mathcal{F}})_{B}^{2} &= \sum_{1\le l}^{L_{o}}\left\{\left[ \sum_{0\le m\le l}\left|\frac{{b_{l}^{m(o)}}}{b_{1}^{0(o)}}-\frac{{b_{l}^{m(f)}}}{b_{1}^{0(f)}} \right|^{2}\right]\right.\\ &\quad\times\left.\left[ \sum_{0\le m\le l}\left|\frac{{b_{l}^{m(o)}}}{b_{1}^{0(o)}}\right|^{2}\right]^{-1}\right\} \end{aligned} $$
at r=r d , i.e., the misfit normalized by the observed field strength for a given degree l. Replacing \({{b_{l}^{m}}}\) by \({\dot b}_{l}^{m}\) in (41), we have \(({\mathcal {O}}\text {-}{\mathcal {F}})_{\dot B}\) of the SV. This modified \(({\mathcal {O}}\text {-}{\mathcal {F}})\) can tell us more accurately how close is the forecast to observation because it eliminates the effect of changes in the magnitude of the individual spectral coefficients. Our modified \(({\mathcal {O}}\text {-}{\mathcal {F}})\) in (41) is different from other quantities used for measuring the difference (or misfit) between the observation y and the forecast H·x f, such as the misfit normalized by the error covariances (e.g., Aubert 2014). Since, in our assimilation, the modeled error covariances are proportional in magnitude to the observed field strength; (41) is actually very similar to the normalized misfit.
Figure 1 is the \(({\mathcal {O}}\text {-}{\mathcal {F}})_{B}\) of Case II (dashed lines) and Case III (solid lines). From this figure, we can observe clearly that their magnitudes in Case III are approximately 30 % smaller than those in Case II over the entire assimilation period, demonstrating a substantial improvement in forecast accuracies with the SV assimilation (21) and (22). Similar improvement can be also observed from \(({\mathcal {O}}\text {-}{\mathcal {F}})_{\dot B}\) shown in Fig. 2.
The SV assimilation also helps accelerate the dynamo model spin-up process. For example, we can observe from Fig. 1 that the time variations of \(({\mathcal {O}}\text {-}{\mathcal {F}})_{B}\) are nearly identical in both cases: they decay nearly monotonically over much of the assimilation period before leveling off in the last 20 years (from 1980 to 2000). But \(({\mathcal {O}}\text {-}{\mathcal {F}})_{\dot B}\), as shown in Fig. 2, are very different in the two cases: in Case II, it increases first from 1900 to 1940 and only starts to decay continuously in the last 20 years. In Case III, however, \(({\mathcal {O}}\text {-}{\mathcal {F}})_{\dot B}\) decays almost monotonically in time, except two small surges around 1940 and 1980. This implies that the dynamo core state x f responds stronger to the SV assimilation. In other words, the SV assimilation helps to accelerate the model spin-up process.
To better understand how do the forecasts \({b_{l}^{m(f)}}\) and \({\dot b}_{l}^{m(f)}\) respond to the observations y b in (17) and \(\mathbf {y}_{\widetilde b}\) in (22), we examine first the \(({\mathcal {O}}\text {-}{\mathcal {F}})\) for individual degrees. In Fig. 3 are \(({\mathcal {O}}\text {-}{\mathcal {F}})_{B}\) for the degrees l≤6. Improvements are clearly shown in all six degrees, as all values are smaller in Case III than those in Case II. But we can also observe different patterns among these degrees. For example, \(({\mathcal {O}}\text {-}{\mathcal {F}})_{B}\) for the odd degrees (l=1,3,5) increase in magnitude again from around 1980 to the end of the assimilation period. But there is no such clear reversing trend in those for the even degrees (l=2,4,6): it decreases monotonically for l=2; for l=4, it oscillates with a damping amplitude; and the \(({\mathcal {O}}\text {-}{\mathcal {F}})\) for l=4 remains nearly constant after the rapid decay in the first two analysis cycles.
The \(({\mathcal {O}}\text {-}{\mathcal {F}})_{B}\) of the first six spherical harmonic degrees in case II (dashed lines) and case III (solid lines)
As shown in Fig. 4, the difference between the odd and even degrees of \(({\mathcal {O}}\text {-}{\mathcal {F}})_{\dot B}\) is even more significant. There is still a strong surge in magnitude for l=3 around 1980 in the both cases. But the reduction for l=5 is minimal. In particular, it does not decay monotonically in time in either case. These differences may indicate potential inconsistencies between the core dynamics of the model and the time variation of the Gauss coefficients. We will discuss this again later in this paper.
Similar to Fig. 3, but for \(({\mathcal {O}}\text {-}{\mathcal {F}})_{\dot B}\)
Responses of the velocity field to SV assimilation
Why does the dynamo model respond faster and stronger in case III than in case II? We can find at least partial answers from the difference between the free-running model solutions x M (case I) and the forecasts x f in cases II and III, in particular the differences in the velocity field v beneath the CMB, because they are the direct consequences of the magnetic induction (12). The knowledge is also very important for obtaining the "effective" observed velocity field (25) for future studies.
Since in our geodynamo model, the CMB is impenetrable and is free-slip, the radial velocity v r =0, and, by (3) and (4), the horizontal velocity v h depends on \(\partial {{v_{l}^{m}}}/\partial r\) and \({{\omega _{l}^{m}}}\) at r=r c . Therefore, it is very convenient to examine the following two variables beneath the CMB:
$${} {\small{\begin{aligned} v_{r}' & \equiv \frac{\partial v_{r}}{\partial r} = \sum_{0\le m\le L}^{L_{M}}\frac{l(l+1)}{{r_{c}^{2}}}\frac{\partial {{v_{l}^{m}}}}{\partial r}{Y_{l}^{m}}\left(\theta,\phi\right) + C.C. \end{aligned}}} $$
$${} {\small{\begin{aligned} \omega_{r} & \equiv \left(\nabla\times\mathbf{v}\right)_{r} = \sum_{0\le m\le L}^{L_{M}}\frac{l(l+1)}{{r_{c}^{2}}}{{\omega_{l}^{m}}} {Y_{l}^{m}}\left(\theta,\phi\right) + C.C., \end{aligned}}} $$
where v r′ is poloidal and describes the up-and-down welling and ω r is toroidal and describes the differential rotation. The rms differences (\({\mathcal {M}}\text {-}{\mathcal {F}}\)) between the two variables of the forecast x f (cases II and III) and of the free-running model x M (case I) can be used to quantify the responses of the core flow to the assimilation of surface observations:
$${} \begin{aligned} (\mathcal{M}\text{-}{\mathcal{F}})_{v_{P}} & \equiv \|v_{r}^{\prime M} - v_{r}^{\prime f}\|_{2}\\ &= \left[ \sum_{0\le m\le l}^{L_{M}}\frac{l^{2}(l+1)^{2}}{{r_{c}^{4}}}\left| \frac{\partial v_{l}^{m(M)}}{\partial r}- \frac{\partial v_{l}^{m(f)}}{\partial r}\right|^{2}\right]^{1/2} \end{aligned} $$
$${} \begin{aligned} ({\mathcal{M}}\text{-}{\mathcal{F}})_{v_{T}} & \equiv \|{\omega_{r}^{M}} - {\omega_{r}^{f}}\|_{2}\\ &= \left[ \sum_{0\le m\le l}^{L_{M}}\frac{l^{2}(l+1)^{2}}{{r_{c}^{4}}}\left| \omega_{l}^{m(M)}-\omega_{l}^{m(f)} \right|^{2}\right]^{1/2} \end{aligned} $$
In the above equations, ∥·∥2 is the L 2−norm (or rms) over the CMB.
In Fig. 5 are the non-dimensional (with the scaling factor 5×10−6 year−1 for dimensional values) ∥v r′∥2 (red) and ∥ω r ∥2 (blue) of the free-running model (case I). As shown in the figure, v r′ increases slightly in magnitude in the assimilation period and ω r remains flat. But, the rms differences \(({\mathcal {M}}\text {-}{\mathcal {F}})_{v_{P}}\) (shown in Fig. 6) and \(({\mathcal {M}}\text {-}{\mathcal {F}})_{v_{T}}\) (shown in Fig. 7) increase in time, i.e., a growing divergence between the forecast state x f and the free-running model state x M.
The non-dimensional ∥v r′∥2 (red) and ∥ω r ∥2/10 (blue) beneath the CMB \(r=r_{c}^{-}\) from the free-running model solutions (case I). The dimensional values can be obtained with the scaling factor 5×10−6 year−1
The \(({\mathcal {M}}\text {-}{\mathcal {F}})\) of the poloidal velocity field v r′ as defined in (44). The dashed lines are the results without SV assimilation (case II) and the solid lines are those with the SV assimilation (case III)
Similar to Fig. 6, but for The \(({\mathcal {M}}\text {-}{\mathcal {F}})\) of the toroidal velocity ω r as defined in (45)
From Figs. 6 and 7, we can also observe that \(({\mathcal {M}}\text {-}{\mathcal {F}})\) of case III (the solid lines) are slightly larger than those of case II (dashed lines), implying that x f moves away from x M faster with the SV assimilation (22), another demonstration of improved model spin-up with the SV assimilation. However, the differences are much less significant than those of the magnetic field. This suggests the need for the effective observed velocity field \({\mathbf {\widetilde {v}}}^{o}\) to increase further \(({\mathcal {M}}\text {-}{\mathcal {F}})\) of the velocity field and thus to expedite the model spin-up process.
To aid the future study of determining the effective core flow from the observed SV via (25), we need to understand better the details of \(({\mathcal {M}}\text {-}{\mathcal {F}})\), e.g., their distributions in the spectral space defined by the spherical harmonic degrees l and orders m. We shall pay special attention to their distributions in l, i.e., the summation of the terms in (44-45) with 0≤m≤l for a given degree l, and their distributions in m, i.e., the summation of the terms in (44-45) with m≤l≤L M for a given order m. Since, as shown in Figs. 6 and 7, the differences between the two cases are very small, we can focus only on case III without loss of generality.
In Fig. 8 is the distribution of \(({\mathcal {M}}\text {-}{\mathcal {F}})_{v_{P}}\) in the degree l and in Fig. 9 is its distribution in the order m. From the figures, we can find that \(({\mathcal {M}}\text {-}{\mathcal {F}})_{v_{P}}\) varies substantially in the spectral spaces. As shown in Fig. 8, the differences for the degrees 15≤l≤35 increase the fastest in time, and their magnitudes are the largest at the end of the assimilation period, with the peak at l=20. The differences are much smaller and grows slower in time for the degrees l≤5 and l≥40. But, as shown in Fig. 9, the distribution in m is more broadband: the differences for 5≤m≤35 increase rapidly in time and reach comparable values in magnitude at the end of the assimilation period. However, \(({\mathcal {M}}\text {-}{\mathcal {F}})\) for m≤4 are very different: they remain small and nearly unchanged throughout the entire assimilation. These suggest that the responses of the poloidal velocity is dominantly non-axisymmetric (m>0).
The distribution of \(({\mathcal {M}}\text {-}{\mathcal {F}})_{v_{P}}\) in spherical harmonic degrees l with the SV assimilation (case III)
Similar to Fig. 8, but for the distribution of \(({\mathcal {M}}\text {-}{\mathcal {F}})_{v_{P}}\) in spherical harmonic orders m
The distribution of \(({\mathcal {M}}\text {-}{\mathcal {F}})_{v_{T}}\) of the toroidal velocity, as shown in Figs. 10 and 11, displays both similar and distinct characteristics. Its distribution in l is very similar to that of \(({\mathcal {M}}\text {-}{\mathcal {F}})_{v_{P}}\), except that it peaks at a higher degree l=30. But its distribution in m (Fig. 11) is very different: the differences for m≤20 remain comparable in both the magnitude and the time increasing rate. But they decay rapidly for larger m. It should be pointed out in particular that, opposite to \(({\mathcal {M}}\text {-}{\mathcal {F}})_{v_{P}}\) (in Fig. 9), \(({\mathcal {M}}\text {-}{\mathcal {F}})_{v_{T}}\) of the axisymmetric toroidal velocity (m=0) remains very large, implying that the axisymmetric toroidal flow is very sensitive to the surface observations.
Similar to Fig. 8, but for \(({\mathcal {M}}\text {-}{\mathcal {F}})_{v_{T}}\)
Discussions and conclusions
In this study, we have examined the consequences of assimilating the observed SV on geomagnetic forecasts and on the responses of the dynamo core state. We argued that,because geomagnetic data sampling frequencies are several orders of magnitude higher than those of the SV, the geomagnetic field and its SV are concurrently measured. We further demonstrated that the observed SV provides unique knowledge of the magnetic field and the velocity field in the core. Thus, assimilations of the observed field and of the observed SV are necessary and are not redundant.
In this study, we incorporate the observed SV into the observation vector y via introducing the effective poloidal field \({{\widetilde {b}}_{l}^{m(o)}}\) (21) in the D ″-layer, which is then used in the sequential assimilation algorithm (15). We designed the three experiments to identify the impact of SV assimilation: a free-running model dynamo simulation (case I), an experiment with the assimilation of the observed field (case II), and an experiment with the assimilation of both the observed field and its SV. The relative \(({\mathcal {O}}\text {-}{\mathcal {F}})\) of the field and SV, defined in (41), at the top of the D ″-layer are used to measure the forecast accuracies; the (\({\mathcal {M}}\text {-}{\mathcal {F}}\)) of the poloidal velocity field (44) and of the toroidal velocity field (45) beneath the CMB are used to characterize the responses of the core state to the SV assimilation.
The results of our experiments demonstrate clearly that the SV assimilation with (21) improves significantly the geomagnetic forecast accuracies since, as shown in Figs. 1 and 2, both \(({\mathcal {O}}\text {-}{\mathcal {F}})_{B}\) and \(({\mathcal {O}}\text {-}{\mathcal {F}})_{\dot B}\) in case III are more than 20 % smaller than those in case II. In particular, the improvements occur to all degrees, as shown in Figs. 3 and 4. The nearly monotonic decay in time of \(({\mathcal {O}}\text {-}{\mathcal {F}})_{\dot B}\) in case III (Fig. 2) shows clearly that the SV assimilation accelerates the spin-up of the dynamo model.
The improvement by the SV assimilation can be also seen from the differences \(({\mathcal {M}}\text {-}{\mathcal {F}})_{v_{P}}\) and \(({\mathcal {M}}\text {-}{\mathcal {F}})_{v_{T}}\) between the free-running model state and those of the assimilations. As shown in Figs. 6 and 7, these differences grow rapidly in time, showing an accelerated departure of the core state with assimilation from the free-running model state. The differences in case III are slightly larger than those in case II, further demonstrating the improvement brought by the SV assimilation, though such increment is less significant that those in the \(({\mathcal {O}}\text {-}{\mathcal {F}})\) of the magnetic field (Figs. 1 and 2).
Our results have further implications. First, even with the help of (21), the dynamo model is still not fully spun up. For example, though \(({\mathcal {O}}\text {-}{\mathcal {F}})_{\dot B}\) decreases monotonically in time, the SV forecast is still very far away from the observations, as \(({\mathcal {O}}\text {-}{\mathcal {F}})_{\dot B} \approx {\mathcal {O}}(1)\) for all degrees (see Fig. 4). This can be shown further by the continuously growing differences \(({\mathcal {M}}\text {-}{\mathcal {F}})_{v_{P}}\) and \({\mathcal {M}}\text {-}{\mathcal {F}})_{v_{T}}\) between the forecast velocity field and that of the free-running model (see Figs. 6 and 7). In addition, the differences at the end of the assimilation period are still very small, approximately 10 % in the magnitude of the velocity field of the free-running model (Fig. 5).
These suggest that much larger velocity differences \(({\mathcal {M}}\text {-}{\mathcal {F}})_{v_{P}}\) and \(({\mathcal {M}}\text {-}{\mathcal {F}})_{v_{T}}\) are needed over a shorter assimilation period for expediting the model spin-up. Assimilation of the effective observed velocity (25) could be an answer. However, as discussed earlier, (25) is an underdetermined system, since both \(\mathbf {v}_{h}^{(o)}\) and \({{\widetilde {b}}_{l}^{m(o)}}\) (more specifically, \(\partial ^{2}{{\widetilde {b}}_{l}^{m(o)}}/\partial r^{2}\)) are unknown beneath the CMB. Thus, the responses of v h to the SV assimilation, e.g., \(({\mathcal {M}}\text {-}{\mathcal {F}})_{v_{P}}\) (in Figs. 8 and 9) and \(({\mathcal {M}}\text {-}{\mathcal {F}})_{v_{T}}\) (in Figs. 10 and 11) are needed to determine \(\mathbf {v}_{h}^{(o)}\). For example, as shown in the two figures, the non-axisymmetric (m>0) poloidal velocity \(v_{l}^{m\prime }\) around the degree l=20 and the toroidal velocity \({\omega _{l}^{m}}\) around the degree l=30 and order m≤20 should be given more attention, as they are most sensitive to the surface observations.
An alternative answer could be the core states inverted from the surface observations and dynamo solutions, such as those of Aubert (2013, 2014). These can be used as the analysis of the assimilation system. But cautions should be taken with this approach. For example, the inverted velocity field beneath the CMB is actually derived with the observed field and SV and the magnetic diffusion of the dynamo state (Aubert 2014). This could potentially lead to dynamical inconsistencies as well as uncertainties in error statistics.
Our results also show several new features that may have implications to field modeling and to core flow inversion. One new knowledge is from the time variation of \(({\mathcal {O}}\text {-}{\mathcal {F}})_{\dot B}\). As shown in Fig. 4, \(({\mathcal {O}}\text {-}{\mathcal {F}})_{\dot B}\) of the odd degrees (l=1,3,5) are significantly different from those of the even degrees (l=2,4,6): the values of the even degrees decay nearly monotonically in time; but those of the odd orders show either spikes (for l=1,3) during the assimilation or even increase over time (for l=5). These even-odd degree disparities suggest inconsistencies between the model and the observations. These inconsistencies could be entirely due to numerical dynamo model which may include a magnetic induction different from those in the Earth's outer core or may include some mechanisms resulting in different symmetry properties of the core state. But the inconsistencies could also come from possible biases in the field models that are not included in the observation error covariances. For example, ionospheric ring current generated field (an external field component) contributes dominantly to the Gauss coefficients of degrees l=1,3,5 and varies on time scales comparable to those of SV, e.g., the solar activity cycles (Sabaka et al. 2015). Model biases exist if this part of the signals is not well separated from those of the core field. This is potentially an area for application of geomagnetic data assimilation.
The core fluid flow responses, i.e., the differences \(({\mathcal {M}}\text {-}{\mathcal {F}})_{v_{P}}\) and \(({\mathcal {M}}\text {-}{\mathcal {F}})_{v_{T}}\) between the forecast velocity field \({\mathbf {v}_{h}^{f}}\) and the \(\mathbf {v}^{M}_{h}\) of the free-running model (see Figs. 8, 9, 10, and 11), from our experiments could also help the inversion of core flow from the observed SV. For example, the different characteristics in \(({\mathcal {M}}\text {-}{\mathcal {F}})_{v_{P}}\) and \(({\mathcal {M}}\text {-}{\mathcal {F}})_{v_{T}}\) suggest that the poloidal velocity field and the toroidal velocity field could be treated separately in the core flow inversion. It should be pointed out that purely toroidal core flow approximations were used in previous studies (e.g., Bloxham et al. 2002; Olsen and Mandea 2008). The strong responses of the high-degree core flow (l≈20 for the poloidal flow and l≈30 for the toroidal flow) to the observed SV (for l≤8) indicate that higher degree velocity field should be included in the core flow inversion. For example, one would normally expect that due to nonlinear effects, i.e., the quadratic terms in Navier-Stokes equation and the induction equation, the core flow up to the degrees twice as much as that of the SV should be sufficient for the core flow inversion (as in Aubert 2013). But our results show that time evolution of the core flow leads to the strongest responses for the degrees more than triple of the maximum degree of the SV. Therefore, inversion of time-dependent core flow from the observed SV (up to degree 13) should include high-degree (l>40) spectral coefficients. Of course, the core flow response may also suggest that the free-running model solutions do not provide a good description of the small-scale dynamical processes in the core. And assimilation of geomagnetic observations could substantially improve our understanding of these small-scale processes. Regardless, our results suggest that unobservable small-scale dynamical processes are very important in interpretation and prediction of fast SV observable at the Earth's surface.
Again, we should point out that our results could be improved in several areas: assimilation with other geodynamo models different from (39) and with more sophisticated assimilation algorithms and field models with more accurate error statistics. The former will build up an ensemble of assimilation results to provide more accurate statistics on the characteristics of the core flow to the observed SV. Improvements from the latter are also obvious. For example, we anticipate more accurate estimation of \(({\mathcal {O}}\)-\({\mathcal {F}})\) for both the field and the SV, and better assessment of the core state responses if a full ensemble approach is used for the covariance P f, and a more appropriate observation error covariance, e.g., those determined by Gillet et al. (2013), than (40) used in this study. Regardless, our assimilation experiments have shown clearly the importance of SV assimilation and the improvements that the SV assimilation brought to forecast accuracies and to model spin-up processes.
CMB:
core-mantle boundary
GDAS:
geomagnetic data assimilation
ICB:
inner core boundary
SV:
secular variation
Aubert, J (2013) Flow through the Earth's core inverted from geomagnetic observations and numerical dynamo models. Geophys. J. Int.192: 537–56.
Aubert, J (2014) Earth's core internal dynamics 1840–2010 imaged by inverse geodynamo modeling. Geophys J Int.197: 1321–34.
Aubert, J, Fournier A (2011) Inferring internal properties of Earth's core dynamics and their evolution from surface observations and a numerical geodynamo. Nonlin Processes Geophys.18: 657–74.
Backus, GE (1968) Kinematics of geomagnetic secular variation in a perfectly conducting core. Phil Trans R Soc LondonA263: 239–66.
Bloxham, J, Zatman S, Dumberry M (2002) The origin of geomagnetic jerks. Nature420: 65–8.
Christensen, U. R, Aubert J, Hulot G (2010) Conditions for Earth-like geodynamo models. Earth Planet Sci Lett.296: 487–96.
Christensen, U. R, Wardingsky I, Aubert J, Lesur V (2012) Timescales of geomagnetic secular acceleration in satellite field models and geodynamo models. Geophys J Int.190: 243–54.
Evensen, G (1994) Sequential data assimilation with a nonlinear quasi-geostrophic model using Monte Carlo methods to forecast error statistics. J Geophys Res.99: 10143–62.
Fournier, A, Eymin C, Alboussiere T (2007) A case for variational geomagnetic data assimilation: insights from a one-dimensional, nonlinear, and sparsely observed MHD system. Nonlinear Process Geophys14: 163–80.
Fournier, A, Hulot G, Jault D, Kuang W, Tangborn A, Gillet N, et al (2010) An introduction to data assimilation and predictability in geomagnetism. Space Sci Rev. doi:10.1007/s11214-010-9669-4.
Fournier, A, Aubert J, Thébault E (2011) Inference on core surface flow from observations and 3-D dynamo modeling. Geophys J Int.186: 118–36.
Fournier, A, Nerger L, Aubert J (2013) An ensemble Kalman filter for the time-dependent analysis of the geomagnetic field. Geochem Geophys Geosys. doi:10.1002/ggge.20252.
Gillet, N, Jault D, Finlay CC, Olsen N (2013) Stochastic modeling of the Earth's magnetic field: inversion for covariances over the observatory era. Geochem Geophys Geosys. doi:10.1002/ggge.20041.
Glatzmaier, GA, Roberts PH (1995) A three-dimensional self-consistent computer simulation of a geomagnetic field reversal. Nature377: 203–6.
Holme, R (2007) Large-scale flow in the core. In: Olson P (ed)Core dynamics, treaties on geophysics.. Elsevier, Amstedam.
Hou, AY, Ledvina DV, da Silva AM, Zhang SQ, Joiner J, Atlas RM, et al (2000) Assimilation of SSM/I-derived surface rainfall and total precipitable water for improving the GEOS analysis for climate studies. Mon Weath Rev.128: 509–37.
Hulot, G, LeMouël JL (1994) A statistical approach to the Earth's main magnetic field. Phys Earth Planet Inter.82: 167–83.
Jackson, A, Jonkers AR, Walker MR (2000) Our centuries of geomagnetic secular variation from historical records. Phil Trans R Soc Lond.A358: 957–90.
Kageyama, A, Sato T (1997) Generation mechanism of a dipole field by a magnetohydrodynamic dynamo. Phys Rev E.55: 4617–26.
Kalnay, E (2003) Atmospheric modeling, data assimilation and predictability. Cambridge University Press, Cambridge, U.K.
Korte, M, Genevey A, Constable CG, Frank U, Schnepp E (2005) Continuous geomagnetic field models for the past 7 millennia: 1. A new data compilation. Geochem Geophys Geosys. doi:10.1029/2004GC000800.
Korte, M, Constable CG, Donadini F, Holme R (2011) Reconstructing the Holocene geomagnetic field. Earth Planet Sci Lett.312: 497–505.
Kuang, W, Bloxham J (1997) An Earth-like numerical dynamo model. Nature389: 371–4.
Kuang, W, Tangborn A (2011) Interpretation of core field models. In: Mandea M Korte M (eds)Geomagnetic observations and models.. Springer, Heidelberg.
Kuang, W, Tangborn A, Jiang W, Liu D, Sun Z, Bloxham J, Wei Z (2008) MoSST −DAS: The first generation geomagnetic data assimilation framework. Comm Comp Phys.3: 85–108.
Kuang, W, Tangborn A, Wei Z, Sabaka TJ (2009) Constraining a numerical geodynamo model with 100 years of surface observations. Geophys J Int.179: 1458–68.
Kuang, W, Wei Z, Holme R, Tangborn A (2010) Prediction of geomagnetic field with data assimilation: a candidate secular variation model for IGRF-11. Earth Planet Space.62: 775–85.
Larmor, J (1919) How could a rotating body such as the Sun become a magnet?Rep Br Assoc87: 159–60.
Li, D, Jackson A, Livermore PW (2011) Variational data assimilation for the initial value dynamo problem. Phys Rev E. doi:10.1103/PhysRevE.84.056321.
Li, D, Jackson A, Livermore PW (2014) Variational data assimilation for a forced, inertia-free magnetohydrodynamic dynamo model. Geophys J Int.199: 1662–76.
Liu, D, Tangborn A, Kuang W (2007) Observing system simulation experiments in geomagnetic data assimilation. J Geophys Res.doi:10.1029/2006JB004691.
Nilsson, A, Holme R, Korte M, Suttie N, Hill M (2014) Reconstructing Holocene geomagnetic field variation: new methods, models and implications. Geophys J Int. doi:10.1093/gji/ggu120.
Olsen, N, Mandea M (2008) Rapidly changing flows in the Earth's core. Nature Geosci.1: 390–4.
Olsen, N, L uhr H, J S. T, Mandea M, Rother M, Tøffner-Clausen L, Choi S (2006) CHAOS-a model of the Earth's magnetic field derived from CHAMP, Ørsted, and SAC-C magnetic satellite data. Geophys J Int.166: 67–75.
Olsen, N, Lühr H, Finlay C. C, Sabaka T. J, Michaelis I, Rauberg J, Tøffner-Clausen L (2014) The CHAOS-4 geomagnetic field model. Geophys J Int.197: 815–27.
Roberts, PH (1992) Geomagnetism. Encyclopedia Earth Syst Sci.2: 277–94.
Roberts, PH, Scott S (1965) On analysis of the secular variation, 1, A hydromagnetic constraint: theory. J Geomag Geoelec.17: 137–51.
Sabaka, TJ, Olsen N, Purucker M (2004) Extending comprehensive models of the Earth's magnetic field with Ørsted and champ data. Geophys J Int.159: 521–47.
Sabaka, T. J, Olsen N, Tyler R. H, Kuvshinov A (2015) CM5, a pre-Swarm comprehensive geomagnetic field model derived from over 12 year of CHAMP, Ørsted, SAC-C and observatory data. Geophys J Int.200: 1596–626.
Stacey, FD (1992) Physics of the Earth. Brookfield Press, Kenmore, Australia.
Sun, Z, Kuang W (2015) An ensemble algorithm based component for geomagnetic data assimilation. Terr Atmos Ocean Sci. doi:10.3319/TAO.2014.08.19.05.
Sun, Z, Tangborn A, Kuang W (2007) Data assimilation in a sparsely observed one-dimensional modeled MHD system. Nonlin Proc Geophys.14: 181–92.
Tangborn, A, Kuang W (2015) Geodynamo model and error parameter estimation using geomagnetic data assimilation. Geophys J Int. doi:10.1093/gji/ggu409.
Zatman, S, Bloxham J (1997) Torsional oscillations and the magnetic field within the Earth's core. Nature388: 760–3.
This research is funded by grants from NASA Earth Surface and Interior Program (NNZ09AK70G) and from NSF (EAR-0757880). We also thank NASA NAS for support on high-end scientific computation.
Planetary Geodynamics LaboratoryNASA Goddard Space Flight Center, 8800 Greenbelt Road, Greenbelt, 20771, MD, USA
Weijia Kuang
Joint Center for Earth System Technologies, University of Maryland, Baltimore County, 1000 Hilltop Circle, Baltimore, 21250, MD, USA
Andrew Tangborn
Search for Weijia Kuang in:
Search for Andrew Tangborn in:
Correspondence to Weijia Kuang.
WK proposed and designed the study and carried out the assimilation experiments. AT provided error covariance models and mathematical derivations for the sequential assimilation algorithm used in the experiments. Both authors read and approved the final manuscript.
Kuang, W., Tangborn, A. Dynamic responses of the Earth's outer core to assimilation of observed geomagnetic secular variation. Prog. in Earth and Planet. Sci. 2, 40 (2015) doi:10.1186/s40645-015-0071-4
Geomagnetic field
Core flow
Multidisciplinary Researches on Deep Interiors of the Earth and Planets


After many years of dedicated work, the community from the Chocolatão slum in Porto Alegre have new homes and for many its brought new livelihoods and opportunities. May 13 2011 marked the inauguration of the new Chocolatão, home to around 200 families or 800 people.
Res idencial Nova Chocolatão was opened by Prefeito Jose Fortunati (left) and the event included the project's many partners, including Mario Porto Fonseca of Solucoes Usimanas, lead partner in the building of the community's new recycling facility. Cities Programme manager, Elizabeth Ryan had the privilege of participating in the event and says the time with the community during this moment of great change was incredibly moving, the results a credit to the collective will of many to bring dignity to the lives of the people of Chocolatão. Community leader, Fernanda (left) melted the crowd with her empassioned speech at the inauguration ceremony. She will manage Nova Chocolatão's beautiful new childcare facility.
Vania Goncalves de Souva from the Prefeitura de Porto Alegre has worked with the community for four years and is committed to continuing this work. Vania is pictured (left) on the last day of Vila Chocolatão with Mario, who has lived at Chocolatão for thirty years. She says "resettlement or settlement is an important, a fundamental big step towards the emancipation of the community …the new house is important because it improves the living conditions. But the point is… to continue capacitating, empowering the community to bring change into their lives to be able to cope with the rights and duties in the new life. We know, even after the resettlement we have to continue working with the community to seek together the true social inclusion". Vania is currently working with researchers from the Cities Programme and Australian University, RMIT to help map the history of the project and develop a model for application in other settings.
For some the improvement has been immediate. Bete – pictured below in front of her shop in Vila Chocolatão and then receiving the keys to her new shop from Prefeito Fortunati – is very pleased with the change in her life. She recently reported delight at the beauty of rain falling on the window of her new shop – a first time experience.
Pictured below is long time Chocolatao resident who now has a wheelchair accessible house and sealed roads at Nova Chocolatão. He is also pictured with Cezar Bussato, Secretary of Political Coordination and Local Governance, Prefeitura de Porto Alegre and project champion. Cezar says Chocolatao's inclusive and consultative model will now be applied to other communities the City is working with. This includes focusing on health, employment and training as well as housing.
The new recycling centre at Nova Chocolatão is a major source of employment for the community who were traditionally catadores, or street pickers. The City has changed its garbage collection and disposal practices so that refuse is delivered to the workers at Chocolatão instead of them collecting from the streets and bringing home the rubbish to sort, a major health and safety risk which resulted in tragic fire a number of years ago. All facilities will be managed and operated by Chocolatão community members.
A network, 'Vila Chocolatão Sustainability Network' was an important governing collaborative force in the years of preparation to move. It included members of health, education and law services, the Chocolatao community and the private sector. Teresa, another community leader is pictured below with friend Mario and Marcio Milleto Mostardeiro, from the City of Porto Alegre in the old community space that was used for network meetings when it was wet (they were otherwise in the adjacent park) and childcare. Teresa is also pictured in her new home at Nova Chocolatão with husband and adopted son on the day of the inauguration.
The relocation has not been without its challenges. A number of families moved to Vila Chocolatao after the project's census date and a range of alternative accommodation was required to be negotiated for these people.


text:
Superior Ave at Public Square Closing Causes Increased Safety & Security Risks, Study Concludes
Commander Harold Pretel at podium [Annie Wu / ideastream]
Closing Superior Avenue through Public Square in downtown Cleveland has caused increased risk to RTA buses, pedestrians and cars in the area. That's the conclusion of a safety study commissioned by the Greater Cleveland Regional Transit Authority after the Mayor's office and RTA agreed to keep the road closed following the Republican National Convention.
A separate traffic study conducted by Parsons Brinckerhoff finds that travel time for buses around Public Square can be reduced by 72 seconds during morning rush hour and 60 seconds during afternoon peak travel time if Superior is open for buses. It would also save RTA more than $800,000 in annual operating costs.
The traffic study also looked at a second scenario: Superior remains closed and traffic lights at 16 intersections would be semi-actuated. That option would reduce bus travel time by 20 seconds in the morning and 27 seconds in the evening as well as $223,000 in RTA annual operating costs.
The safety study conducted by K & J Consulting also concludes Public Square is vulnerable to terrorism whether Superior is open or closed and recommends further study by local and federal agencies.
"I think the one major point that's important was if there was a major event around Public Square like the tree lighting, like the Fourth of July – the plan always was that it would be at that point in time closed to RTA buses and secured on the perimeters," says RTA General Manager Joe Calabrese.
RTA provided a summary of conclusions from each study. Cleveland city officials say they have yet to receive the full safety study, and they are not convinced the city's security and mitigation concerns have been adequately addressed.
Commander Harold Pretel is Cleveland's in-house homeland security expert. He's scheduled to meet with RTA and the study consultants on Wednesday afternoon.
"It would be irresponsible to move forward and just open the Square without addressing any of those potential threats and our vulnerability and without discussing any mitigation strategies," says Pretel.
RTA could owe the federal government $12 million dollars for violating the terms of a transit grant if Superior Avenue is not re-opened by February 21st.

redpajama_set_name: RedPajamaCommonCrawl


text:
Sam Kouvaris: At 90, Tony Trabert savors tennis memories
Sam Kouvaris
As sports fans we're all getting used to the different timelines that are happening in 2020. The NBA, NHL, college and pro football have all had jumbled schedules. In golf they played the U.S. Open in September and the Masters is in November.
In tennis they cancelled Wimbledon but the French Open will conclude today, five months after its normal spot on the sports calendar.
I called my friend Tony Trabert, who lives in Ponte Vedra, to talk about what was going on at Roland Garros in Paris this week and he casually mentioned, "I won that 70 years ago."
Whoa. "Seventy years ago," I said to myself. "That's a big number."
Trabert turned 90 in August and his Hall of Fame career has been well documented: Instrumental in America's Davis Cup success and eventually the Davis Cup captain, Trabert won ten Grand Slam titles, back-to-back French Championships and three of the four majors in 1955. He was on the cover of Sports Illustrated.
But an offer and a little help from a fellow Cincinnatian in 1950 really got things started.
Already an accomplished junior player, Trabert won the Ohio state high school singles championships three straight years, He was at the University of Cincinnati when his idol, Bill Talbert, called to offer him a chance to play as his doubles partner in Europe that spring and summer.
"Talbert had taken me under his wing when I was twelve," Trabert recalled. "He saw me swinging at volleys and came out on the court to help."
Trabert asked both the president of the University of Cincinnati and his parents about going. The school president said, "You'll learn more in three months there than you would in three months going over books at school."
Talbert told the national tennis governing body — then the United States Lawn Tennis Association — that he was taking the fledgling star, Trabert, to Europe as his doubles partner and asked for some financial assistance.
"They said no," Trabert recalled with a laugh. "They didn't think I was going to be any good."
Talbert was determined. As one of the top players in the world, he had never played in Europe and told Tony he'd figure it out. So Trabert, Talbert and his wife Nancy traveled through Europe in the spring and summer of 1950 playing in all of the big tournaments as doubles partners.
And winning everything.
They flew to Nice and played at Monte Carlo, winning in five sets in the final.
"We played Jaroslav Drobny, a lefthander from Czechoslovakia, and his partner in the final," Trabert remembered in vivid detail." He'd kick his serve to my forehand, and I'd never seen that before. I finally figured that out and we won in the fifth set."
They played in Monte Carlo and won in April. They went on the road and won in Nice, in Rome at the Italian Championships,(where Talbert lost to Drobny in the singles finals) the Paris City Championships, the French, and at Queens Club in London on grass.
The only tournament they didn't play as doubles partners was at Wimbledon.
Talbert felt that he owed it to his regular doubles partner, Gardnar Mulloy, to play at Wimbledon. Talbert had won the U.S. Championships (now the US Open) with Mulloy four times before 1950. At Wimbledon that year, Talbert and Mulloy lost in the third round. Meanwhile, Trabert and his partner, Budge Patty, the eventual singles champion, made it to the semifinals losing to the eventual champions.
After Wimbledon, they came back to the States and Trabert returned to the University of Cincinnati.
"That was a big jump start," Trabert said. "By far the most important thing, tennis wise, to happen to me. I saw more silverware in front of me at dinner one night in Europe than we had in our whole house. I learned so much, and not just about tennis. Writing thank-you notes, how to dress, all of that."
"Each tournament gave you something for your travel and per diem," he added. "When I finally made a little money I told Bill I wanted to repay him. He said just do something nice for another tennis player, which was typical of Talbert."
Back to the University of Cincinnati in the fall of 1950, Trabert won the NCAA singles championship in '51 and played guard on the Bearcat basketball team that played in the NIT at Madison Square Garden that year. The NIT was much bigger than the NCAAs at the time.
Trabert also beat Talbert for the first time in singles in '51 at the Cincinnati International Championships, now known as the Cincinnati Masters.
"I had played him four or five times before that," Trabert said. "He was better than I was, he was ten years older than me. But he was my idol and it was hard for me to want to beat him. I have a picture of us shaking hands after the first time I beat him and he has the most beautiful smile on his face."
Drafted by the Navy in 1951 during the Korean War, Trabert reported to Bainbridge, Maryland for boot camp that fall.
"They joked I'd be playing tennis with Admirals soon. I got orders to Norfolk next and found myself hanging off the deck of the aircraft carrier Coral Sea, chipping paint off the sides of the ship."
Looking around for something different as a Seaman Apprentice, Trabert was offered a job as a quartermaster aboard the ship.
"I figured anything was better than what I was doing," he explained. "I was on the bridge a lot as the 'Captains Talker' on air defense." He was on the carrier for sixteen months, including a six month cruise in the Mediterranean.
The Coral Sea was deployed to "the Med" in the spring of '52, so there was a chance Trabert could get some leave and possibly play at the French and Wimbledon once again.
"I got liberty and went to play in the French championships. I won a couple of matches but lost to Felicisimo Ampon from the Philippines in the fourth round. He got everything back."
Another trip to play at Wimbledon the next month on liberty seemed to be in the plans but it wasn't to be.
"The Navy department sent a message to the captain that it was OK for me to go play at Wimbledon, but the captain turned me down," he explained. "I was on the bridge and the captain said, 'I guess you're mad at me.' He had already made admiral. 'I said, 'No sir,'. Who was I, as just a seaman, to be mad at an admiral? 'Nothing to do with you,' he told me. 'I just wanted to let them know who was in charge of this ship.'"
The next year Trabert was out of the Navy and playing at Wimbledon. The admiral contacted him to see if he could get him tickets to the tournament.
"I wish I could, but I didn't have any juice. I couldn't help him."
Trabert played singles in three of the four Grand Slam tournaments that year, only missing the Australian Championships. During his stint in the Navy he missed nine of the twelve major tournaments he was eligible for, but finished 1953 with a victory at the U.S. Championships (now the U.S. Open) at Forest Hills. In the next two years he won four of the eight majors, was a semi-finalist in two and a quarterfinalist in another.
His 1955 season is considered one of the best ever in tennis history. After a famous, and tight Davis Cup victory over the Australians in December of 1954, Trabert stayed to play in the Australian Open three weeks later. A loss to eventual champion Ken Rosewall in the semi-finals was one of the few losses he'd have that year.
"We went to Australia to get acclimated in November," he recalled. "It was fifty hours of flying to get there. Five, ten-hour legs. We stopped off in Hawaii for six hours or so to get gas. Vic (Seixas) and I went and played in an exhibition at the Royal Hawaiian hotel, in Honolulu, jumped in the ocean and got back on the plane."
Trabert went on to win Wimbledon, the French and the U.S. Championships in 1955 among his 18 tournament victories.
He helped Rene Lacoste develop what became the 'T2000' racket and promote the now-famous Lacoste "alligator" logo shirts. In Paris for three years to develop professional tennis in Europe, Asia and Africa, Trabert signed Rod Laver to his first professional contract in 1963. He served as president of the International Tennis Hall of Fame and had more than a three decades long career as a television commentator for CBS here in the U.S. and Channel 9 in Australia.
While equipment has changed the game in the last seventy years, Trabert says the top players of any era could compete against each other. The difference is the depth of the draw.
"Now guys have to be ready to play in the first round," he said. "When I was playing if you didn't go to sleep it was hard to lose before the quarterfinals. Now you have to be ready on the first Monday. Players are taller and bigger, better servers. Players come from all over these days. It's about as international a sport as you can get. It's great for the game."


Fashion photographer Lindsay Adler has risen to the top of her industry as both a photographer and educator. Based in New York City, her fashion editorials have appeared in numerous fashion and photography publications including Marie Claire, InStyle, Noise Magazine, Zink Magazine, Elle, Rangefinder, Professional Photographer and dozens more. As a photographic educator, she is one of the most sought after speakers internationally, teaching on the industry's largest platforms and most prestigious events, having been named one of the top 10 best fashion photographers in the world. Lindsay has worked with some of the top brands in the photographic and related industries such as Canon, Adobe and Microsoft.
A clean, bold, and graphic style has become the hallmark of her work, whether shooting advertising campaigns, designer look books, jewelry, hair campaigns, fashion editorials, or professional athletes. Lindsay is renowned for her creativity and collaborating with designers and stylists in order to create fresh looks.


What better place to start then the original source – the Black Pudding Producers. A great opportunity to go direct to those who know most about the product. Use the filters or advance search to find your Regional favourite, a local producer or a source of fresh blood black pudding.
If you are a Producer, upgrade your listing to share why you believe you produce the best black pudding!
Annascaul Black Pudding Co. offers home made black & white pudding, traditional sausages and dry cured rashers.
Appleton's is a traditional pork butcher. We've been selling quality meat from our shop in the heart of Ripon since 1867.
We take great pride in the quality of the meat we butcher, which will always be the basis for our business.
Sarah Pettegree is Bray's Cottage Pork Pie Queen. Supplying the finest independent retailers & pubs UK wide. Available online.
The only company making Bury Black Puddings in the town of Bury offers you the very highest quality, award winning Black Puddings!
Selfridges specialist butcher. Royal Warrant holders providing a fine selection of meat and game.
We sell Michelin star quality meat, fish and deli and we pride ourselves in being the best online family butcher in Scotland.
Chadwicks Quality Family Run Organic Butchers. Maker of Hobbit Pies for 'Lord Of The Rings' Premier!
Made to the original 70 year old recipe, our award winning Stornoway Black Puddings are legendary.
Purveyors of Blackpudding, made to a secret recipe since the 1880s.
A family business which has been manufacturing high quality and great tasting sausage and black pudding for nearly 20 years.
We're the leading independent foodservice brand. Let's bring our products to life – be a part of our Inspiration Population.
Purveyors of the finest traditional Irish Black and White Pudding, made to a secret family recipe.
Traditional butchers of the finest quality produce.
Using the finest ingredients, we make an old style black and white pudding, handed down through the years in the De Roiste family.


text:
Online preregistration
Urban Studies Institute
Hannah Arendt Institute - Urban Studies Institute
UAntwerp
Hannah Arendt Institute
Hannah Arendt Institute focuses on diversity, urbanity and citizenship
On Wednesday 27 May 2020, the Vrije Universiteit Brussel (VUB) and the University of Antwerp have launched the Hannah Arendt Institute for Diversity, Urbanity and Citizenship. The institute is an initiative of the two universities and is supported by the city of Mechelen and Bart Somers, Flemish Minister for Home Affairs, Public Governance, Civic Integration and Equal Opportunities.
In December, the VUB and the UAntwerp jointly launched the knowledge platform Diversity, Urbanity and Citizenship, a unique partnership between the two universities. In addition to this scientific platform will now have an educational centre that focuses on translating academic insights to professional target groups and the general public. These two pillars form the Hannah Arendt Institute for Diversity, Urbanity and Citizenship.
The new Hannah Arendt Institute will address social challenges such as migration, diversity and socioeconomic and political deprivation, as they are most strongly expressed in an urban context, and the broad societal impact of globalisation, ecological disruption, digitisation and the pressure on human rights and democratic values. In that context, the centre wants to make a fundamental contribution to the social debate and living together in the city. Christophe Busch will take on the daily management of the institute. As its home base, the institute will have a side wing of the Mechelen city hall, which will also contain space for workshops.
USI chair Stijn Oosterlynck (right) is one of the initiators of the Hannah Arendt Institute, here holding the banner with Pieter Ballon (director imec-SMIT, VUB).
Who was Hannah Arendt?
The German-Jewish political philosopher Hannah Arendt (1906-1975) was one of the most important political thinkers of the twentieth century. She emphasized the importance of connectedness, conversation and plurality for a democracy. Her studies on the nature of power, the public sphere and the political, direct democracy, freedom, civil disobedience, violence, authority, and totalitarianism still contain particularly relevant lessons for today's challenges.
The University of Antwerp website uses cookies and similar technologies to ensure the basic functionality of the site and for statistical and optimisation purposes. It also uses cookies to display content such as YouTube videos and for marketing purposes. This last category consists of tracking cookies: these make it possible for your online behaviour to be tracked. You consent to this by clicking on Accept. Also read our Privacy statement
Manage your personal preferences
Cookie prreferences
They are necessary for the proper functioning of the website and can not be switched off. They are generally defined only as a result of actions you carry out and that enable you to make use of the services we offer (definition of your privacy preferences, connection, filling in forms, etc.). You can configure your browser so that these cookies are blocked or you are warned about these cookies, but then some parts of the site will no longer work.
They allow us to measure surfing patterns and determine which sections of our website were visited, and to measure habits and improve our services. They show us which pages are the most and the least popular and how visitors move around the site. If you do not allow these cookies, we do not know when you visited our site and we can not check their habits.
They enable improved functionality and customization of the site. They can be activated by us or by third parties whose services we have added to our pages. If you do not allow these cookies, some or all of the services may not work properly.
They can be set up by our advertising partners on our site to build a profile of your interests and show you relevant ads on other sites. They do not store personal information directly, but are based on the unique identification of your browser and device. If you do not allow these cookies, you will see less targeted ads.


text:
Q: Can I get the google geocoder to tell me when a location can't be found, instead of just picking the center of a city? I'm working with the google maps API, and testing a section which geocodes all the addresses in a database which do not have coordinates, then generates a report of which were successful, and which failed.
If any of the locations couldn't be found, it should put a form at the top of the page to allow the user to either modify the address or specify coordinates.
I know that the geocoder can call back my code with the G_GEO_UNKNOWN_ADDRESS status code, but it never uses it. I put fake addresses into the database and got this:
The coordinates of 1337 Rawrmander Rd.
Were successfully set to:
    latitude:   49.170171
    longtitude: -123.136579

The coordinates of 5678 Imaginary Lane
Were successfully set to:
    latitude:   49.170171
    longtitude: -123.136579
It just gave back the center of the city (and continues to do so no matter what the address is) instead of calling it an unknown address.
This is how I call the geocoder:
//in a for loop
setTimeout( "geocoder.getLocations(\"" +
    //address to send to google
    values[i]['house'] + " " +
    values[i]['street'] + " " +
    values[i]['streetType'] + ", " +
    "Richmond, BC, Canada\", " +
    //google returns to this function asynchronously
    "function(reply) { handleReply(processItem++, reply); } );"
    ,
    callTimeout);
callTimeout += 220;

I tried changing "Richmond, BC, Canada" to just "Richmond" in hopes that it would not know which city to pick and then return an unknown address, but instead it just picked the middle of Richmond, VA.
Does anyone know a way to "hint" at google to give me an unknown address instead of just picking anything?

A: I think what you are looking for is the Accuracy Value.  When you run a geocode, look at the accuracy value that comes back.  Or are we talking about different things?
For example:
Using the address you provide (5678 Imaginary Lane) I get the following:
200,4,44.9774820,-93.2643510

versus using a real address (701 1st Ave N) the result is:
200,8,44.9782860,-93.2760830

Note the second value there (4=TownLevel vs 8=AddressLevel).  Looking for Imaginary,MN I'm able to get the Accuracy down to 2=RegionLevel.  I was able to get a 602=G_GEO_UNKNOWN_ADDRESS response (accuracy of 0=Unknown) when I looked for Imaginary,XX.
602,0,0,0

A: Check the coordinates of the address you have on hand, eg: "1234 Imaginary St, Vancouver" and then check the coordinates of just the city: "Vancouver". If they're the same, there's a problem. Unless they are living right in the very centre of town, that'd be a good indication the address doesn't work.

A: The Geocoder will return the status code of 602 when it really can't find a location.  More often than not though - it will just take a shot in the dark and give you a result.  The accuracy value will not tell you that you just got a shot in the dark.  It only tells you whether your result is a country, city, address, etc.  Unfortunately, there doesn't seem to be a way around this.

redpajama_set_name: RedPajamaStackExchange


Set of 2 double prong robe hooks.
Recently featured in a DIY mug rack on the popular DIY blog My Repurposed Life as seen in the additional images.
The additional images are from the popular DIY and repurposing blog My Repurposed Life.
And the wonderful DIY rustic hook rack photos comes from The Boondocks Blog!


This is a project to make new types of musical instruments by collaborating with other various fields: sound artists, sculptors, engineers, makers, and architects, etc. As mentioned on the project's name, the aim of this project was to search the boundaries of musical instruments. As a participant within various fields, I started this project asking how we can interpret musical instruments by our own way, depending on our own fields.
Not all space is an architecture. In order to be termed as an architecture, space needs to be organized and designed with a specific purpose. Likewise, it can't be said that all sound-producing devices are musical instruments. Not all sound is termed as a music. Therefore, with delicate observations and interests on ignored sound, I could finally change the sound into a music which contains the purpose of performer.
The Ball is a typological modification of a reed organ that the sound is produced from vibration of its thin metal plates called 'Reeds'. Comparing with a pair of bellows in a reed organ presses air through the reeds, each air pocket produces air movement through each reed in 12 modules of the Ball.
The ball-shaped instrument made of 3D printed modules produces sound in reaction to light with assistance of a MCU(Micro Controller Unit, Arduino).
The project ended with joint performance of 9 other instruments created. A documentary, which includes all the construction process, final performance and Q&A with experts from each fields, was filmed by renown educational channel of Korea, EBS, and broadcasted on September 2014.