Geophysical Journal International

SummaryThe aim of this study is to assess the potential of rotational and strain measurements to provide complementary information on seismic wave scattering, in addition to the conventional seismological observables. We begin by evaluating the accuracy of numerical solutions to the elastic wave equation, solved via the Spectral Element Method, for modeling wave propagation in 3D complex heterogeneous media. These simulations are benchmarked against predictions from the Radiative Transfer Equation (RTE), which models energy transport in scattering media. The comparison focuses on key scattering parameters: mean free path, diffusion onset, and temporal evolution of P/S energy partitioning. Three levels of velocity heterogeneity (10%, 17%, and 25%) are tested in both full-space and half-space configurations. The analysis highlights how scattering strength, numerical accuracy, and theoretical assumptions, such as those underlying the Born approximation, affect the agreement between the two modeling approaches. This comparison helps define the conditions under which RTE and wave equation-based simulations produce consistent results. Following this assessment, we analyze the energy envelopes of the displacement wavefield and its spatial gradients. The results demonstrate that rotational measurements preserve source-induced polarization longer than other observables. This persistence can provide valuable information for better constraining the source mechanism. Furthermore, analysis of the rotational components can provide complementary constraints on the medium’s elastic and scattering properties.

SummaryThis study develops a Bayesian optimisation method of coseismic slip distributions and the correlation lengths of the von Kármán autocorrelation function, which can realize more reliable regularisation grounded in geophysical analysis (VKR). To validate VKR and its dependence on the observation, synthetic tests using data from inhomogeneously distributed Global Navigation Satellite System (GNSS) stations and Synthetic Aperture Radar (SAR) were conducted. When observation stations surround the source fault, the assumed slip distribution and the correlation length were well recovered with some artificially extended slip due to the inhomogeneity of the observation network. Moreover, when observation stations exist only on one side of the fault, the method recovered the slip and other parameters with accuracy comparable to that of the surrounding-case scenario, albeit with slightly increased uncertainty of the parameters. This highlights the importance of uncertainty evaluation for slip and correlation length parameters, especially under biased observation networks. Applying existing Laplacian smoothing methods to the same experiments produced models globally consistent with VKR, demonstrating that the proposed method, despite additional non-linear parameters, achieves comparable estimation accuracy. Existing method showed isotropic correlations of slip variables, whereas VKR exhibited correlations elongated along the strike direction, reflecting its ability to independently regularise along strike and dip via correlation length parameters. These correlation patterns were most pronounced in deeper fault regions, where regularisation dominated over observational constraints. Cluster analysis of the Markov samples revealed that VKR captured a more diverse set of slip distribution models, with cluster differences most evident in deeper fault regions. These analyses underscore the importance of regularisation choice and its impact when evaluating and interpreting slip distributions. We also applied VKR to the actual data of the 2024 Noto Peninsula earthquake in Japan, observed by dense GNSS networks and SAR. The estimated slip distribution featured multiple slip areas in the eastern and western peninsula, consistent with previous studies and Laplacian smoothing-based results. The correspondence between estimated slip distribution and correlation length was confirmed, although it did not match previous empirical findings. To further examine the effect of observation network configuration, we performed additional inversions using a reduced set of actual data. The resulting slip distribution was smoother, and the estimated correlation lengths were larger and more uncertain. Our results reconfirm that simultaneous estimation of slip distribution and correlation length produces mutually consistent values, both of which depend on the spatial distribution of observation points. The results also demonstrate the strong dependence of estimated slip distribution and correlation lengths on the observation network. Therefore, future studies of coseismic slip self-similarity, using observed data, should incorporate parameter analyses that account for network resolution effects.

SummaryThe precise extraction of tidal signals from non-stationary gravity observations is a central challenge in geophysics, where accuracy is often limited by mode mixing in data preprocessing algorithms. This study evaluates the performance of the Time-Varying Filter-based Empirical Mode Decomposition (TVF-EMD) method to address this issue. We employed a progressive validation pipeline: the method was first verified on simulated signals, then rigorously tested against a high-fidelity benchmark from a superconducting gravimeter (SG), and finally applied to one month of continuous data from an Atom Gravimeter (AG) at the Yilan station. Results demonstrate that TVF-EMD dramatically suppresses mode mixing, with the energy of transient spikes in its residual being an order of magnitude lower than that from the conventional Ensemble Empirical Mode Decomposition (EEMD) method. The tidal signal reconstructed by TVF-EMD achieved the highest cross-correlation coefficient and the smallest root mean square error when compared to the theoretical gravity tide. Subsequent harmonic analysis confirmed that TVF-EMD yielded the lowest errors across all major tidal constituents. These findings validate TVF-EMD as a superior preprocessing framework for tidal analysis, particularly for enhancing the reliability of geophysical parameter inversion from short-duration records obtained with next-generation quantum sensors.

SummaryStructural boundaries are often the features of most interest geologically, but imaging them can be difficult due to wavefield scattering and interference caused by the sharp velocity contrasts. One example of this is the apparent Rayleigh-wave anisotropy (1-psi anisotropy) that has been observed near major structural boundaries using seismic arrays. The cause of the apparent anisotropy is the interference between the incident surface wave and waves scattered from velocity discontinuities. In this study, we first investigate the sensitivity of apparent anisotropy measurements to lateral boundary sharpness through 2D full waveform simulations. We demonstrate that 1-psi anisotropy can vary based on boundary sharpness, station spacing, and period of surface waves. We show that a misfit defined using triple-difference travel times, i.e. the difference in double-difference travel times between station pairs with opposite propagation directions, well characterizes the apparent anisotropy. The sensitivity kernel for this triple-difference misfit can be constructed using the adjoint method. We show that triple-difference travel times are mainly sensitive to velocity contrasts rather than absolute velocities, in contrast to double-difference travel times. With sensitivity kernels constructed, we demonstrate how triple-difference travel times can be combined with double-difference travel times into a tomography inversion. We show that by including triple-difference travel times, seismic inversions converge faster and resolve boundary and average structure better in early iterations, compared to using double-difference travel times alone. Recent advancements in dense array experiments could facilitate the application of this method to better delineate tectonic and basin structural boundaries.

SummaryWe present a probabilistic framework for evaluating earthquake forecasting models that use an alarm-based approach. In this approach, alarms are triggered by specific precursor signals. In a previous paper we compared such models and two ensemble models combining them in additive and multiplicative mode, with the ETAS (Epidemic Type Aftershock Sequence) forecasting model, which is defined in a probability-based approach, by making the latter to issue an alarm when the expected rate exceeds a predefined threshold. In this work we compare the alarm-based models with the ETAS and with another probability-based model, EEPAS (Every Earthquake a Precursor According to Scale) previously applied to Italy, using the testing procedures developed for probability-based models within the Collaboratory Study for Earthquake Predictability (CSEP) initiative. To do that, for the four alarm-based models, we compute empirical probabilities (frequencies) of Mw ≥ 5.0 earthquakes in Italy, inside and outside alarm time intervals issued by such models from 1990 to 2011. We then compare pseudo-prospectively the forecasting ability of all six models, by applying the CSEP tests on the time interval from 2012 to 2023. We found that the evaluation method used has a strong impact on the ranking of model performance. Probabilistic models like ETAS and EEPAS tend to score better under the CSEP testing framework whereas alarm-based models generally outperform probability-based ones when assessed using alarm-based metrics.

SummaryMonitoring subsurface fluid transport processes using electrical resistivity tomography (ERT) gained attention due to its sensitivity to variations in fluid content and temperature. Capturing these dynamic processes, such as tracer transport or pollutant dispersion, requires not only advanced imaging techniques but also efficient data acquisition strategies to ensure high-resolution and cost-effective monitoring. This study advances Optimal Experimental Design (OED) for transient flow monitoring using ERT. We developed three different strategies based on the Compare-R algorithm, and present as well as compare their results. By optimizing the selection of measurement configurations, these strategies aim to maximize the information content of the acquired data while minimizing redundancy and acquisition costs. The data-driven approach uses data from prior time steps to generate focusing masks. The model-driven approach incorporates transport models to target critical areas, and the hybrid approach combines both methods, iteratively refining the transport model based on the acquired data. We use synthetic studies to demonstrate the model-driven strategy”s advantage in identifying transport-affected areas, especially for rapid processes or longer monitoring intervals. The data-driven approach performs well when masked and process-affected areas overlap, but struggles with spatial delays otherwise. The hybrid OED effectively combines the strengths of both methods, detecting discrepancies between modeled and observed fluid flow, and refining transport simulations for subsequent steps, thereby ensuring reliable optimization for dynamic monitoring scenarios.

SummaryShear wave splitting, caused by seismic anisotropy, provides insights into convective processes in the mantle. While upper mantle anisotropy beneath the continents is well-resolved with high lateral resolution, its characterization beneath the oceans remains limited due to the paucity of seismic stations. In this study, we leverage the sensitivity of P waves that convert to S upon reflection at the surface (PS) to infer seismic anisotropy near their bounce point. We measure PS splitting automatically from a global dataset that includes all earthquakes with magnitudes ≥5.9 from January 1995 to present, collected from 25 datacenters, totaling approximately 5,800 events and 18 million three-component seismograms. After careful quality control, we obtain 889 PS splitting measurements from this dataset, mostly for regions without other splitting constraints. These regions include the South China Sea, the Philippine Sea, the western Pacific Ocean, the western Atlantic Ocean near the Caribbean, and the South Atlantic Ocean. Where independent SKS constraints are available, these are in general agreement with our PS results. Comparisons of PS fast directions to azimuthal seismic anisotropy derived from surface waves and the direction of absolute plate motion (APM) show a dominant contribution of asthenospheric deformation to the measurements. This agreement is particularly pronounced beneath the oceans, away from subduction zones, for instance, across much of the Pacific basin. In contrast, splitting patterns near subduction zones are generally more complex: fast directions are most often trench-parallel or oblique, and only rarely trench-perpendicular, suggesting that a trench-parallel component of shear deformation is commonly generated in subduction settings, as commonly inferred from SKS splitting. We make all measurements publicly available as a data product, along with detailed metadata to enable future work.

SummaryMass redistribution during earthquake rupture, along with subsequent wave propagation, perturbs Earth’s gravity field, generating so-called ‘prompt elasto-gravitational signals’ (PEGS) and ‘prompt gravity-strain signals’ (PGSS). These signals are detectable on accelerometers and gradiometers before P-wave arrivals, and therefore offer the potential for early-warning systems and rapid assessment of event magnitude and tsunamigenic risk. Despite their significance, numerical modelling of PEGS and PGSS has been restricted to 1D Earth models, assuming negligible effects of 3D heterogeneities.In this study, we utilise a Spectral-Infinite-Element method to compute these PEGS and PGSS for both 1D and 3D Earth models. Kernels of sensitivity to model heterogeneity are investigated using the adjoint method, which is re-formulated to write the self-gravitating adjoint equations directly in terms of density and elastic perturbations only. We find that PEGS and PGSS display the greatest sensitivity to model perturbations below the source and receiver; however, sensitivity is generally weak, validating previous assumptions that 1D models are sufficient.

SummaryThe strength of the subduction interface plays a major role in controlling subduction dynamics on both local and global scales. While previous studies have primarily examined interface strength in two-dimensional models, natural subduction zones are inherently three-dimensional, with interface strength varying along-strike due to spatial differences in factors such as sediment input. Here, we use the geodynamic code ASPECT to conduct fully dynamic 3D subduction models in which interface strength varies along-strike. We find that the interaction between strong and weak segments of the interface leads to a narrower range of convergence velocities while broadening the range of viable interface stresses compared to 2D or homogeneous 3D models. Stronger segments, when adjacent to weaker ones, exhibit increased convergence velocities. This promotes higher interface stresses and facilitates the subduction of otherwise stagnant strong segments. We find that the interface viscosity of the strong segment controls the baseline stress, whereas the viscosity contrast along-strike controls the magnitude of amplification of the stress due to velocity increases. The elevated interface stresses at strong segments also generate greater compressional forces in the overriding plate than expected from 2D models. Combined with along-strike variations in convergence velocity, this results in trench migration, with stronger segments displaying more advanced trench positions relative to weaker segments. Possible natural analogs include the Bolivian Orocline in the central Andes and the Lesser Antilles, both of which show enhanced overriding plate compression and trench advance in areas of reduced sediment supply.

SummaryRealistic models of earthquake sequences can be simulated by assuming faults governed by rate-and-state friction embedded in an elastic medium. Exploring the possibility of using such models for earthquake forecasting is challenging due to the difficulty of integrating Partial Differential Equation (PDE) models with sparse, low-resolution observational data. This paper presents a machine-learning-based reduced-order model (ROM) for earthquake sequences that addresses this limitation. The proposed ROM captures the slow/fast chaotic dynamics of earthquake sequences using a low-dimensional representation, enabling computational efficiency and robustness to high-frequency noise in observational data. The ROM’s efficiency facilitates effective data assimilation using the Ensemble Kalman Filter (EnKF), even with low-resolution, noisy observations. Results demonstrate the ROM’s ability to replicate key scaling properties of the sequence -namely the magnitude-frequency, moment-duration, and moment-area relationships- and to estimate the distributions of fault slip rate and state variable, enabling predictions of large events in time and space with uncertainty quantification. These findings underscore the ROM’s potential for forecasting and for addressing challenges in inverse problems for nonlinear geophysical systems.

SummaryThis study introduces a novel inversion approach to resolve the layered anisotropic structure of the Earth’s crust and upper mantle using harmonic patterns observed in Ps receiver functions (RFs). Designed for dense seismic networks, our method effectively captures the complexity and deformation of subsurface structures by leveraging the decomposition of RF patterns into five harmonic functions representing distinct terms in a harmonic regression model. Our approach combines residual weighting at individual and coherency weighting between stations to reduce noise and enhance structural signal, followed by a direct inversion of harmonic patterns alongside a multi-phase stack. This process is optimized using a Markov chain Monte Carlo framework to iteratively explore model space and define posterior distributions, which yields robust model parameter estimates with quantified uncertainties. We validated our approach with both synthetic models and real data from the dense SEISConn seismic transect in Northern Connecticut. The synthetic test highlighted the reliability of coherency weighting in reducing noise. Real data analysis revealed anisotropic and structural features consistent with geological expectations and previous studies, including shallow anisotropy in western Connecticut, a shallowing Moho towards the Hartford basin, and indications for a west-dipping layer beneath the Moho in the lithospheric mantle. Our approach offers a promising tool for revealing the details of anisotropic features beneath dense seismic profiles, facilitating insights into tectonic history and lithospheric deformation.

SummaryLarge-scale geodetic monitoring of volcanic earthquakes is essential for understanding the physical mechanisms governing volcanic activity and magma migration. Recently, a significant earthquake swarm occurred around the Scotia plate. Geodetic data of 14 permanent GNSS stations on the Antarctic Peninsula, King George Island, the South Sandwich Islands, and South America were collected and analyzed, to monitor the crustal deformation of King George Island, the expansion of Bransfield Strait, the drift of Scotia plate, and sea level anomalies. Tidal data from four permanent tide stations were analyzed to monitor sea level anomalies. Results showed that after earthquakes the King George Island’s movement speed increased tenfold and its direction altered by 90 degrees. Land surface fluctuations in southeast King George Island were observed a year before the earthquakes, followed by continuous uplift. A combinatorial model including a point pressure source and expanding dike fit well with new geodetic monitoring data, revealing the impact of volcanic activity on this region. Geodetic monitoring and modeling quantitatively depicted the pre-seismic, co-seismic, and post-seismic phases of geological changes, providing new evidence and insights into the complex geological structures.

SummaryA 220 km long active seismic profile crossing the Saronic and Corinth Gulfs was performed using 35 4C Ocean Bottom Seismographs (OBS) and 4 3C stand-alone land stations. We recorded shots fired along line at every 120 m from a 48-l airgun array of 51 bar-m power. The velocity model was developed by first break tomographic inversion, followed by kinematic and dynamic ray tracing. The final velocity model was used to prestack depth migrate the Common Receiver Gathers. Moho depth below the central Corinth Basin is located at 32 km, thinning to 22 km below the Lechaio Gulf, at the transition to the Saronic Basin. In the Saronic domain Moho is found at 19 to 21 km depth, whereas below we identified a low velocity upper mantle of Vp 7.5 to 7.6 km/s, extending to 34 km depth. This is a low velocity asthenosphere wedge that intruded from the Cyclades region below the Saronic Gulf, driving the volcanic activity. It does not extend in the Corinth Rift domain to the west, where Pn is 8.0 km/s. Two major extensional detachments were mapped along the profile: one to the NW in the central Corinth Gulf, east of Galaxidi, corresponding to the onshore Itea-Amfissa Detachment, having a throw of more than 2000 m; the other to the SE, west of Agios Georgios Island, separates the Cycladic Metamorphic Core Complex from the non-metamorphic internal nappes of the Hellenides. Thrusting to the NW observed in the SE Saronic Basin has doubled the thickness of the high velocity limestones. At 20 km SE of Agios Georgios Island a major dextral strike slip fault was mapped. Normal faults with more than 1 km throw are observed around the Isthmus of Corinth trending E-W and WNW-ESE. Maximum thickness of 2000 m of Middle Miocene-Quaternary sediments (Vp 2.1–2.8 km/s and 3.4 km/s) is observed in the Corinth Gulf and minimum of 200 m of Quaternary sediments above the western Cyclades. Beneath the volcanoes of Aegina, Methana, Paphsanias, and Sousaki the crust has a lateral Vp velocity increase of nearly 3 per cent, whereas depth migrated data show high reflectivity, indicating magma intrusions through the crust. Mesozoic limestones with Vp 5.2 to 5.8 km/s and 1800 to 2000 m thickness occur along the profile, corresponding probably to the Tripolis external carbonate platform, whereas limestones and flysch with Vp 4.3 to 5.5 km/s are overlying, probably corresponding to pelagic sequences like the Pindos nappe. The Saronic domain is characterized by arc-parallel NW-SE structures, hosting the volcanic arc, and is tectonically controlled by magma uplift and thermally triggered deformation. The Corinth domain is affected by E-W transverse-oblique structures of a rapidly developing continental rift, and deformation driven by fracturing of a brittle crust of the Central Hellenic Shear Zone.

SummaryTo better understand Eocene ophiolite emplacement and subduction initiation in northeastern Zealandia, we analysed ambient noise to image shallow (0–3 km) shear wave velocity structures of and beneath an ophiolite nappe in southern Grande Terre, New Caledonia. We assessed the uncertainties of each dispersion curve to obtain stable dispersion curves at short periods (<1 s) from a network of 17 seismic stations, whose average interstation distance is ∼15 km. We obtained 1D velocity profiles and interpolated them to generate 2D transects with a lateral resolution <2 km. Two velocity discontinuities were imaged at depths of 100 m and 400–700 m, representing surface regolith and the base of the ophiolite nappe, respectively. The ophiolite nappe is underlain by continental basement rocks in the centre of the island and sedimentary rocks near the east and west coasts. The base of the nappe shallows at ∼2.5° westward to the surface at its southwest flank. Based on the geometry of the ophiolite nappe, we suggest a down-going gravity-driven emplacement mechanism, and note similarities to allochthons in Reinga Basin and Raukumara Basin of northern New Zealand. The ophiolite nappe and underlying bedrock are more fractured on their east flank due to syn/post emplacement deformation, isostatic adjustment and present flexural bending.

SummarySeismological inversion traditionally targets either source parameters, such as location and moment tensor, or structural parameters, such as velocity and anisotropy. However, the natural formulation of Full-Waveform Inversion, often used for high-resolution structural model estimation, is to jointly invert for source and structural parameters. The common practice of holding source parameters, after initial estimation, fixed throughout the inversion inherently leads to biased solutions of the structural model, and vice versa. Whereas a joint inversion suffers from severe non-uniqueness, we demonstrate that leveraging the large amounts of data available from Distributed Acoustic Sensing (DAS) can yield robust and unbiased estimations of source and structural parameters, provided an appropriate misfit function and optimisation scheme are used. We show how the size of the data space and eventual convergence can be improved by supplementing the phase misfit objective function with amplitude information. To this end, we formulate a new misfit function, the normalised envelope. To support native DAS data implementations, we calculate the adjoint sources for the new misfit function when defined directly on strain or strain-rate data. We also show how a new approach to preconditioning as part of the L-BFGS optimisation scheme allows for effective updates of all parameters in the same iteration, despite enormous differences in their relative importance. We test our approach in a challenging synthetic noisy 2D scenario, showing a considerable reduction in source parameter errors and an improved S-wave velocity model. We also show a 3D synthetic case with an idealised DAS recording array, demonstrating a significant reduction of source parameter errors using realistic initial estimates and structural model errors. We argue that the proposed methodology can be used to improve the quality of earthquake catalogues and high-resolution structural models in seismically active regions, especially at the local-to-regional scale. None the less, computational cost remains a major challenge of the method.

SummaryTo understand the melt source of hotlines with asynchronous volcanoes, we investigate the lithospheric structure of the Cameroon Volcanic Line (CVL), an intraplate hotline without age progression stretching from the Atlantic Ocean into Central Africa. We analyze Bouguer gravity anomalies from the World Gravity Model 2012 using the 2‐D power spectrum techniques and 2-D forward modeling to estimate the crustal and lithospheric thickness. We find: (1) thin crust (20–30 km) beneath the oceanic CVL; (2) thick crust (30–43 km) beneath the continental CVL and the Oubanguides Belt, and thicker crust (43–50 km) beneath the Congo Craton; (3) thin lithosphere (90–120 km) beneath the oceanic CVL and thinner lithosphere (75–90 km) beneath the continental CVL; and (4) thicker lithosphere (150–234 km) beneath the Congo Craton. Our seismically constrained forward models reveal a delaminated body beneath the continental CVL and a sharp transition from thick lithosphere beneath the Congo Craton to thin lithosphere beneath the Oubanguides Belt. We interpret that the thin lithosphere beneath the continental CVL is a result of lithospheric delamination. The delaminated body in the uppermost mantle deflects rising mantle plume material, resulting in the Y-shaped distribution of continental volcanoes. Edge-Driven Convection (EDC) resulting from the sharp gradient in lithospheric thickness between the Congo Craton and the Oubanguides Belt focuses the plume material beneath thin lithosphere, producing the continental CVL. The southern volcanoes of the continental CVL are formed from the southward deflection of plume material by the delaminated body, with melt ascent facilitated by the lithospheric-scale Central African Shear Zone. The northward-directed plume material forms the distinct Biu Plateau, and the eastward-deflected plume material forms the Adamawa Plateau. With a continuous influx of plume material beneath the thin continental lithosphere, for mass to be conserved, part of the plume material defiles the gradient of the thicker oceanic lithosphere adjacent to the Congo Craton to flow oceanward. The oceanward flow of plume material is modulated by upwellings from EDC, producing the oceanic CVL, which explains the oceanward decrease in the timing of the onset of volcanism. We therefore conclude that only the continental CVL lacks age progression resulting from the complex interaction of the rising plume with the delaminated body and the lithospheric architecture.

SummaryInferring the spatio-temporal distribution of slip during earthquakes remains a significant challenge due to the high dimensionality and ill-posed nature of the inverse problem. As a result, finite-source inversions typically rely on simplified assumptions. Moreover, in the absence of ground-truth measurements, the performance of inversion methods can only be evaluated through synthetic tests. Laboratory earthquakes offer a valuable alternative by providing “simulated real data” and ground truth observations under controlled conditions, enabling a more reliable evaluation of source inversion procedures. In this study, we present static and quasi-static slip inversion results from data recorded during laboratory earthquakes. Each event is instrumented with 20 accelerometers along the fault, and the recorded acceleration data are used to invert for the slip history. We consider two different types of Green’s functions (GF): simplistic GF assuming a homogeneous elastic half-space and realistic GF computed by finite element modeling of the experimental setup. The inversion results are then compared to direct observations of fault slip and rupture velocity obtained independently during the experiments. Our results show that, regardless of the GF used, the inversions fit well with the data and result in small formal uncertainties of model parameters. However, only the inversion with realistic GF yields slip distributions consistent with the true fault slip measurements and successfully recovers the distribution of rupture velocity along the fault. These findings emphasize the critical role of GF selection in accurately resolving slip dynamics and highlight an important distinction in Bayesian inversion: while posterior uncertainty quantification is essential, it does not guarantee accuracy, especially if forward modeling uncertainties are not properly accounted for. Thus, confidence in inversion results must be paired with careful modeling choices to ensure physical reliability.

SummaryThe regional seismic travel-time (RSTT) model predicts travel times of regional seismic phases accounting for three-dimensional structure of the crust and the upper mantle on a global scale. Previous versions of the RSTT model have been implemented using nodes separated by ∼1° spacing across the globe. A regional-scale study using regional Pn and Pg travel times across Israel and the Middle East demonstrated that data driven, systematic grid refinement reduces travel time residuals and enhances resolution of smaller tectonic features in regions having dense ray coverage. High density Pn ray coverage in the western US, Europe, Middle East, and East Asia can likewise provide the resolution that allows systematic global grid refinement of the RSTT model. In this study, we use a large number of Pn ray paths originating from events located with an epicentral location uncertainty of 25 km (GT25) or better. We conduct targeted grid refinements at 1.0°, 0.5°, 0.25°, and 0.125° on a global scale, producing a refined RSTT model that yields a 21.6% reduction in median event location error in Europe and the Middle East, when compared with the original global RSTT model presented in Begnaud et al. (2021a). The new model also resolves finer tectonic structures in regions with high Pn ray density.

SummaryWe present and apply a pseudo trans-dimensional inversion method for 3D geometrical gravity inversion, in which the number of rock units, their geometry, and their density can vary during sampling. The method is designed for efficient exploration of the model space and to infer the presence and properties of units not directly observable but detectable with geophysical data. Sampling relies on a non-reversible Metropolis-Hastings algorithm, during which rock units can be added or removed from the model, interface geometries are perturbed using random fields, and densities are sampled from distributions informed by prior information. To visualise the space of sampled models and to aid interpretation, a workflow is proposed that combines dimensionality reduction with the clustering of models in families. The capabilities of the inversion method are evaluated using two synthetic cases. The first is a motivating example aimed at recovering an intrusion missing from the prior model. It features a horizontal layer-cake where fixed-dimensional inversion fails to adequately fit the data and sample models close to the true model, while the proposed pseudo trans-dimensional approach is much more successful. The second case investigates the recovery of two missing units and the capability to overcome prior model biases. Results show the potential of our method to infer the presence of unseen geological features such as intrusions. However, they suggest that with biased prior geological modelling, it may be challenging to infer with certainty the presence of more than two previously unknown rock units at depth.

SummaryWhen highly viscous fluids are present in a rock medium, the viscous effect of such fluids cannot be neglected in the propagation of elastic waves. In this paper, a fractional thermoporoelastic theory is newly proposed, which is a further improvement of the two-temperature generalized thermoporoelastic theory. Firstly, by introducing the Kelvin-Voigt model into the stress-strain constitutive equation, the viscous effect of highly viscous fluids is considered. Then, fractional derivatives are introduced into the heat conduction equations of the solid and fluid phase to consider the anomalous heat conduction caused by the viscous effect in rock media. Plane wave analysis method is adopted to obtain the phase velocity and attenuation factor of four longitudinal waves (P1, P2, T1, T2). Numerical results show that the introduction of fluid viscosity leads to the appearance of new relaxation peaks in the P wave at high frequencies, and the introduction of fractional derivatives causes a decrease in the phase velocity and attenuation factor of T waves. The results provide a reference for further research on the wave propagation in rock media containing highly viscous fluids.