ALBERT

Description: 〈span〉〈div〉SUMMARY〈/div〉The 2017 July 20, 〈span〉M〈/span〉〈sub〉w〈/sub〉6.6 Bodrum–Kos earthquake occurred in the Gulf of Gökova in the SE Aegean, a region characterized by N–S extension in the backarc of the easternmost Hellenic Trench. The dip direction of the fault that ruptured during the earthquake has been a matter of controversy where both north- and south-dipping fault planes were used to model the coseismic slip in previous studies. Here, we use seismic (seismicity, main shock modelling, aftershock relocations and aftershock mechanisms using regional body and surface waves), geodetic (GPS, InSAR) and structural observations to estimate the location, and the dip direction of the fault that ruptured during the 2017 earthquake, and the relationship of this event to regional tectonics. We consider both dip directions and systematically search for the best-fitting locations for the north- and south-dipping fault planes. Comparing the best-fitting planes for both dip directions in terms of their misfit to the geodetic data, proximity to the hypocenter location and Coulomb stress changes at the aftershock locations, we conclude that the 2017 earthquake ruptured a north-dipping fault. We find that the earthquake occurred on a 20–25 km long, ∼E–W striking, 40° north-dipping, pure normal fault with slip primarily confined between 3 and 15 km depth, and the largest slip exceeding 2 m between depths of 4 and 10 km. The coseismic fault, not mapped previously, projects to the surface within the western Gulf, and partly serves both to widen the Gulf and separate Kos Island from the Bodrum Peninsula of SW Anatolia. The coseismic fault may be an extension of a mapped, north-dipping normal fault along the south side of the Gulf of Gökova. While all of the larger aftershocks are consistent with N–S extension, their spatially dispersed pattern attests to the high degree of crustal fracturing within the basin, due to rapid trenchward extension and anticlockwise rotation within the southeastern Aegean.〈/span〉

Description: 〈span〉〈div〉SUMMARY〈/div〉The harmonic electromagnetic noise produced by anthropic electrical structures is a critical component of the global noise affecting geophysical signals and increasing data uncertainty. It is composed of a series of harmonic signals whose frequencies are multiple integers of the fundamental frequency specific to the electrical noise source. To date, most model-based noise removal strategies assume that the fundamental frequency constraining the harmonic noise is single and constant over the duration of the geophysical record. In this paper, we demonstrate that classical harmonic processing methods lose efficacy when these assumptions are not valid. We present several surface nuclear magnetic resonance field data sets, which testify the increasing probability of recording the harmonic noise with such multiple or unstable frequency content. For each case (multiple frequencies or unstable frequency) we propose new processing strategies, namely, the 〈span〉2-D grid-search〈/span〉 and the 〈span〉segmentation〈/span〉 approach, respectively, which efficiently manage to remove the harmonic noise in these difficult conditions. In the process, we also apply a fast frequency estimator called the Nyman, Gaiser and Saucier estimation method, which shows equivalent performance as classical estimators while allowing a reduction of the computing time by a factor of 2.5.〈/span〉

Description: 〈span〉〈div〉SUMMARY〈/div〉Combinatorial methods are used to determine the spatial distribution of earthquake magnitudes on a fault whose slip rate varies along strike. Input to the problem is a finite sample of earthquake magnitudes that span 5 kyr drawn from a truncated Pareto distribution. The primary constraints to the problem are maximum and minimum values around the target slip-rate function indicating where feasible solutions can occur. Two methods are used to determine the spatial distribution of earthquakes: integer programming and the greedy-sequential algorithm. For the integer-programming method, the binary decision vector includes all possible locations along the fault where each earthquake can occur. Once a set of solutions that satisfy the constraints is found, the cumulative slip misfit on the fault is globally minimized relative to the target slip-rate function. The greedy algorithm sequentially places earthquakes to locally optimize slip accumulation. As a case study, we calculate how earthquakes are distributed along the megathrust of the Nankai subduction zone, in which the slip rate varies significantly along strike. For both methods, the spatial distribution of magnitudes depends on slip rate, except for the largest magnitude earthquakes that span multiple sections of the fault. The greedy-sequential algorithm, previously applied to this fault (Parsons et al., 2012), tends to produce smoother spatial distributions and fewer lower magnitude earthquakes in the low slip-rate section of the fault compared to the integer-programming method. Differences in results from the two methods relate to how much emphasis is placed on minimizing the misfit to the target slip rate (integer programming) compared to finding a solution within the slip-rate constraints (greedy sequential). Specifics of the spatial distribution of magnitudes also depend on the shape of the target slip-rate function: that is, stepped at the section boundaries versus a smooth function. This study isolates the effects of slip-rate variation along a single fault in determining the spatial distribution of earthquake magnitudes, helping to better interpret results from more complex, interconnected fault systems.〈/span〉

Description: 〈span〉〈div〉SUMMARY〈/div〉The southcentral Hikurangi subduction margin (North Island, New Zealand) has a wide, low-taper accretionary wedge that is frontally accreting a 〉3-km-thick layer of sediments, with deformation currently focused near the toe of the wedge. We use a geological model based on a depth-converted seismic section, together with physically realistic parameters for fluid pressure, and sediment and décollement friction based on laboratory experiments, to investigate the present-day force balance in the wedge. Numerical models are used to establish the range of physical parameters compatible with the present-day wedge geometry and mechanics. Our analysis shows that the accretionary wedge stability and taper angle require either high to moderate fluid pressure on the plate interface, and/or weak frictional strength along the décollement. The décollement beneath the outer wedge requires a relatively weaker effective strength than beneath the inner (consolidated) wedge. Increasing density and cohesion with depth make it easier to attain a stable taper within the inner wedge, while anything that weakens the wedge—such as high fluid pressures and weak faults—make it harder. Our results allow a near-hydrostatic wedge fluid pressure, sublithostatic fluid overpressure at the subduction interface, and friction coefficients compatible with measurements from laboratory experiments on weak clay minerals.〈/span〉

Description: 〈span〉〈div〉SUMMARY〈/div〉We present a new methodology to compute the gravitational fields generated by tesseroids (spherical prisms) whose density varies with depth according to an arbitrary continuous function. It approximates the gravitational fields through the Gauss–Legendre Quadrature along with two discretization algorithms that automatically control its accuracy by adaptively dividing the tesseroid into smaller ones. The first one is a preexisting 2-D adaptive discretization algorithm that reduces the errors due to the distance between the tesseroid and the computation point. The second is a new density-based discretization algorithm that decreases the errors introduced by the variation of the density function with depth. The amount of divisions made by each algorithm is indirectly controlled by two parameters: the distance-size ratio and the delta ratio. We have obtained analytical solutions for a spherical shell with radially variable density and compared them to the results of the numerical model for linear, exponential, and sinusoidal density functions. The heavily oscillating density functions are intended only to test the algorithm to its limits and not to emulate a real world case. These comparisons allowed us to obtain optimal values for the distance-size and delta ratios that yield an accuracy of 0.1 per cent of the analytical solutions. The resulting optimal values of distance-size ratio for the gravitational potential and its gradient are 1 and 2.5, respectively. The density-based discretization algorithm produces no discretizations in the linear density case, but a delta ratio of 0.1 is needed for the exponential and most sinusoidal density functions. These values can be extrapolated to cover most common use cases, which are simpler than oscillating density profiles. However, the distance-size and delta ratios can be configured by the user to increase the accuracy of the results at the expense of computational speed. Finally, we apply this new methodology to model the Neuquén Basin, a foreland basin in Argentina with a maximum depth of over 5000 m, using an exponential density function.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉The dynamics of dyke emplacement are typically modeled by assuming an elastic rheology for the host rock. However, the resulting stress field predicts significant shear failure in the region surrounding the dyke tip. Here, we model the dyking process in an elastic-perfectly plastic host rock in order to simulate distributed shear fracturing and subsequent frictional slip on the fracture surfaces. The fluid mechanical aspects of the magma are neglected as we are interested only in the fracture mechanics of the process. Magma overpressure in dykes is typically of the same order of magnitude as the yield stress of the host rock in shear, especially when the pressure effect of volatiles exsolving from the magma is taken into account. Under these conditions, the plastic deformation zone has spatial dimensions that approach the length of the dyke itself, and concepts based on linear elastic fracture mechanics (LEFM) no longer apply. As incremental plasticity is path dependent, we describe two geologically meaningful endmember cases, namely dyke propagation at constant driving pressure, and gradual inflation of a pre-existing crack. For both models, we find that plastic deformation surrounding the fracture tip enhances dyke opening, and thus increases the energy input into the system due to pressure work integrated over the fracture wall. At the same time, energy is dissipated by plastic deformation. Dissipation in the propagation model is greater by about an order of magnitude than it is in the inflation model because the propagating dyke tip leaves behind it a broad halo of deformation due to plastic bending and unbending in the relict process zone. The net effect is that plastic deformation impedes dyke growth in the propagation model, while it enhances dyke growth in the inflation model. The results show that, when the plastic failure zone is large, a single parameter such as fracture toughness is unable to capture the physics that underpin the resistance of a fracture or dyke against propagation. In these cases, plastic failure has to be modeled explicitly for the given conditions. We provide analytical approximations for the propagation forces and the maximum dyke aperture for the two endmember cases, that is, the propagating dyke and the dyke formed by inflation of a crack. Furthermore, we show that the effect of plasticity on dyke energetics, together with an overestimate of magma pressure when interpreting dyke aspect ratios using elastic host rock models, offers a possible explanation for the long-standing paradox that laboratory measurements of fracture toughness of rocks consistently indicate values about two orders of magnitude lower than those derived from dyke observations.〈/span〉

Description: 〈span〉〈div〉SUMMARY〈/div〉Building geomechanical models for induced seismicity in complex reservoirs poses a major challenge, in particular if many faults need to be included. We developed a novel way of calculating induced stress changes and associated seismic moment response for structurally complex reservoirs with tens to hundreds of faults. Our specific target was to improve the predictive capability of stress evolution along multiple faults, and to use the calculations to enhance physics-based understanding of the reservoir seismicity. Our methodology deploys a mesh-free numerical and analytical approach for both the stress calculation and the seismic moment calculation. We introduce a high-performance computational method for high-resolution induced Coulomb stress changes along faults, based on a Green's function for the stress response to a nucleus of strain. One key ingredient is the deployment of an octree representation and calculation scheme for the nuclei of strain, based on the topology and spatial variability of the mesh of the reservoir flow model. Once the induced stress changes are evaluated along multiple faults, we calculate potential seismic moment release in a fault system supposing an initial stress field. The capability of the approach, dubbed as MACRIS (〈strong〉M〈/strong〉echanical 〈strong〉A〈/strong〉nalysis of 〈strong〉C〈/strong〉omplex 〈strong〉R〈/strong〉eservoirs for 〈strong〉I〈/strong〉nduced 〈strong〉S〈/strong〉eismicity) is proven through comparisons with finite element models. Computational performance and suitability for probabilistic assessment of seismic hazards are demonstrated though the use of the complex, heavily faulted Gullfaks field.〈/span〉

Description: 〈span〉〈div〉SUMMARY〈/div〉We test the feasibility of GPS-based rapid centroid moment tensor (GPS CMT) methods for Taiwan, one of the most earthquake prone areas in the world. In recent years, Taiwan has become a leading developer of seismometer-based earthquake early warning systems, which have successfully been applied to several large events. The rapid determination of earthquake magnitude and focal mechanism, important for a number of rapid response applications, including tsunami warning, is still challenging because of the limitations of near-field inertial recordings. This instrumental issue can be solved by an entirely different observation system: a GPS network. Taiwan is well posed to take advantage of GPS because in the last decade it has developed a very dense network. Thus, in this research, we explore the suitability of the GPS CMT inversion for Taiwan. We retrospectively investigate six moderate to large (〈span〉M〈/span〉〈sub〉w〈/sub〉6.0 ∼ 7.0) earthquakes and propose a resolution test for our model, we find that the minimum resolvable earthquake magnitude of this system is ∼〈span〉M〈/span〉〈sub〉w〈/sub〉5.5 (at 5 km depth). Our tests also suggest that the finite fault complexity, often challenging for the near-field methodology, can be ignored under such good station coverage and thus, can provide a fast and robust solution for large earthquake directly from the near field. Our findings help to understand and quantify how the proposed methodology could be implemented in real time and what its contributions could be to the overall earthquake monitoring system.〈/span〉

Description: 〈span〉〈div〉SUMMARY〈/div〉Time-domain processing of seismic reflection data has always been an important engine that is routinely utilized to produce seismic images and to expeditiously construct subsurface models. The conventional procedure involves analysing parameters related to the derivatives of reflection traveltime with respect to offset including normal moveout (NMO) velocities (second-order derivatives) and quartic coefficients (fourth-order derivatives). In this study, we propose to go beyond the typical assumption of 1-D laterally homogeneous medium when relating those ‘processing’ parameters to the subsurface medium parameters and take into account the additional influences from lateral heterogeneity including curved interfaces and smoothly variable velocities. We fill in the theoretical gap from previous studies and develop a general framework for such connection in layered anisotropic media. We show that in general, the influences of lateral heterogeneity get accumulated from all layers via a recursive relationship according to the Fermat’s principle and can be approximately quantified in terms of the lateral derivatives of the layer interface surfaces and velocities. Based on the same general principle, we show that our approach can also be used to study the lateral heterogeneity effects on diffraction traveltime and its second-order derivative related to time-migration velocity. In this paper, we explicitly specify expressions for NMO and time-migration velocities with the influences from both types of heterogeneity suitable for 2-D data sets and also discuss possible extensions of the proposed theory to 3-D data sets and to parameters related to higher-order traveltime derivatives. Using numerical examples, we demonstrate that the proposed theory can lead to more accurate reflection and diffraction traveltime predictions in comparison with those obtained based on the 1-D assumption. Both the proposed theoretical framework and its numerical testing for forward traveltime computation presented in this study aid in understanding the effects from lateral heterogeneity on time-processing parameters and also serve as an important basis for designing an efficient technique to separate those influences in important processes such as Dix inversion for a more accurate subsurface model in the future.〈/span〉

Description: 〈span〉〈div〉SUMMARY〈/div〉Estimating shear wave velocity with depth from Rayleigh-wave dispersion data is limited by the accuracy of fundamental and higher mode identification and characterization. In many cases, the fundamental mode signal propagates exclusively in retrograde motion, while higher modes propagate in prograde motion. It has previously been shown that differences in particle motion can be identified with multicomponent recordings and used to separate prograde from retrograde signals. Here we explore the domain of existence of prograde motion of the fundamental mode, arising from a combination of two conditions: (1) a shallow, high-impedance contrast and (2) a high Poisson ratio material. We present solutions to isolate fundamental and higher mode signals using multicomponent recordings. Previously, a time-domain polarity mute was used with limited success due to the overlap in the time domain of fundamental and higher mode signals at low frequencies. We present several new approaches to overcome this low-frequency obstacle, all of which utilize the different particle motions of retrograde and prograde signals. First, the Hilbert transform is used to phase shift one component by 90° prior to summation or subtraction of the other component. This enhances either retrograde or prograde motion and can increase the mode amplitude. Secondly, we present a new time–frequency domain polarity mute to separate retrograde and prograde signals. We demonstrate these methods with synthetic and field data to highlight the improvements to dispersion images and the resulting dispersion curve extraction.〈/span〉