Resolving coupled physical processes in porous rocks: From linear quasi-static and dynamic phenomena to non-linear instabilities

La majorité des processus physiques sur Terre sont couplés. Un processus physique peut en induire un autre, ce qui est le cas d’une onde se propageant dans une roche poreuse saturée de fluide, ce qui induit un écoulement de fluide. Dans un milieu biphasique, l’interaction entre les phases solides et fluides conduit à des effets physiques qui ne sont pas observés dans un milieu monophasique. Ainsi, une description satisfaisante de systèmes physiques complexes nécessite un traitement particulier. L’application des théories décrivant les processus physiques couplés dans les roches poreuses fracturées est d’une grande importance dans les scénarios impliquant la séquestration géologique du CO2, l’élimination des déchets nucléaires, l’exploration et la production d’énergie géothermique et les hydrocarbures. Le développe- ment de méthodes géophysiques non invasives de détection et de surveillance de ces formations géologiques est crucial. La recherche scientifique vise à faire progresser la description quantitative et qualitative des processus physiques couplés dans les roches poreuses. Conformément à ces objectifs, les contributions présentées ici sont réparties dans quatre disciplines différentes : la micromécanique, la géophysique, la mécanique computationnelle et la poroélasticité computationnelle, et la théorie des instabilités non linéaires. Différentes méthodes analytiques et numériques sont utilisées pour résoudre la physique aux micro- et macro-échelles. Cela comprend l’étude des processus linéaires quasi-statiques et dynamiques. Ce travail de recherche contient également des résultats basés sur la théorie des processus physiques non linéaires. À l’échelle microscopique (ou à l’échelle des pores), une roche est constituée de grains, de pores et de fractures individuels. De petites déformations causées par une onde sismique se propageant à travers la roche induisent des gradients de pression dans les fractures conformes. En conséquence, un écoulement de fluide (appelé “écoulement de jet”) a lieu jusqu’à ce que la pression interstitielle s’équilibre. En raison de la viscosité du fluide et du frottement visqueux associé, un tel écoulement de fluide provoque une forte dissipation d’énergie. Des simulations numériques tridimensionnelles de l’écoulement de jet à l’aide d’une approche par éléments finis sont effectuées et les résultats sont comparés à un modèle analytique publié. A partir de cette comparaison, de nombreuses limitations de la solution analytique publiée sont quantifiées et décrites. Par la suite, un nouveau modèle analytique pour la dispersion sismique et l’atténuation associées à l’écoulement de jet est présenté pour une géométrie de pores qui a été classiquement utilisée pour expliquer ce mécanisme d’expulsion de l’eau. Ensuite, ce modèle analytique est étendu pour traiter des géométries plus complexes de l’espace poreux, beaucoup plus représentatives de celle d’une roche. Les paramètres clés de l’espace poreux qui contrôlent la fréquence caractéristique (à laquelle se produit le maximum d’atténuation) sont redéfinis. De plus, des expressions analytiques pour calculer les propriétés de rigidité effective d’un modèle de roche, dont l’espace poreux est décrit par une fracture reliée à un pore ou à plusieurs pores, sont fournies. À l’échelle macroscopique, une roche poreuse peut être décrite par un ensemble de propriétés macro- scopiques, par exemple, les modules élastiques effectifs, la perméabilité, etc. Les équations de Biot décrivent un système couplé hydro-mécaniquement et établissent la théorie largement reconnue de la poroélasticité. Le milieu biphasique est représenté par une matrice poreuse solide élastique saturée d’un fluide visqueux compressible. La réponse dynamique d’un tel milieu biphasique et isotrope se traduit par deux ondes longitudinales et une onde de cisaillement, comme prédit par Yakov Frenkel. La modélisation numérique efficace et précise des équations de la poroélasticité de Biot nécessite la connaissance des conditions exactes de stabilité. Cette recherche présente les résultats de l’analyse de stabilité de von Neumann des équations de Biot discrétisées et de l’équation d’onde amortie linéaire discrétisée. Les conditions exactes de stabilité pour un certain nombre de schémas implicites et explicites sont dérivées. De plus, un solveur numérique d’unités de traitement multi-graphiques (GPU) est développé pour résoudre les équations élastodynamiques anisotropes de Biot afin de simuler, en quelques secondes, des champs d’ondes pour des domaines spatiaux impliquant plus de 4,5 milliards de mailles. Une analyse dimensionnelle complète est présentée, réduisant ainsi le nombre de paramètres matériels nécessaires pour les expériences numériques de dix à quatre. Une analyse de dispersion en fonction de paramètres adimensionnels est effectuée, conduisant à des relations de dispersion simples et transparentes. La haute efficacité de notre implémentation numérique la rend facilement accessible pour étudier des scénarios tridimensionnels et à haute résolution d’applications pratiques. Dans le cadre de la théorie des instabilités non linéaires, une nouvelle théorie de la nucléation sismique est présentée. La rhéologie visco-plastique ou élasto-plastique la plus simple permet de modéliser la nucléation sismique spontanée. En augmentant lentement la contrainte dans le milieu, elle atteint la limite plastique, produisant ainsi la localisation de la déformation et entraînant le développement lent de bandes de cisaillement fractales. Au fil du temps, ces dernières se développent spontanément et des chutes de contrainte se produisent dans le milieu. Une chute de contrainte correspond à un nouveau modèle particulier de localisation de déformation, qui agit alors comme source sismique et déclenche la propagation des ondes sismiques. Cette nouvelle approche de modélisation est basée sur des lois de conservation sans aucune relation constitutive dérivée expérimentalement. -- The majority of the physical processes on the Earth are coupled. A physical process might induce a different one, which is the case of a wave propagating in a fluid-saturated porous rock and inducing fluid flow. In a two-phase medium, the interaction between solid and fluid phases leads to physical effects, that are not observed in a single-phase medium. Thus, a successful description of complex physical systems requires special treatment. The applications of theories describing coupled physical processes in cracked porous rocks are of great importance in scenarios involving CO2 geological sequestration, nuclear waste disposal, the exploration and production of geothermal energy, and hydrocarbons. Developing non-invasive geophysical detection and monitoring methods for these geological formations is crucial. Scientific research aims to advance the quantitative and qualitative description of coupled physical processes in porous rocks. In line with these objectives, the contributions presented here are distributed across four different disciplines: micromechanics, geophysics, computational mechanics and computational poroelasticity, and the theory of non-linear instabilities. Different analytical and numerical methods are used to resolve the physics at the micro- and macro-scales. It includes the study of linear quasi-static and dynamic processes. This research work also contains some results based on the theory of non-linear physical processes. At the micro-scale (or pore-scale), a rock consists of individual grains, pores, and cracks. Small deformations caused by a passing seismic wave propagating through the rock induce pressure gradients in compliant cracks. As a result, fluid flow (so-called squirt flow) takes place until the pore pressure equilibrates. Due to the fluid viscosity and the associated viscous friction, such fluid flow causes strong energy dissipation. Three- dimensional numerical simulations of squirt flow using a finite-element approach are performed and the results are compared against a published analytical model. From this comparison, many limitations of the published analytical solution are quantified and described. Subsequently, a new analytical model for squirt flow associated seismic dispersion and attenuation is presented for a pore geometry that has been classically used to explain squirt flow. Then, this analytical model is extended to deal with more complex geometries of the pore space, which are much more closely representative of that of a rock. The key parameters of the pore space which control the characteristic frequency (at which the maximum of attenuation occurs) are re-defined. Additionally, closed-form analytical expressions to calculate the effective stiffness properties of a rock model whose pore space is described by a crack connected to a pore or multiple pores are provided. At the macro-scale, a porous rock can be described by a set of macroscopic properties, e.g., effective elastic moduli, permeability, etc. Biot’s equations describe a hydro-mechanically coupled system and establish the widely recognized theory of poroelasticity. The two-phase medium is represented by an elastic solid porous matrix saturated with a compressible viscous fluid. The dynamic response of such an isotropic two-phase medium results in two longitudinal waves and one shear wave, as predicted by Yakov Frenkel. The efficient and accurate numerical modeling of Biot’s equations of poroelasticity requires the knowledge of the exact stability conditions. This research presents the results of the von Neumann stability analysis of the discretized Biot’s equations and the discretized linear damped wave equation. The exact stability conditions for several implicit and explicit schemes are derived. Additionally, a multi-graphical processing units (GPU) numerical solver is developed to resolve the anisotropic elastodynamic Biot’s equations to simulate, in a few seconds, wave fields for spatial domains involving more than 4.5 billion grid cells. A comprehensive dimensional analysis is presented reducing the number of material parameters needed for the numerical experiments from ten to four. A dispersion analysis as a function of dimensionless parameters is performed leading to simple and transparent dispersion relations. The high efficiency of our numerical implementation makes it readily accessible to investigate three-dimensional and high-resolution scenarios of practical applications. As a part of the theory of non-linear instabilities, a new theory for earthquake nucleation is presented. The simplest visco-plastic or elasto-plastic rheology allows us to model spontaneous earthquake nucleation. By slowly increasing the stress in the medium, it reaches the yield surface, strain localization occurs resulting in the slow development of fractal shear bands. As time evolves, shear bands grow spontaneously, and stress drops take place in the medium. A stress drop corresponds to a particular new strain localization pattern which acts as seismic source and triggers seismic wave propagation. This new modeling approach is based on conservation laws without any experimentally derived constitutive relations