JGR–Solid Earth

We demonstrate the effectiveness of seismic dense array surface wave focal spot imaging using USArray data from the western-central United States. We study dispersion in the 60–310 s period range and assess the image quality of fundamental mode Rayleigh wave phase velocity maps. We apply isotropic spatial autocorrelation models to the time domain zero lag noise correlation wavefield data at distances of about one wavelength. Local estimates of the phase velocity, its uncertainty, and the regression quality imply overall better ZZ relative to ZR or RZ results. The extension of the depth resolution compared to passive surface wave tomography is demonstrated by the inversion of three clustered dispersion curves from different tectonic units. We observe anisotropic surface wave energy flux and the influence of body wave energy, but sensitivity tests at 60 s targeting the data range, correlation component, and processing choices show that the ZZ focal spots yield consistent high-quality images compared to regional tomography results in the 60–150 s period range. In contrast, at 200–300 s the comparatively small scales of the imaged structures and the imperfect agreement with low-resolution global tomography results highlight the persistent challenge to reconcile imaging results based on different data sources, theories, and techniques. Our study shows that surface wave focal spot imaging is an accurate, robust, local imaging approach. Better control over clean autocorrelation fields can further improve applications of this seismic imaging tool for increased resolution of the elastic structure below dense seismic arrays.

No abstract is available for this article.

The Mongol-Okhotsk suture and the Adaatsag ophiolite belt are associated with the closure of the Mongol-Okhotsk paleo-ocean and are located within the Central Asian Orogenic Belt (CAOB) and Mongolia. The suture zone is flanked by volcanic-plutonic belts that host significant metallogenic zones, containing deposits of copper and gold. The tectonic evolution of this region is not fully understood and the lithospheric structure has been poorly studied. We analyze magnetotelluric data and generate a model of the electrical resistivity distribution across this region. Whereas the northern segment has a sharp transition from a high-resistivity upper crust to a low-resistivity lower crust, as observed beneath the Hangai Dome, the southern segment does not show this transition. A wide, low-resistivity zone (1–100 Ωm) imaged in the crust and lithospheric mantle is coincident with the Mongol-Okhotsk suture and ophiolite, revealing a clear and significant lithospheric-scale feature. Across the profile, numerous narrow, vertically oriented, low-resistivity features (1–100 Ωm) are spatially associated remarkably well with the proposed boundaries of tectonic domains. These results confirm ideas about the development of the CAOB. Some of these low-resistivity features are beneath the surface locations of large mineral zones, and likely represent fossil fluid pathways. We show congruent seismic velocity models for comparison and the results show a large-scale low-velocity anomaly (decrease of 2%–3%) that correlates with the location of the low-resistivity anomaly below the Mongol-Okhotsk suture. The geophysical results, combined with geological and geochemical data, provide insights into the structure of this region and help shed light on unanswered questions.

Vehicle-induced seismic waves, generated as vehicles traverse the ground surface, carry valuable information for imaging the underlying near-surface structure. These waves propagate differently in the subsurface depending on soil properties at various spatial locations. By leveraging wave propagation characteristics, such as surface-wave velocity and attenuation, this study presents a novel method for near-surface monitoring. Our method employs passing vehicles as active, non-dedicated seismic sources and leverages pre-existing telecommunication fibers as large-scale and cost-effective roadside sensors empowered by Distributed Acoustic Sensing (DAS) technology. A specialized Kalman filter algorithm is integrated for automated DAS-based traffic monitoring to accurately determine vehicles' location and speed. Then, our approach uniquely leverages vehicle trajectories to isolate space-time windows containing high-quality surface waves. With known vehicle (i.e., seismic source) locations, we can effectively mitigate artifacts associated with suboptimal distribution of sources in conventional ambient noise interferometry. Compared to ambient noise interferometry, our approach enables the synthesis of virtual shot gathers with a high signal-to-noise ratio and spatiotemporal resolution at reduced computational costs. We validate the effectiveness of our method using the Stanford DAS-2 array, with a focus on capturing spatial heterogeneity and monitoring temporal variations in soil seismic properties during rainfall events. Specifically, in non-built-up areas, we observed an evident decrease in phase velocity and group velocity and an increase in attenuation due to the rainfall. Our findings illustrate our method's sensitivity and resolution in discerning variations across different spatial locations and demonstrate that our method is a promising advancement for high-resolution near-surface imaging in urban settings.

The entire editorial board of the Journal of Geophysical Research-Solid Earth would like to sincerely thank all our colleagues who reviewed manuscripts for us in 2023. The hours they spent reading in order to provide insightful comments on manuscripts not only help improve the quality of these manuscripts but also ensure the scientific rigor of our reviewing process and eventually, of the research published in the field of Solid Earth Geophysics by our journal. With the advent of open science and AGU's data policy, the reviewing process now also encompasses checking the accessibility and availability of data and developed software. This is a key objective of AGU's FAIR (Findable, Accessible, Interoperable and Reusable) policy, for which many reviewers have provided suggestions that helped to improve the data presentation and availability, and which also fed the editorial board's reflection on the matter. Of course, we particularly appreciate timely reviews, particularly in light of the growing demands imposed by the increase of manuscripts submitted to Journal of Geophysical Research-Solid Earth. We received 1,869 submissions in 2023, and 1,472 reviewers contributed to their evaluation by providing 2,237 reviews in total. We are deeply thankful for all of their contributions. The editorial board of Journal of Geophysical Research-Solid Earth: Rachel Abercrombie, Yves Bernabé, Michael Bostock (former editor), Mark Dekkers, Anke Friedrich, Shin-Chan Han, Satoshi Ide, Isabelle Manighetti (former EIC), Fenglin Niu, Douglas R. Schmitt, Alexandre Schubnel (EIC), Jun Tsuchiya, and all the associate editors of JGR-SE.

Continental topography is dominantly controlled by a combination of crustal thickness and density variations. Nevertheless, it is clear that some additional topographic component is supported by the buoyancy structure of the underlying lithospheric and convecting mantle. Isolating these secondary sources is not straightforward, but provides valuable information about mantle dynamics. Here, we estimate and correct for the component of topographic elevation that is crustally supported to obtain residual topographic anomalies for the major continents, excluding Antarctica. Crustal thickness variations are identified by assembling a global inventory of 26,725 continental crustal thickness estimates from local seismological data sets (e.g., wide-angle/refraction surveys, calibrated reflection profiles, receiver functions). In order to convert crustal seismic velocity into density, we develop a parametrization that is based upon a database of 1,136 laboratory measurements of seismic velocity as a function of density and pressure. In this way, 4,120 new measurements of continental residual topography are obtained. Observed residual topography mostly varies between ±1 and 2 km on wavelengths of 1,000–5,000 km. Our results are generally consistent with the pattern of residual depth anomalies observed throughout the oceanic realm, with long-wavelength free-air gravity anomalies, and with the distribution of upper mantle seismic velocity anomalies. They are also corroborated by spot measurements of emergent marine strata and by the global distribution of intraplate magmatism that is younger than 10 Ma. We infer that a significant component of residual topography is generated and maintained by a combination of lithospheric thickness variation and sub-plate mantle convection. Lithospheric composition could play an important secondary role, especially within cratonic regions.

Changbaishan volcano (China/North Korea) is one of the most active and hazardous volcanic centers in Northeast Asia. Despite decades of intensive research, the eruption history of this stratovolcano remains poorly constrained. One of the major puzzles is the timing of the eruption that produced the Tianwen Yellow Pumice (TYP) deposit at the caldera rim. Here we identify a new cryptotephra layer in sediment core 13PT-P4 from the East Sea. Grain-specific major, minor, and trace element analyses of glass shards allow a clear correlation of this distal tephra to the proximal TYP deposit of Changbaishan. Age-depth modeling using radiocarbon (14C) dates of sediment bulk organic fractions and other tephrochronological markers from the sediment sequence constraints the age of the cryptotephra and thus the TYP eruption to 29,948–29,625 cal yr BP (95.4% confidence interval). Our findings lead to a revision of the history of Changbaishan explosive activity, and show that the volcano has been particularly active during ca. 51–24 ka BP in the last 100 ka. Using high resolution palaeo-proxy records, we find the TYP tephra almost coeval with regional to hemispheric-scale climatic changes known as Heinrich Event 3 (H3). With its precise age determination and wide geographic dispersion, the tephra offers a key isochron for dating records of past climatic changes and addressing the phasing relationships in environmental response to H3 across East Asia.

Crystal-hosted melt embayments and melt inclusions partially record magmatic processes at depth, but it is not always obvious how to interpret this record. One impediment is our incomplete understanding of how embayments and melt inclusions form. In this study, we investigate the formation mechanism of embayments and melt inclusions during quartz growth to quantify the relationship between the compositions of the entrapped and average melt. We study the growth of embayments and inclusions through direct numerical simulations that couple the growth of a crystal surface with the evolution of the concentrations of incompatible components in the surrounding melt. We find that H2O is more enriched in the interior of defects on crystal surface compared to the exterior. The resultant lower disequilibrium in the defect interior causes lower growth rate than in the exterior, elongating the defect into an embayment. If crystal growth stops, the composition in the embayment equilibrates with the average melt within days to months. If crystal growth continues until the embayment neck closes, a melt inclusion forms. The melt entrapped by both embayments and melt inclusions is enriched in incompatible components, such as H2O and CO2. In addition to inclusion size, the enrichment of incompatible components in melt inclusions also depends on component diffusivity and the crystal growth regime. High-diffusivity components like H2O have similar enrichment levels in all scenarios, while lower-diffusivity components like CO2 are more enriched in melt inclusions with smaller sizes or formed in continuous crystal growth.

Paleomagnetic records of middle Neoproterozoic (820 to 780 Ma) rocks display high amplitude directional variations that lead to large discrepancies in paleogeographic reconstructions. Hypotheses to explain these data include rapid true polar wander (TPW), a geomagnetic field geometry that deviates from a predominantly axial dipole field, a hyper-reversing field (>10 reversals/Ma), and/or undiagnosed remagnetization. To test these hypotheses, we collected 1,057 oriented cores over a 85 m stratigraphic succession in the Laoshanya Formation (Yangjiaping, Hunan, China). High precision U-Pb dating of two intercalated tuff layers constrain the age of the sediments between 809 and 804 Ma. Thermal demagnetization isolates three magnetization components residing in hematite which are not time-progressive but conflated throughout the section. All samples possess a north and downward directed component in geographic coordinates at temperatures up to 660°C that is ascribed to a Cretaceous overprint. Two components isolated above 660°C reveal distinct directional clusters: one is interpreted as a depositional remanence, while the other appears to be the result of a mid-Paleozoic (460 to 420 Ma) remagnetization, which is likely widespread throughout South China. The high-temperature directions are subtly dependent on lithology; microscopic and rock magnetic analyses identify multiple generations of hematite that vary in concentration and distinguish the magnetization components. A comparison with other middle Neoproterozoic paleomagnetic studies in the region indicates that the sudden changes in paleomagnetic directions, used elsewhere to support the rapid TPW hypothesis (ca. 805 Ma), are better explained by mixtures of primary and remagnetized components, and/or vertical axis rotations.

Constraining strain localization and the growth of shear fabrics within brittle fault zones at sub-seismic slip rates is important for understanding fault strength and frictional stability. We conducted direct shear experiments on simulated sandstone-derived fault gouges at an effective normal stress of 40 MPa, a pore pressure of 15 MPa, and a temperature of 100°C. Using a passive strain marker and X-ray Computed Tomography, we analyzed the spatial distribution of deformation in gouges deformed in the strain-hardening, subsequent strain-softening, and then steady-state regimes at displacement rates of 1, 30, and 1,000 µm/s. We developed a machine-learning-based automatic boundary detection method to recognize the shear fabrics and quantify displacement partitioning between each fabric element. Our results show fabrics oriented along R1 and Y (including boundary) shears are the two major fabric elements. At rates of 1 and 30 µm/s, the relative amount of displacement on R1 shears is displacement dependent, increasing to ∼20% of the total displacement up to the strain-softening stage, then decreasing to ∼10%–18% at the steady state. This trend is absent at the high rate where ∼18% of the displacement occurs on R1 shears throughout all investigated stages. At all rates, the relative amount of displacement on Y shears increases linearly with displacement to a total of larger than 50% at the steady state. Our study provides constraints on the development of the active slip zone, which is an important factor controlling heating and weakening associated with small-magnitude earthquakes with limited displacement (mm-dm), such as induced seismicity.

We present a new geodetic strain rate and rotation rate model for Greece that has been derived using a highly dense GPS velocity field. The spatial distribution and the resolved rates of the various velocity gradient tensor quantities provided updated constraints on the present-day upper crustal deformation in the region and revealed new details not reported previously. The spatial distribution of the second invariant demonstrated that the overall magnitude of strain rates is highest across two well-defined provinces. The first follows the North Anatolian Fault and its two branches within the north Aegean, crosses central Greece and through the Gulf of Corinth it terminates in western Greece, while the second encompasses the extensional province of western Turkey and the eastern Aegean Sea islands. Our estimates revealed that shearing affects some of the fault-bounded grabens of central Greece that lie to the SW of the North Aegean Basin implying considerable oblique extension. We identified a narrow region of counterclockwise rotation whose location and kinematics have been induced by the net effect across the intersection of the clockwise rotating domains of western and central Greece. The Aegean microplate and the Anatolian plate are separated by a wide transition zone which accommodates the curved stretching of the entire plate system. In both edges of the Hellenic forearc the dominant mode of crustal strain is E-W extension. We found that earthquakes of M ≥ 5.6 are spatially well-correlated with high-strain areas, indicating that strain rate mapping could be used to inform future probabilistic seismic hazard analyses.

The CsCl-type (B2) phase of FeSi (B2-FeSi) has been proposed as a candidate phase in the ultralow-velocity zones (ULVZs) at the base of the lower mantle and in the Earth's inner core. However, the elastic properties of B2-FeSi under relevant conditions remain unclear. Here we determine the density, elastic constants, and velocities of B2-FeSi at high pressures (90–390 GPa) and temperatures (3,000–6,000 K) relevant to the Earth's lower most mantle and the inner core, using first-principles molecular dynamics simulations. At the base of the lower mantle, B2-FeSi shows significantly lower velocities and a higher density than those of the ambient mantle. Mechanical mixing models suggest the presence of ∼27–39 vol% B2-FeSi in the silicate mantle is consistent with the reduced velocities and the elevated density of ULVZs observed seismically. On the other hand, the hcp-Fe and B2-FeSi mixture exhibits higher bulk sound velocity compared to the PREM under inner core conditions. Adding superionic H in the interstitial sites of B2-FeSi lowers its density but has little effect on the bulk sound velocity of B2-FeSi, precluding H-bearing B2-FeSi as a major component in the Earth's inner core.

Fault zones exhibit 3D variable thickness, a feature that remains inadequately explored, particularly with regard to the impact on fluid flow. Upon analyzing an analytic solution, we examine 3D thermal-hydraulic (TH) dynamical models through a benchmark experiment, which incorporates a fault zone with thickness variations corresponding to realistic orders of magnitude. The findings emphasize an area of interest where vigorous convection drives fluid flow, resulting in a temperature increase to 150°C at a shallow depth of 2.7 km in the thickest sections of the fault zone. Moreover, by considering various tectonic regimes (compressional, extensional, and strike-slip) within 3D thermal-hydraulic-mechanical (THM) models and comparing them to the benchmark experiment, we observe variations in fluid pressure induced by poroelastic forces acting on fluid flow within the area of interest. These tectonic-induced pressure changes influence the thermal distribution of the region and the intensity of temperature anomalies. Outcomes of this study emphasize the impact of poroelasticity-driven forces on transfer processes and highlight the importance of addressing fault geometry as a crucial parameter in future investigations of fluid flow in fractured systems. Such research has relevant applications in geothermal energy, CO2 storage, and mineral deposits.

We monitored the time history of the velocity change (dv/v) from 2002 to 2022 to investigate temporal changes in the physical state near the Parkfield Region of the San Andreas Fault throughout the interseismic period. Following the coseismic decrease in dv/v caused by the 2003 San Simeon (SS) and the 2004 Parkfield earthquakes, the dv/v heals logarithmically and shows a net long-term increase in which the current dv/v level is equivalent to, or exceeding, the value before the 2003 SS earthquake. We investigated this long-term trend by fitting the model accounting for the environmental and coseismic effects to the channel-weighted dv/v time series. We confirmed with the metrics of Akaike information criterion and Bayesian information criterion that the additional term of either a linear trend term, or a residual healing term for the case where the healing had not been completed before the SS earthquake occurred, robustly improved the fit to the data. We eventually evaluated the sensitivity of the dv/v time history to the GNSS-derived strain field around the fault. The cumulative dilatational strain spatially averaged around the seismic stations shows a slight extension, which is opposite to what would be expected for an increase in dv/v. However, the cumulative rotated axial strain shows compression in a range near the maximum contractional horizontal strain (azimuth of N35°W to N45°E), suggesting that the closing of pre-existing microcracks aligned perpendicular to the axial contractional strains would be a candidate to cause the long-term increase observed in the multiple station pairs.

We have extended the distributional finite difference method (DFDM) to simulate the seismic-wave propagation in 3D regional earth models. DFDM shares similarities to the discontinuous finite element method on a global scale and to the finite difference method locally. Instead of using linear staggered finite-difference operators, we employ DFDM operators based on B-splines and a definition of derivatives in the sense of distributions, to obtain accurate spatial derivatives. The weighted residuals method used in DFDM's locally weak formalism of spatial derivatives accurately and naturally accounts for the free surface, curvilinear meshing, and solid-fluid coupling, for which it only requires setting the shear modulus and the continuity condition to zero. The computational efficiency of DFDM is comparable to the spectral element method (SEM) due to the more accurate mass matrix and a new band-diagonal mass factorization. Numerical examples show that the superior accuracy of the band-diagonal mass and stiffness matrices in DFDM enables fewer points per wavelength, approaching the spectral limit, and compensates for the increased computational burden due to four Lebedev staggered grids. Specifically, DFDM needs 2.5 points per wavelength, compared to the five points per wavelength required in SEM for 0.5% waveform error in a homogeneous model. Notably, while maintaining the same accuracy in the solid domain, DFDM demonstrates superior accuracy in the fluid domain compared to SEM. To validate its accuracy and flexibility, we present various 3D benchmarks involving homogeneous and heterogeneous regional elastic models and solid-fluid coupling in both Cartesian and spherical settings.

Large block-bounding faults on the Tibetan plateau are significant geological structures that accommodate tectonic movements and accumulate stress, leading to large earthquakes. Quantifying the interseismic slip deficit rate helps to better assess the earthquake potential. We combine available InSAR (2015–2020) and interseismic GPS data to determine fault coupling along 14 major block-bounding faults. Spatially dense InSAR measurements remarkably improve the resolution of the coupling model. Combined with a GPS-constrained block model, we examine the performance of the inversion approach with the stress constraint and the common Laplacian smoothing based on both synthetic tests and real data. We suggest that, for continental strike-slip faults, adding the stress constraint can mitigate unphysical coupling distributions due to unreasonable assumptions or modeling artifacts, reducing the model uncertainty and improving the accuracy of the coupling model. This is particularly useful for segments featured by a highly heterogeneous distribution of coupling along the transition zone from locking to creeping region, partially-coupling segment, and junction zone between main and subsidiary faults. We present a large-scale fault coupling map along the major block-bounding faults on the northeastern and eastern Tibetan plateau, highlighting the distinct degrees of fault coupling and lateral variations. The collage of coupling maps along different faults demonstrates the kinematic features over a broad time scale during earthquake cycles ranging from early to late interseismic phases, such as the segments ruptured during the 2001 Kokoxili earthquake and the 1920 Haiyuan earthquake.

Surface topography is an important yet largely neglected aspect of the early evolution of cratons. The lateral accretion of cratonic nuclei inevitably forms orogenic belts that subsequently provide a sediment source for large, resource-rich intracratonic basins, but to date, geodynamic models have focused exclusively on lithospheric root processes. Here we use two-dimensional thermal-mechanical models to study the topography and lithospheric deformation during 50 Myr of compression of a cratonic nucleus, to simulate the lateral accretion phase of craton growth in the Neoarchean. Although the cratonic nucleus thickens slightly during the compression phase, most of the deformation occurs in the regions adjacent to the nucleus that have weaker lithosphere. Here, crustal thickness triples developing high topography in excess of 10 km without active erosion. Models with different initial rheological parameters will have different final topography and lithosphere geometry, but in general it is difficult to shorten and deform the depleted cratonic nucleus, unless there are significantly weak heterogeneities in the mantle lithosphere. We apply two quantitative analysis techniques to objectively evaluate a multitude of model outputs. Cross-correlation clustering (CCC) measures the degree of similarity between topography profiles and categorizes models based on the general topographic character. Six different topography families are possible in the context of our models and crustal strength is the most important parameter affecting the shape. From principal component analysis (PCA) we identify four dominant lithosphere geometries. When used together, these two methods provide distinct yet complementary information about the surface and subsurface deformation features in our models.

Clay-rich rocks are integral to subduction zone dynamics and of practical importance, for example, as barriers in nuclear waste and CO2 repositories. While the effects of swelling strain on the self-sealing capabilities of these rocks are relatively well-established, the implications of polar fluids interacting with charged clay particles on the frictional behavior, and the role of swelling stress in initiating slip in critically stressed faults, remain ambiguous. To address these uncertainties, we conducted triaxial friction experiments using saw-cut samples, with the upper half composed of Opalinus claystone (OPA) and the lower half of Berea sandstone (BER). The frictional strength of the non-wetted OPA-BER interface was estimated based on experiments at confining pressures of 4–25 MPa and constant axial loading rate (0.1 mm/min). Fluid injection friction experiments were performed using decane (non-polar fluid) or deionized water (polar fluid) at 10 and 25 MPa confining pressures and constant piston displacement control. Macroscopic mechanical data were complemented by distributed strain sensing on the sample surface. Compared to decane, the frictional strength of the OPA-BER interface tended to decrease when injecting water, which is attributed to phyllosilicate lubrication and the transition of the OPA from a solid rock to an incohesive mud near the saw-cut surface. When injecting water, slip was initiated during initial hydration of the OPA-BER interface, although the apparent stress state was below the yield stress. To explain this behavior, we propose that the swelling stress is a crucial factor that should be integrated into the effective stress model.

Crack interactions leading to shear localization were quantified using microstructural analysis for brittle faults and high-temperature ductile faults formed during experiments on quartz sandstone. In both faulting regimes, the nucleation of macroscopic faults results from the interactions of microfractures at two length scales in ensemble. Brittle faults nucleate when the longest mesoscale shear fractures and transgranular tensile cracks critically interact. In contrast, ductile faults nucleate when the longest mesoscale shear fractures and multi-grain scale intergranular shear cracks critically interact. For both faulting regimes, we conclude the interaction and coalescence of the longest mesoscale shear fractures is the fundamental process responsible for fault nucleation. Hence, mesoscale shear fractures, which accommodate the majority of axial strain prior to shear localization in both faulting regimes, also serve as the nucleus of macroscopic faults. Locally, the growth of the mesoscale shear fractures is promoted by the interaction and coalescence of the multi-grain scale cracks in both faulting regimes. We hypothesize that attainment of a critical microstructure for shear localization (i.e., local clustering of the longest microfractures) requires a characteristic amount of plastic axial strain, which depends on deformation conditions. In brittle faulting, distributed microfracturing is confined within limited regions of the rock volume, which expedites crack clustering and fault nucleation at low characteristic strains. In ductile faulting, distributed microfracturing occurs more uniformly throughout the rock volume, delaying shear localization to high characteristic strains. Accurate prediction of shear localization requires models that describe crack interactions of the largest flaws that account for crack clustering.

This paper aims to update our understanding of the carbon cycle in the Himalayas, the most intense collisional orogeny globally, by providing new insight into its impact on Cenozoic climate cooling through the use of isotopic variations in both organic and inorganic carbon and an isotopic mass balance model. Our results from 20 selected hot springs show that the relative contributions of dissolved carbon from the mantle, metamorphic decarbonization, aqueous dissolution, and soil organic matter are approximately 2%, 82%, 6%, and 10%, respectively. Approximately 87% ± 5% of CO2 generated in the deep crust precipitates as calcite, while approximately 5.5% ± 1% of this carbon is converted to biomass through microbial chemosynthesis at depths less than 2 km. Our random forests approach yielded a metamorphic carbon flux from the entire Himalayan orogenic belt of approximately 2.7 ∼ 4.5 × 1012 mol/yr. The minor CO2 released into the atmosphere (2.5 ∼ 4.2 × 1011 mol/yr) is comparable to the carbon consumption driven by Himalayan weathering. This paper provides new insights into deep carbon cycling, notably that approximately 93% of deeply sourced carbon is trapped in the shallow crust, rendering orogenic processes carbon neutral and possibly acting as one of the major triggers of long-term climate cooling in the Cenozoic.