Geophysical Journal International

SummaryIn Distributed Acoustic Sensing (DAS), a fiber-optic cable is used as a distributed seismic sensor, with channels representing successive short sections of the fiber, spaced at defined intervals along the 1D fiber axis. Typically, the positions of these channels are assumed to be a line projection along the cable's position. In reality, a fiber-optic cable contains many fibers that are not perfectly straight and are thus longer than the cable itself. Consequently, the real channel positions may not correspond to a simple interpolation along the cable axis. Moreover, the precise cable coordinates are often sensitive information and may not be provided to the end user who uses the cable for sensing applications. On land, a tap test is usually carried out before the start of a DAS acquisition to determine the exact channel locations. DAS with marine horizontal cables has recently been used for various offshore applications, including seismic imaging. To avoid errors in the seismic image, a precise receiver location is required. In this paper, we propose a traveltime-based inversion workflow to determine a more accurate channel position on the seafloor. Moreover, we show that we can resolve an unknown time shift between the acquisition and the recording system, in addition to the fiber position.

SummaryTransform faults are one of the major tectonic plate boundaries offsetting the global mid-oceanic ridge system. The topographic features within these transform faults provide crucial evidence for tectono-magmatic processes and crustal accretion in transform fault zones. These interesting features include median ridges, which are major bathymetric anomalies found within both slow-slipping and fast-slipping transform faults, often associated with exposures of ultramafic rocks on the seafloor. To explain the origin of median ridges, previous studies have invoked multiple processes such as serpentinite diapirism, thermal uplift at ridge-transform intersections, or transpressive uplift induced by global plate reorganization, without any knowledge of the seismic structure. Here, we present results from 2D travel time tomography of downward-continued multi-channel seismic data along and across an ∼80 km long median ridge that lies within the eastern end of the slow-slipping (∼3.4 cm/yr) Chain transform fault in the equatorial Atlantic Ocean. The data were acquired during the 2018 ILAB-SPARC survey using a 6-km long streamer. Our high-resolution P-wave velocity model of the median ridge shows distinct high and low velocities ranging from 2.5 to 5 km/s within 500 m below the seafloor, on either side of the presently active strike-slip fault trace that cuts through the ridge. The low velocity on the eastern side of the ridge could be due to the presence of highly fractured basalt (with porosity in the range of 28 to 36 per cent) due to transform fault motion, whereas the high velocity on the western flank could be due to the presence of gabbro or highly serpentinised peridotite. The basaltic origin of the median ridge is supported by the observation of a seismic triplication event, which we call the T-event. The depth at which the T-event maps is shallow (200–500 m below seafloor) in high-velocity regions and deeper (600–1400 m) in low-velocity regions. We also find that the currently active strike-slip fault has been active since at least 0.26 Ma and has sliced the ridge. We image low-velocity pockets at the northern and southern limits of the median ridge that could represent the expression of the currently less active strike-slip faults.

SummaryThe northeastern Tibetan Plateau is bounded by the left-lateral Altyn Tagh and Haiyuan faults. How crustal motion along these fault systems transitions to crustal shortening and uplift is key for deciphering the geodynamic link between the escape tectonics and the growth of the Tibetan Plateau. Here, we use the PS-InSAR observations, combined with GNSS and leveling data, to obtain a high-resolution 3D model of the present-day crustal motion in the northeastern Tibetan Plateau. The resolved deformation field covers the entire northeastern Tibetan Plateau with a spatial resolution of approximately 0.01° × 0.01 °. Our analysis of slip rates and strain partitioning reveals that crustal motion along the Altyn Tagh fault gradually diminishes eastward and is absorbed by thrusting and uplift in the Qilianshan orogenic belt within the plateau. A similar tectonic transition occurs between the Haiyuan fault and the Liupanshan orogen on the eastern margin of the plateau. Some of the eastward crustal motion is accommodated by the younger Xiangshan-Tianjingshan fault system to the north of the Haiyuan fault, indicating the ongoing northward expansion of the Tibetan Plateau. Our results align with geological evidence of crustal deformation in the past few million years, highlighting the continuing tectonic transition from eastward crustal motion along the left-lateral strike-slip faults to the growth of the Tibetan Plateau.

SummaryThe continuous NE-SW compression due to the Indo-Asian collision creates active and complex deformation in the northeastern (NE) Tibetan plateau. How the lithosphere of the NE Tibetan plateau deforms both vertically and laterally in response to the ongoing collision is still a question. Further investigations with refined lithospheric structure are required. Here we present a high-resolution radially anisotropic model of the lithosphere beneath the NE Tibetan plateau and surrounding areas that was constrained by the joint analysis of Rayleigh and Love wave dispersion at periods from 6 to 100 s using the methods of ambient noise cross-correlation for short periods and earthquake two-plane-wave for long periods. Results show relatively small regions of significant slow shear wave velocity and positive radial anisotropy (Vsh > Vsv) in the middle crust beneath the Qilian orogenic belt, suggesting the existence of partial melting and probably limited channel flow. Considering the variable lateral strength of shear wave velocity and radial anisotropy in the middle crust with large parts mechanically strong enough to pass the strain, vertical coherent lithospheric deformation could still work in the Qilian orogenic belt. Extensive low shear wave velocity anomaly in the uppermost mantle extend from the Qilian orogenic belt northward to the Alxa block and eastward to the southwestern Ordos block, implying a hot and weak mantle lithosphere. The observed negative radial anisotropy (Vsh < Vsv) in such warm mantle lithosphere beneath the Qilian orogenic belt, Alxa block and the southwestern Ordos block is ascribed to vertical deformation fabrics arising from the convergence between Indian and Asian plates. These observations imply that the lithosphere of Qilian orogenic belt, Alxa block and southwestern Ordos block deform coherently and the NE Tibetan plateau is expanding towards Alxa block and the southwestern Ordos block.

SummaryWe develop a new method for estimating the autocorrelation function (ACF) from segmented data with the assumption of stochastic stationarity. The ACF of a signal is represented as the summation of the cross terms of sub-segments of arbitrary length. To successfully remove undesired transients in the data, this method introduces a correction for the amplitude bias associated with the removal of sub-segments, based on the comparison between the expected stationary signal and the measured signal. The method reconstructs and accesses later lag times, provides finer frequency resolution, obtains a better signal-to-noise ratio, which enables the extraction of detailed temporal or spectral structures from noisy data sets. As an application, we successfully retrieved a spectrum of the Earth’s seismic hum on the vertical component with fine frequency resolution and compared it to synthetic autocorrelation for spatially isotropic and homogeneous excitation by random shear traction at the ocean bottom and random pressure at the Earth’s surface. Although both models can explain the observed fundamental spheroidal modes, shear traction is better at explaining the observed overtones above 3 mHz. From 2 to 3 mHz, the pressure source also contributes to the excitation of the overtones, and the shear traction becomes dominant again below 2 mHz. This new method is anticipated to be effective in extracting valuable information from rare records within the context of extraterrestrial seismology.

SummaryNortheast China, with its complicated regional tectonic evolution, situated within the eastern Central Asian Orogenic Belt, is a key region for understanding lithospheric deformation and mantle dynamics. However, the ongoing debate surrounding its lithospheric structure and evolutionary processes remains, largely attributed to data limitations and methodological constraints. In this study, we integrate topography, geoid height, surface heat flow, and Rayleigh wave phase velocity dispersion curves to conduct a detailed imaging of the lithospheric thermal and compositional structure in Northeast China. We find a significant east-west gradient in lithospheric thickness, ranging from approximately 60 km in the east to 140 km in the west, and a compositional transition in the lithospheric mantle from fertile peridotite in the east to refractory peridotite in the west. By integrating analyses of upper mantle anisotropy and the spatiotemporal distribution of Mesozoic basalts, we argue that the lithospheric delamination and mantle upwelling may have combined to cause the lithospheric thinning in the region. This study highlights the significance of joint inversion of multiple datasets and integrated multidisciplinary analysis.

SummaryTraditional methods for measuring geopotential difference using optical fiber frequency transfer or satellite-based time-frequency transfer, based on general relativity, require the use of two clocks and the calibration of these clocks. Here we present a simplified clock transportation experiment using a single hydrogen clock to measure the geopotential difference between two time-frequency stations, separated by 129 km with a height difference of 1,245 m, by GNSS precise point positioning time-frequency transfer. Taking the reference clock of International GNSS Service (IGS) time as a ‘bridge’, we extract the gravity frequency shift between the two stations by comparing the fractional frequency differences between the hydrogen clock and the ‘bridge’ before and after clock transportation. The determined geopotential difference between the two stations is 12075.9$\pm $118.5 m²/s², which closely aligns with the value computed by the EIGEN-6C4 global gravity field model, with a difference of -78.7 m²/s². These results validate the feasibility of geopotential difference measurements with a single clock and highlight several advantages compared to the dual-clock method: elimination of inter-clock calibration, low operational complexity and equipment cost, high data utilization efficiency but similar precision of geopotential difference measurement. Furthermore, this method can be extended to other similar techniques to measure geopotential differences, provided that they enable users to connect to a stable time-frequency reference.

SummaryShale gas extraction could produce underground stress perturbation and local seismicity, which could put a threat human casualty. The Weiyuan area in the Sichuan province, China, underwent massive gas production and a significant increase of earthquake since 2015. In this study, we focus on human-induced subsurface hydrofracturing, calculate cumulative underground Coulomb-stress changes using a 3D numerical model, and probe the main cause of recent seismic activity in the Weiyuan area based on continuous regional stress/displacement loading. The simulation reveals the regional extent of positive Coulomb stress change with fractures matches the distribution of the moderate and micro seismicity in the past ten years. Background regional tectonic stress in the vicinity of the active fault likely resulted in earthquake preparation within and around the active faults; hydraulic fracturing changes mainly the displacement and stress pattern in the vicinity of the fracturing wells, and enhanced fracturing intensity (fracture volume-to-model volume ratio (θ), causes more obvious difference; faults may be locked prior to fracturing, and even small fracturing intensity may trigger the earthquakes near faults and fracturing wells; the seismic risk will be significantly increased near the two faults and fracturing wells in the next 50 years.

SummaryIn this study, we compare the usability of a simplified microtremor-based empirical method and a conventional microtremor method based on an inversion analysis of a subsurface velocity structure model for constructing a map of average S-wave velocity (AVS) values. In the simplified (empirical) method, the phase velocities of Rayleigh waves, which can be obtained by processing a microtremor array, at wavelengths of 13, 25, and 40 m are regarded as AVS values from the ground surface to depths of 10, 20, and 30 m (${\overline {Vs} }_{10},{\overline {Vs} }_{20},\ {\rm{and}}\ {\overline {Vs} }_{30}$), respectively. Microtremor array surveys were conducted at 173 observation points within a 15 km × 17 km area east of Aso caldera, Kyushu, Japan (target area). AVS values are obtained by applying the empirical method to the phase velocities obtained at each observation point. The AVS values at an observation point (located near the centre of the target area) with velocity logging data are verified by a comparison with those based on the velocity logging data (i.e. overestimations by 6 per cent at maximum). It is found that for the entire target area, the spatial distribution of the obtained AVS values is consistent with the geological distribution. The AVS values within areas of the Aso-3 ignimbrite are 30–40 per cent larger than those within areas of thick soil and tephra on the strongly consolidated Aso-4 ignimbrite. In addition, the AVS values of the Aso-3 deposits are more than 10 per cent larger than those of the Aso-4 deposits and about 10 per cent smaller than those of geological units older than the Aso-3 deposits. We also apply a conventional (i.e. inversion) method to the phase velocity data at each observation point to obtain a one-dimensional S-wave velocity (Vs) structure model from which we deduce AVS values. The deduced AVS values at the velocity logging point are underestimated by -8 per cent, with differences from the AVS values obtained using the empirical method reaching 13 per cent. The average systematic difference between the two methods is 15 per cent, as determined from a statistical analysis. None the less, a strong correlation is found between the methods, with an average correlation coefficient of 0.94, with no evidence showing that either method is more accurate. The empirical method can be used to construct an AVS map if overestimation is carefully considered. This analysis also reveals that the average maximum survey depths of the one-dimensional Vs structures based on the inversion method are only 23±10 m, making them often insufficient to map ${\overline {Vs} }_{20}$ and ${\overline {Vs} }_{30}$ (the ratios of the available to total numbers of data points are only 60 and 21 per cent, respectively). In contrast, the empirical method can determine ${\overline {Vs} }_{10},{\overline {Vs} }_{20},\ \ {\rm{and}}\ {\overline {Vs} }_{30}$ at more than 80 per cent of all sites. The construction of AVS maps using the empirical method is effective in terms of the simplicity and reliability of planning, observational efficiency, and simplicity of data processing, which support a practical and objective approach to seismic assessments.

SummaryIn this study we obtain 35 903 high-quality P-wave receiver functions from 1737 teleseismic events recorded at 120 dense broadband TanluArray temporary stations deployed in and around the Tanlu fault zone (TLFZ). After station azimuth and sediment correction are made, a detailed Moho depth distribution is obtained by CCP stacking. Our results show a sharp change in the Moho depth across the TLFZ from the west to east, which well corresponds to the surface geological structure. The deepest Moho (38.0 ∼ 40.0 km) occurs beneath the Dabie orogenic belt and the Sulu orogenic belt. The Moho beneath the Luxi uplift, Jiangnan orogenic belt and Jiaodong uplift is deeper (36.0 ∼ 37.0 km), whereas the Subei basin and the southern basin of the South Yellow Sea have a shallow Moho (28.0 ∼ 30.0 km). There is an obvious Moho uplift near Weifang, which corresponds to the Changle ancient volcano on the surface and may be a channel for upwelling of hot mantle material. The Moho is unclear under the fault zone near Tancheng, which is speculated to be a channel for upwelling of hot mantle material. It may be related to upwelling of hot and wet flows in the big mantle wedge above the subducted Pacific slab that is stagnant in the mantle transition zone beneath East Asia, which is a possible cause of the 1668 M8.5 Tancheng earthquake.

SummaryThe Limpopo transform margin offshore southern Mozambique results from the separation of Gondwana along the East Africa continental margin. Over the last three decades, more than thirty different reconstruction models have been proposed, sometimes contradicting each other. Here, we present results from the travel-time tomography of wide-angle seismic data acquired during the second China-Mozambique Joint Cruise, allowing the interpretation of the crustal structure and magmatism in the Limpopo Corridor and the Mozambique Basin. Using these results, we determine the extent of the Continent Ocean Transition and the location of the Continent Ocean Boundary on the southern Mozambique margin. The seismic profile is 442-km long, extending from the eastern part of the North Natal Valley in the west and crossing the Limpopo Corridor and the Mozambique Basin to the east. Based on the tomographic velocity model, we delineated three distinct domains from west to east along the profile: (1) a western transitional domain with anomalous or mixed crust, bounded by the Mozambique Fracture Zone to the east, where the crust gradually thins eastward from ∼14 km at distance 45 km to ∼10.8 km at distance 140 km; (2) a domain of thickened oceanic crust resulting from enhanced magmatism, where the crust thins eastward, from ∼10.8 km to ∼8.5 km over ∼100 km distance; and (3) an eastern domain of normal oceanic crust, where the average crustal thickness is ∼8 km. We suggest that (1) the western transitional domain roughly corresponds to the Limpopo Corridor and is of continental crustal origin but was affected and modified by strike-slip motion and magmatic activity, resulting in anomalous or mixed crust. The eastern Continent Ocean Boundary of the Limpopo Margin is close to the Mozambique Fracture Zone; (2) The thickened oceanic domain thins eastward, and the crustal velocity and thickness change dramatically compared to the oceanic domain. This domain seems to have strongly interacted and contaminated by the Limpopo Corridor during the opening of Mozambique Basin and seafloor spreading; (3) The eastern oceanic domain shows a relatively uniform oceanic crust of ∼8 km and high velocity up to 7.4 km/s in the lower crust, suggestive of a hotter mantle that produces more MgO-rich melts probably due to the influence of a thermal mantle anomaly.

SummaryThe downward continuation of gravity field can provide valuable information for 3-D gravity-field modeling, shallow-layer geological interpretation, source depth estimation, and so on. However, downward continuation is ill-posed, and traditional approaches often suffer from computational instability, poor noise resistance, and limited continuation depth, making it a longstanding challenge in gravity data processing. We present a new approach for fast, stable and large-depth downward continuation of gravity anomalies by using frequency-domain 3-D imaging. First, we utilize the frequency-domain 3-D imaging approach to invert the gravity anomalies at the original observational plane to quickly obtain the equivalent density model in the subsurface. Then, we apply the optimized strategy of frequency-domain 3-D forward calculation on the equivalent density model to rapidly obtain high-precision gravity anomalies at the downward-continuation plane. The synthetic data tests prove the effectiveness of our approach, and demonstrate that our approach enables fast, stable, robust noise resistance and large-depth downward continuation of large-scale gravity anomalies data, and has superior performance compared to the traditional regularized filtering approach and spatial-domain equivalent-source approach. The real data test of the free-air gravity anomalies data in the central South China Sea also verifies the fast, stable and reliable downward continuations of large depths by our approach. The 3-D gravity-field model built by our approach will provide significant support for the tectonic studies and resource exploration in this area.

SummaryModelling and inversion of controlled-source electromagnetic data requires elaborate numerical tools. The major challenge is the high computational cost of computing solutions to numerous forward problems (for the forward responses as well as the sensitivity matrix). Forward modelling is accomplished using either a direct or an iterative solver. Current modelling suites predominantly employ direct solution methods in the forward operator since multiple solutions are easily accessible using inexpensive and quick forward-backward substitution after an initial resource-demanding matrix factorisation step. Iterative techniques, on the other hand, require little resources for single forward solutions, and are yet very time consuming if solutions for many right-hand sides are to be computed. Evaluations of different solution techniques for modelling and inverse problems are only sparsely investigated. In light of this, we integrated an iterative solver as alternative in the forward and and inversion operators of the open-source software custEM and pyGIMLi. In particular, we implemented a two-level iterative scheme where the outer solver employs a generalised conjugate residual algorithm preconditioned with a highly efficient block-based preconditioner for square blocks. The inner-level solver is either of the same type as the outer solver, but preconditioned with the auxiliary-space Maxwell preconditioner, or may alternatively be a direct solver. In this paper, we evaluate the described iterative forward operator for forward modelling tasks for the Marlim R3D model for a single as well as numerous right-hand side vectors and compare the performance to the direct solver MUMPS. We further investigate the solver’s applicability on small and medium-sized computing platforms. We then examine the iterative solver for inversions of synthetic land-based and semi-airborne data in terms of computational requirements. Our results demonstrate that forward modelling tasks are best performed using an iterative approach for single source problems. Moreover, simulations of large and complex problems are accessible on even on small computing platforms such as laptops in very reasonable time. For inversions, the iterative forward operator, in particular the mixed iterative-direct-based one, performs equally well in terms of time as the direct one while reducing the memory demands for the computations of the forward responses and the data sensitivities.

SummaryIn a region of complex geology, we examine the influence of spatial resolution of conductivity models on Geomagnetically Induced Currents (GICs) estimations. We focus on the southern region of Portugal mainland, for which magnetotelluric (MT) sounding measurements have been obtained with lower noise from human activity. Using two conductivity models inverted from sets of MT soundings with different sampling distance, we look for an interpretation of the differences in GIC estimations at substation grounding resistances. We make use of two different proxies, the Local Effective Field (LEF) and the Regional Electromotive Source (RES), built from the electric induced field at each substation site and the sum of electromotive forces along all transmission lines connected to that substation, respectively. We compare different time signals associated to GICs using a parameter that combines Pearson correlation and linear regression slope, the Correlation Regression Coefficient (CRC). Our main conclusion is that spatially detailed information on lateral heterogeneities of the conductivity associated to complex geology is crucial for a rigorous assessment of GIC hazard, leading to relative differences in GIC standard deviation and in GIC peak values that can amount to more than 100% in certain cases. Additionally, using LEF and RES, we emphasise the non-locality of GIC drivers and bring new input concerning the choice of proxies used to monitor and forecast this kind of hazard.

SummaryThe Scandinavian Peninsula and its vicinity comprise highly tectonically diverse blocks, including the Baltic Shield, the continental margin, and the North Sea Basin. The crustal rheology is a critical constraint to understanding the tectonic evolution in this region. Based on 19 416 Lg waveforms from 233 earthquakes and 560 broadband digital stations, using an inversion method combining both single- and two-station ray paths, we constructed a broadband (0.05 and 10.0 Hz) Lg wave attenuation model in the study region, with the resolution approaches to 110 km (∼1°) or higher in areas with dense ray path coverages. The QLg distributions correlate well with regional geological features. The Baltic Shield exhibits the highest QLg, consistent with its thick Precambrian crust and high rheological rigidity developed through Archean Svecofennian orogeny. In contrast, passive margins with crustal thinning, magmatic modification, and thick sedimentary sequences exhibit strong attenuation, reflecting a reduction in rheological strength resulting from interactions with mantle plumes and extensional tectonics. The North Sea Basin exhibits the lowest QLg values and the presence of hydrocarbon-bearing sediments. The extremely high QLg distribution reveals the ancient cratonic core of the Baltic Shield, particularly in areas where the surface rock dating sample cannot be collected due to seawater coverage.

SummaryWe analysed infrasound waves associated with the Gyeongju earthquake (ML 5.8) that occurred on September 12, 2016, in the southeastern Korean Peninsula. For infrasound wave detection, the Progressive Multi-channel Correlation method was applied to the infrasound dataset recorded at 7 arrays operating in South Korea at epicentral distances ranging from 178 to 472 km. Based on the back-projection method constrained by array-dependent celerity and azimuth deviation models, the source regions were identified in both the epicentral and nonepicentral regions. Remarkably, the nonepicentral secondary sources of this earthquake were located in regions with shallow water depths: i) the western coastal area in the Yellow Sea and ii) the shallow ocean basin and bank in the East Sea. The location results obtained from the earthquake could be corroborated through its foreshock (ML 5.1), yielding location results consistent with those of the mainshock. The generation of infrasound waves over shallow water depths was fortuitously validated by direct recordings of dominant single-frequency (∼0.3 Hz) infrasound waves at close range via temporary sensors near the ocean basin and bank. We interpreted that low-frequency infrasound signals could be generated from interactions among the ocean floor, shallow seawater, and atmosphere. We performed numerical simulations of seismoacoustic fields to predict ground motions on the seafloor and acoustic transmission efficiency between the water and air interface. The simulations quantified the energy transfer through different media and clarified our observational results. We found that because this solid Earth‒water‒atmosphere coupled air wave has a relatively low frequency (∼0.3 Hz), it can survive propagation over long distances compared with high-frequency infrasound waves generated in inland and mountain regions. In this study, we extend our understanding of water‒atmosphere coupling and the monitoring framework for earthquake-associated nonepicentral infrasound waves, encompassing not only inland ground shaking but also shallow sea regions located far from the epicentre.

SUMMARYIn recent years, machine learning (ML) techniques have emerged as a powerful tool in seismology, enabling the detection of small-magnitude seismic events that typically go unnoticed by traditional methods. Here, we apply ML-based methods to improve the characterization of normal fault systems and aftershock activity in the Central Apennines, using data from the 2009 L'Aquila seismic sequence. By processing data from both permanent and temporary seismic stations, we identified approximately 191 000 events—with a local magnitude range of -1.83 and 5.96, recorded during January-December 2009—nearly ten times more than the standard catalog maintained by INGV. These events were relocated using a combination of absolute and relative location techniques, resulting in a high-resolution catalog of 148 000 earthquakes. This catalog is distinguished by an increased number of S-wave pickings, which significantly reduces localization errors and enhances the accuracy of fault geometry reconstruction. Compared to an existing semi-automatic catalog, we observe a full recovery of seismic events, and a significant improvement of new events identified and well-located by the ML-approach, with a marked increase in the quality and quantity of P- and S-wave arrivals. The refined seismic catalog not only provides a more detailed and accurate definition of the fault architecture but also offers new insights into the distribution of aftershocks, unrolling the complex pattern of faulting that normally remains masked during standard analyses. This work highlights the potential of ML methods in advancing our understanding of complex fault systems and seismic sequences.

AbstractProgressively denser mapping of ocean-floor magnetization has led to detailed reconstructions of past plate motions in the Cenozoic. These reconstructions often reveal rapid kinematic changes that provide crucial information for identifying geodynamic mechanisms that may have caused them, and for quantifying force budgets upon plates. In parallel to these advances, the notion of thin, low-viscosity asthenosphere beneath tectonic plates that facilitates their motions has emerged and consolidated. This weak, mobile layer promotes the formation of the pressure-driven Poiseuille flow that, in turn, generates basal shearing upon plates. In addition, it can be linked to dynamic topography variations due to pulsing plume activity. In this study, we use publicly available finite-rotation compilations of the North American plate (NA) to investigate its kinematic history since Oligocene time. After removing data that are possibly impacted by significant noise, we find that NA experienced a westward speedup near 27 Ma. Next, we explore the role that asthenospheric Poiseuille-type flow caused by increased Canary plume activity may have had in generating this kinematic change. Such plume activity is inferred from the combination of anomalously shallow residual bathymetry and records of past ocean-floor magmatism offshore northwestern Africa. We compare estimates of torque variation upon NA that are (i) required to explain the reconstructed kinematic change, and (ii) predicted by the Poiseuille-type flow associated with the Canary plume activity. Our results indicate that these two torque-variations estimates are in agreement with each other, both in terms of direction and magnitude. This inference suggests that the increased Canary plume activity is a geodynamically-plausible process to explain the Oligocene plate-motion change of NA.

SummarySeismic traveltime tomography represents a popular and useful tool for unravelling the structure of the subsurface across the scales. In this work we address the case where the forward model is represented by the eikonal equation and derive a formalism to solve the inverse problem where gradients are calculated efficiently using the discrete adjoint state method. Our approach provides gradients with respect to both velocity structure and source locations, allowing us to perform a consistent joint inversion. The forward problem is solved using a second-order fast-marching method, which provides a strategy to efficiently solve the adjoint problem. We allow for arbitrary positions of both sources and receivers and for a refined grid around the source region to reduce errors in computed traveltimes. We show how gradients computed using the discrete adjoint method can be employed to perform either deterministic inversion, i.e., solving an optimization problem, or for a probabilistic (Bayesian) approach, i.e., obtaining a posterior probability density function. We show applications of our methodology on a set of synthetic examples both in 2D and 3D using the L-BFGS algorithm for the deterministic case and the Hamiltonian Monte Carlo algorithm for the probabilistic case.

SummaryHigh-quality maps of subsurface temperature and the geothermal gradient are useful when assessing the geothermal potential of a region. However, determining geothermal potential is a challenge when direct measurements of in-situ temperature and thermal property information are sparse and indirect geophysical methods are sensitive to a range of parameters, not just temperature. Here, we produce subsurface temperature maps of Ireland using a joint geophysical-petrological inversion, where seismic and other geophysical and petrophysical data are inverted directly for temperature in 1D columns and are collated into a pseudo 3D temperature volume. Additionally, the inversion produces new models for Moho and LAB depth and for the average crustal radiogenic heat production.To assess the robustness of the resulting temperature model, an uncertainty analysis has been performed by inverting all of the 1D columns for a range of reasonable input parameters applicable to the Irish crust (rather than the ‘best’ input parameters). The resulting uncertainty model suggests temperature estimates at 2 km depth in our model could vary by ± 2 to 5°C with an average of 3.5°C in most locations. The uncertainty model can be used to assess confidence in different regions of the temperature model. In addition, 3D forward modelling was performed to assess the lateral heat flow variations when compared to the purely 1D inversion. The upper-crustal geothermal gradient ranges from 20 to 40°C/km indicating a higher geothermal gradient for Ireland than previously reported with subsurface temperatures at 2 km depth > 60°C everywhere, sufficient for residential and industrial heating purposes. The temperature gradient is typically higher in areas with thinner lithosphere. However, in some locations, the observed geotherms are elevated further due to high radiogenic heat production in granitic rocks. In Northern Ireland, a thin lithosphere, coupled with a weakly conductive basalt layer overlying warm crust, results in elevated temperatures. These are the first temperature maps for Ireland that include uncertainty estimates, providing ranges for the subsurface temperature values, and demonstrate that the maps are comparable to direct independent borehole temperature measurements, which are observed to fall within the model uncertainty. Our new methodology provides workflows for determining the geothermal potential in areas with limited direct temperature measurements. The final temperature model with uncertainty provides useful constraints for geothermal exploration and utilisation on the island of Ireland.