[1] A model is developed which couples fully saturated porous compaction to the viscous-plastic deformation of the skeleton matrix. The Darcy fluid flow during compaction is described by an advection-diffusion equation for the excess pressure with two source/sink terms that depend on the mechanical compressibility and viscous compaction of the pore space, the latter representing the effect of pressure solution. The incompressible deformation of the composite medium is described by a force balance equation and its rheology can be viscous, plastic, or viscoplastic (Bingham material). For the plastic and viscoplastic cases, the coupling between the compacting and plastically deforming parts of the system is through the Drucker-Prager frictional-plastic yield criterion modified by Terzaghi's principle, so that the yield strength depends on the effective dynamical pressure. The coupled system is solved using a two-dimensional (2-D) finite element method. Two problems are solved to demonstrate the behavior of our theory. The first considers compaction of a uniform sediment layer. The numerical results agree with the predictions of the nondimensional control parameters and previously published results. The second problem concerns 2-D kinematic progradation of deltaic sediments. Substratum and delta sediments have the same compaction properties and a Bingham rheology during deviatoric deformation, such that the delta undergoes linear postyield viscous flow. For certain depositional regimes, overpressure is generated. When pore pressures approach critical values, yielding occurs and the delta front fails and becomes unstable, spreading gravitationally under its own weight. The flow velocity is limited to geological rates by the Bingham viscosity. For the range of parameter values considered, pressure solution is the most effective mechanism for generating nearlithostatic fluid pressures that lead to initial failure, and it appears that mechanical compaction hardly contributes to the fluid overpressure at this stage.Citation: Morency, C., R. S. Huismans, C. Beaumont, and P. Fullsack (2007), A numerical model for coupled fluid flow and matrix deformation with applications to disequilibrium compaction and delta stability,

[1] Sudden uplift, extension, and increased igneous activity are often explained by rapid mechanical thinning of the lithospheric mantle. Two main thinning mechanisms have been proposed, convective removal of a thickened lithospheric root and delamination of the mantle lithosphere along the Moho. In the latter case, the whole mantle lithosphere peels away from the crust by the propagation of a localized shear zone and sinks into the mantle. To study this mechanism, we perform two-dimensional (2-D) numerical simulations of convection using a viscoplastic rheology with an effective viscosity depending strongly on temperature, depth, composition (crust/mantle), and stress. The simulations develop in four steps. (1) We first obtain ''classical'' sublithospheric convection for a long time period ($300 Myr), yielding a slightly heterogeneous lithospheric temperature structure.(2) At some time, in some simulations, a strong thinning of the mantle occurs progressively in a small area ($100 km wide). This process puts the asthenosphere in direct contact with the lower crust. (3) Large pieces of mantle lithosphere then quickly sink into the mantle by the horizontal propagation of a detachment level away from the ''asthenospheric conduit'' or by progressive erosion on the flanks of the delaminated area. (4) Delamination pauses or stops when the lithospheric mantle part detaches or when small-scale convection on the flanks of the delaminated area is counterbalanced by heat diffusion. We determine the parameters (crustal thicknesses, activation energies, and friction coefficients) leading to delamination initiation (step 2). We find that delamination initiates where the Moho temperature is the highest, as soon as the crust and mantle viscosities are sufficiently low. Delamination should occur on Earth when the Moho temperature exceeds $800°C. This condition can be reached by thermal relaxation in a thickened crust in orogenic setting or by corner flow lithospheric erosion in the overriding lithosphere of subduction zones.

No abstract

S U M M A R YWe present a derivation of the equations describing wave propagation in porous media based upon an averaging technique which accommodates the transition from the microscopic to the macroscopic scale. We demonstrate that the governing macroscopic equations determined by Biot remain valid for media with gradients in porosity. In such media, the well-known expression for the change in porosity, or the change in the fluid content of the pores, acquires two extra terms involving the porosity gradient. One fundamental result of Biot's theory is the prediction of a second compressional wave, often referred to as 'type II' or 'Biot's slow compressional wave', in addition to the classical fast compressional and shear waves. We present a numerical implementation of the Biot equations for 2-D problems based upon the spectralelement method (SEM) that clearly illustrates the existence of these three types of waves as well as their interactions at discontinuities. As in the elastic and acoustic cases, poroelastic wave propagation based upon the SEM involves a diagonal mass matrix, which leads to explicit time integration schemes that are well suited to simulations on parallel computers. Effects associated with physical dispersion and attenuation and frequency-dependent viscous resistance are accommodated based upon a memory variable approach. We perform various benchmarks involving poroelastic wave propagation and acoustic-poroelastic and poroelastic-poroelastic discontinuities, and we discuss the boundary conditions used to deal with these discontinuities based upon domain decomposition. We show potential applications of the method related to wave propagation in compacted sediments, as one encounters in the petroleum industry, and to detect the seismic signature of buried landmines and unexploded ordnance.

Computational physics has become an essential research and interpretation tool in many fields. Particularly, in reservoir geophysics, ultrasonic and seismic modeling in porous media is used to study the properties of rocks and characterize the seismic response of geological formations. Here, we give a brief overview of the most common numerical methods used to solve the partial differential equations describing wave propagation in fluid-saturated rocks, namely finite-difference, pseudospectral and finite-element methods, including the spectral-element technique. The modeling is based on Biot-type theories of dynamic poroelasticity, which constitute a general framework to describe the physics of wave propagation. We provide a review of the various techniques and discuss numerical implementation aspects for application to seismic modeling and rock physics, as for instance the role of the Biot diffusion wave as a loss mechanism and interface waves in porous media.

We have drawn connections between imaging in exploration seismology, adjoint methods, and emerging finite-frequency tomography. All of these techniques rely on spatial and temporal constructive interference between observed and simulated waveforms to map locations of structural anomalies. Modern numerical methods and computers have facilitated the accurate and efficient simulation of 3D acoustic, (an)elastic, and poroelastic wave propagation. Using a 2D cross section of the SEG/EAGE salt model, we have determined how such waveform simulations might be harnessed to improve onshore and offshore seismic imaging strategies and capabilities. We have found that the density sensitivity kernel in adjoint tomography is related closely to the imaging principle in exploration seismology, and that in elastic modeling the impedance kernel actually is a better diagnostic tool for reflector identification. The shear- and compressional-wave speed sensitivity kernels in adjoint tomography are related closely to finite-frequency banana-doughnut kernels, and these kernels are well suited for mapping larger-scale structure, i.e., for transmission tomography. These ideas have been substantiated by addressing problems in subsalt time-lapse migration.

No abstract

Summary SPiRaL is a joint global-scale model of wave speeds (P and S) and anisotropy (vertical transverse isotropy) variations in the crust and mantle. The model is comprised of >2.1 million nodes with five parameters at each node that capture velocity variations for P- and S-waves traveling at arbitrary directions in transversely isotropic media with a vertical symmetry axis (VTI). The crust (including ice, water, sediments, crystalline layers) is directly incorporated into the model. The default node spacing is approximately 2° in the lower mantle and 1° in the crust and upper mantle. The grid is refined with ∼0.25° minimum node spacing in highly sampled regions of the crust and upper mantle throughout North America and Eurasia. The data considered in the construction of SPiRaL includes millions of body wave travel times (crustal, regional, and teleseismic phases with multiples) and surface wave (Rayleigh and Love) dispersion. A multi-resolution inversion approach is employed to capture long-wavelength heterogeneities commonly depicted in global-scale tomography images as well as more localized details that are typically resolved in more focused regional-scale studies. Our previous work has demonstrated that such global-scale models with regional-scale detail can accurately predict both teleseismic and regional body wave travel times, which is necessary for more accurate location of small seismic events that may have limited signal at teleseismic distances. SPiRaL was constructed to predict travel times for event location and long-period waveform dispersion for seismic source inversion applications in regions without sufficiently tuned models. SPiRaL may also serve as a starting model for full-waveform inversion (FWI) with the goal of fitting waves with periods 10–50 seconds over multiple broad regions (thousands of kilometers) and potentially the globe. To gain insight to this possibility, we simulated waveforms using SPiRaL and independent waveform-based models for comparison. The performance of the travel-time-based SPiRaL model is shown to be on par with regional 3-D waveform-based models in three regions (western United States, Middle East, Korean Peninsula) suggesting SPiRaL may serve as a starting model for FWI over broad regions.