Baseline monitoring using aircraft laser ranging

The use of aircraft laser ranging for the determination of baselines between ground based retroreflectors was investigated via simulations and with tests at Wallops Flight Center using the Airborne Oceanographic Lidar (AOL) on the Wallops C-54 aircraft ranging to a reflector array deployed around one of the Wallops runways. The aircraft altitude and reflector spacing were chosen on the basis of scaled down modeling of spacecraft tracking from 1000 km of reflectors separated by some 52 km, or of high altitude (10 km) aircraft tracking of reflectors separated by some 500 m. Aircraft altitudes flown for different passes across the runway reflector array varied from 800 m to 1350 m, with 32 reflectors deployed over an approximtely 300 m x 500 m ground pattern. The AOL transmitted 400 pulses/sec with a scan rate of 5/sec in a near circular pattern, so that the majority of the pulses were reflected by the runway surface or its environs rather than by retroreflectors. The return pulse characteristics clearly showed the high reflectivity of portions of the runway, with several returns indistinguishable in amplitude from reflector returns. For each pass across the reflector field, typically six to ten reflector hits were identified, consistent with that predicted by simulations and the observed transmitted elliptical pulse size

3-D multiobservable probabilistic inversion for the compositional and thermal structure of the lithosphere and upper mantle: III. Thermochemical tomography in the Western-Central U.S.

Acknowledgments We are indebted to F. Darbyshire and J. von Hunen for useful comments on earlier versions of this work. This manuscript benefited from thorough and constructive reviews by W. Levandowski and an anonymous reviewer. We also thank J. Connolly, M. Sambridge, B. Kennett, S. Lebedev, B. Shan, U. Faul, and M. Qashqai for insightful discussions about, and contributions to, some of the concepts presented in this paper. The work of J.C.A. has been supported by two Australian Research Council Discovery grants (DP120102372 and DP110104145). Seismic data are from the IRIS DMS. D.L.S. acknowledges support from NSF grant EAR-135866. This is contribution 848 from the ARC Centre of Excellence for Core to Crust Fluid Systems (http://www.ccfs.mq.edu.au) and 1106 in the GEMOC Key Centre (http://www.gemoc.mq.edu.au).Peer reviewedPublisher PD

A new free-surface stabilization algorithm for geodynamical modelling:Theory and numerical tests

The surface of the solid Earth is effectively stress free in its subaerial portions, and hydrostatic beneath the oceans. Unfortunately, this type of boundary condition is difficult to treat computationally, and for computational convenience, numerical models have often used simpler approximations that do not involve a normal stress-loaded, shear-stress free top surface that is free to move. Viscous flow models with a computational free surface typically confront stability problems when the time step is bigger than the viscous relaxation time. The small time step required for stability (<2. Kyr) makes this type of model computationally intensive, so there remains a need to develop strategies that mitigate the stability problem by making larger (at least ~10 Kyr) time steps stable and accurate. Here we present a new free-surface stabilization algorithm for finite element codes which solves the stability problem by adding to the Stokes formulation an intrinsic penalization term equivalent to a portion of the future load at the surface nodes. Our algorithm is straightforward to implement and can be used with both Eulerian or Lagrangian grids. It includes α and β parameters to respectively control both the vertical and the horizontal slope-dependent penalization terms, and uses Uzawa-like iterations to solve the resulting system at a cost comparable to a non-stress free surface formulation. Four tests were carried out in order to study the accuracy and the stability of the algorithm: (1) a decaying first-order sinusoidal topography test, (2) a decaying high-order sinusoidal topography test, (3) a Rayleigh-Taylor instability test, and (4) a steep-slope test. For these tests, we investigate which α and β parameters give the best results in terms of both accuracy and stability. We also compare the accuracy and the stability of our algorithm with a similar implicit approach recently developed by Kaus et al. (2010). We find that our algorithm is slightly more accurate and stable for steep slopes, and also conclude that, for longer time steps, the optimal α controlling factor for both approaches is ~2/3, instead of the 1/2 Crank-Nicolson parameter inferred from a linearized accuracy analysis. This more-implicit value coincides with the velocity factor for a Galerkin time discretization applied to our penalization term using linear shape functions in time

The morphology, evolution and seismic visibility of partial melt at the core-mantle boundary: Implications for ULVZs

SUMMARY Seismic observations indicate that the lowermost mantle above the core–mantle boundary (CMB) is strongly heterogeneous. Body waves reveal a variety of ultra-low velocity zones (ULVZs), which extend not more than 100 km above the CMB and have shear velocity reductions of up to 30 per cent. While the nature and origin of these ULVZs remain uncertain, some have suggested they are evidence of partial melting at the base of mantle plumes. Here we use coupled geodynamic/thermodynamic modelling to explore the hypothesis that present-day deep mantle melting creates ULVZs and introduces compositional heterogeneity in the mantle. Our models explore the generation and migration of melt in a deforming and compacting host rock at the base of a plume in the lowermost mantle. We test whether the balance of gravitational and viscous forces can generate partially molten zones that are consistent with the seismic observations. We find that for a wide range of plausible melt densities, permeabilities and viscosities, lower mantle melt is too dense to be stirred into convective flow and instead sinks down to form a completely molten layer, which is inconsistent with observations of ULVZs. Only if melt is less dense or at most ca. 1 per cent more dense than the solid, or if melt pockets are trapped within the solid, can melt remain suspended in the partial melt zone. In these cases, seismic velocities would be reduced in a cone at the base of the plume. Generally, we find partial melt alone does not explain the observed ULVZ morphologies and solid-state compositional variation is required to explain the anomalies. Our findings provide a framework for testing whether seismically observed ULVZ shapes are consistent with a partial melt origin, which is an important step towards constraining the nature of the heterogeneities in the lowermost mantle and their influence on the thermal, compositional and dynamic evolution of the Earth.ER

Unraveling the physics of the Yellowstone magmatic system using geodynamic simulations

The Yellowstone magmatic system is one of the largest magmatic systems on Earth, and thus an ideal location to study magmatic processes. Whereas previous seismic tomography results could only image a shallow magma reservoir, a recent study using more seismometers showed that a second and massive partially molten mush reservoir exists above the Moho (Huang et al., 2015). To understand the measurable surface response of this system to visco-elasto-plastic deformation, it is thus important to take the whole system from the mantle plume up to the shallow magma reservoirs into account. Here, we employ lithospheric-scale 3D visco-elasto-plastic geodynamic models to test the influence of parameters such as the connectivity of the reservoirs and rheology of the lithosphere on the dynamics of the system. A gravity inversion is used to constrain the effective density of the magma reservoirs, and an adjoint modeling approach reveals the key model parameters affecting the surface velocity. Model results show that a combination of connected reservoirs with plastic rheology can explain the recorded slow vertical surface uplift rates of around 1.2 cm/year, as representing a long term background signal. A geodynamic inversion to fit the model to observed GPS surface velocities reveals that the magnitude of surface uplift varies strongly with the viscosity difference between the reservoirs and the crust. Even though stress directions have not been used as inversion parameters, modeled stress orientations are consistent with observations. However, phases of larger uplift velocities can also result from magma reservoir inflation which is a short term effect. We consider two approaches: (1) overpressure in the magma reservoir in the asthenosphere and (2) inflation of the uppermost reservoir prescribed by an internal kinematic boundary condition. We demonstrate that the asthenosphere inflation has a smaller effect on the surface velocities in comparison with the uppermost reservoir inflation. We show that the pure buoyant uplift of magma bodies in combination with magma reservoir inflation can explain (varying) observed uplift rates at the example of the Yellowstone volcanic system

GEM-Selektor geochemical modeling package: revised algorithm and GEMS3K numerical kernel for coupled simulation codes

Reactive mass transport (RMT) simulation is a powerful numerical tool to advance our understanding of complex geochemical processes and their feedbacks in relevant subsurface systems. Thermodynamic equilibrium defines the baseline for solubility, chemical kinetics, and RMT in general. Efficient RMT simulations can be based on the operator-splitting approach, where the solver of chemical equilibria is called by the mass transport part for each control volume whose composition, temperature, or pressure has changed. Modeling of complex natural systems requires consideration of multiphase-multicomponent geochemical models that include nonideal solutions (aqueous electrolytes, fluids, gases, solid solutions, and melts). Direct Gibbs energy minimization (GEM) methods have numerous advantages for the realistic geochemical modeling of such fluid-rock systems. Substantial improvements and extensions to the revised GEM interior point method algorithm based on Karpov's convex programming approach are described, as implemented in the GEMS3K C/C+ + code, which is also the numerical kernel of GEM-Selektor v.3 package ( http://gems.web.psi.ch ). GEMS3K is presented in the context of the essential criteria of chemical plausibility, robustness of results, mass balance accuracy, numerical stability, speed, and portability to high-performance computing systems. The stand-alone GEMS3K code can treat very complex chemical systems with many nonideal solution phases accurately. It is fast, delivering chemically plausible and accurate results with the same or better mass balance precision as that of conventional speciation codes. GEMS3K is already used in several coupled RMT codes (e.g., OpenGeoSys-GEMS) capable of high-performance computin

Entwicklung und Anwendung von Hochleistungs-Software für Mantelkonvektionssimulationen

The Earth mantle convects on a global scale, coupling the stress field at every point to every other location at an instant. This way, any change in the buoyancy field has an immediate impact on the convection patterns worldwide. At the same time, mantle convection couples to processes at scales of a few kilometers or even a few hundred meters. Dynamic topography and the geoid are examples of such small-scale expressions of mantle convection. Also, the depth of phase transitions varies locally, with strong influences on the buoyancy, and thus the global stress field. In order to understand these processes dynamically it is essential to resolve the whole mantle at very high numerical resolutions. At the same time, geodynamicists are trying to answer new questions with their models, for example about the rheology of the mantle, which is most likely highly nonlinear. Also, due to the extremely long timescales we cannot observe past mantle states, which calls for simulations backwards in time. All these issues lead to an extreme demand in computing power. To cater to those needs, the physical models of the mantle have to be matched with efficient solvers and fast algorithms, such that we can efficiently exploit the enormous computing power of current and future high performance systems. Here, we first give an extensive overview over the physical models and introduce some numerical concepts to solve the equations. We present a new two-dimensional software as a testbed and elaborate on the implications of realistic mineralogic models for efficient mantle convection simulations. We find that phase transitions present a major challenge and suggest some procedures to incorporate them into mantle convection modeling. Then we give an introduction to the high-performance mantle convection prototype HHG, a multigrid-based software framework that scales to some of the fastest computers currently available. We adapt this framework to a spherical geometry and present first application examples to answer geodynamic questions. In particular, we show that a very thin and very weak asthenosphere is dynamically plausible and consistent with direct and indirect geological observations.Englische Übersetzung des Titels: Development and application of high performance software for mantle convection modelin