La détermination des paramètres d'un milieu poreux, notamment de la porosité, de la perméabilité et des propriétés du fluide saturant est un enjeu important pour des problèmes hydrologiques, pétroliers ou de risques naturels. L'objectif principal de ce travail est d'estimer ces propriétés à partir des ondes sismiques réfléchies.<br />Les théories poro-élastiques (Biot,1956) nécessitent de nombreux paramètres pour caractériser les milieux poreux et conduisent à des propriétés particulières des ondes sismiques (deux ondes de compression, atténuation intrinsèquement définies,...). Ces équations sont résolues pour un milieu poreux stratifié plan saturé par un fluide homogène par une méthode de réflectivité associée à une intégration en nombre d'ondes discrets. Ce programme de simulation est tout d'abord utilisé pour estimer la sensibilité des ondes réfléchies à la localisation et à la concentration du dioxyde de carbone dans le cas d'un stockage dans un aquifère marin profond. La sensibilité de la réponse sismique aux différents paramètres du milieu poreux est ensuite établie de manière plus systématique par le calcul analytique des dérivées de Fréchet des sismogrammes et leur mise en oeuvre numérique. Les applications numériques réalisées indiquent que les paramètres primordiaux à déterminer sont la porosité et la consolidation.<br />Ces opérateurs de sensibilité ont ensuite été intégrés dans un code d'inversion de formes d'ondes complètes (algorithme de Quasi-Newton). Les calculs d'inversion réalisés à partir de données synthétiques indiquent que les distributions de porosité et les paramètres caractérisant le solide et le fluide (densité et modules mécaniques) peuvent être correctement reconstruits lorsque les autres paramètres sont bien déterminés.<br />Cependant, l'inversion de plusieurs paramètres reste un problème difficile du fait des couplages sismiques existant entre eux. Il est cependant possible de résoudre des problèmes complexes en ne considérant qu'un seul paramètre pour le fluide (saturation) et un pour les minéraux (lithologie), ou en effectuant des inversions différentielles pour suivre des variations du fluide.<br />La méthode d'inversion est finalement appliquée à un jeu de données réelles acquis sur le site côtier de Maguelonne dans l'Hérault. Les variations du milieu peuvent être reconstruites en utilisant de l'information a priori venant de forages.

La propagation des ondes sismiques dans les milieux poreux multiphasiques présente des enjeux nombreux, tant sur le plan environnemental (risques naturels, géotechnique, pollutions de nappes...) que pour les réservoirs (aquifères, hydrocarbures, stockages de CO2...). L'utilisation des ondes sismiques pour étudier ces milieux se justifie par le fait qu'en se propageant, les ondes sont déformées par le milieu qu'elles traversent et contiennent ainsi des informations aux capteurs sur les phases fluides et solides et sur le squelette poreux. Ce travail de thèse s'intéresse aux caractéristiques des ondes sismiques dans les milieux multiphasiques (plusieurs phases fluides et solides), depuis la description physique jusqu'à la caractérisation des paramètres constitutifs par inversion, en passant par la modélisation numérique 2D de la propagation. La première partie du travail a consisté à décrire la physique des milieux multiphasiques (phase par phase et leurs intéractions dynamiques) en utilisant des méthodes d'homogénéisation pour se ramener à un milieu équivalent défini par sept paramètres. Ainsi, dans des milieux simple porosité saturés et dans des milieux plus complexes (double porosité, partiellement saturés ou visco-poroélastiques), je peux calculer la propagation des ondes sismiques sans approximation. En effet, j'utilise une méthode numérique dans le domaine fréquence-espace qui permet de prendre en compte tous les termes qui dépendent de la fréquence sans approximation. La discrétisation spatiale utilise une méthode d'éléments finis discontinus (Galerkin discontinu) qui permet de considérer des milieux hétérogènes.Je montre notamment que les attributs sismiques (vitesses et atténuations) des milieux poreux complexes sont fortement dispersifs et les formes d'ondes complètes, calculées sans approximation, sont fortement dépendantes de la description physique du milieu. La caractérisation des paramètres poroélastiques s'effectue par inversion. Une méthode en deux étapes a été proposée : la première consiste en une inversion ``classique`` (tomographie, inversion des formes d'ondes complètes) des données (sismogrammes) pour obtenir des paramètres macro-échelles (attributs sismiques). La seconde étape permet de reconstruire, à partir des paramètres macro-échelles, les paramètres poroélastiques micro-échelles. Cette étape d'inversion utilise une méthode d'optimisation semi-globale (algorithme de voisinage). Une analyse de sensibilité montre qu'en connaissant a-priori certains paramètres, on peut inverser avec précision les paramètres du squelette poroélastique ou retrouver la nature du fluide saturant, à partir des vitesses de propagation. En revanche, pour retrouver la saturation en fluide, il est préférable de connaître les atténuations. Deux applications réalistes (monitoring de réservoir et hydrogéophysique) mettent en oeuvre ce type d'inversion en deux étapes et démontrent qu'à partir de données estimées par des méthodes classiques d'imagerie, on peut remonter à certains paramètres poroélastiques constitutifs. / Seismic wave propagation in multiphasic porous media have various environmental (natural risks, geotechnics, groundwater pollutions...) and ressources (aquifers, oil and gas, CO2 storage...) issues. When seismic waves are crossing a given material, they are distorted and thus contain information on fluid and solid phases. This work focuses on the characteristics of seismic waves propagating in multiphasic media, from the physical complex description to the parameter characterisation by inversion, including 2D numerical modelling of the wave propagation. The first part consists in the description of the physics of multiphasic media (each phase and their interactions), using several upscaling methods, in order to obtain an equivalent mesoscale medium defined by seven parameters. Thus, in simple porosity saturated media and in complex media (double porosity, patchy saturation, visco-poroelasticity), I can compute seismic wave propagation without any approximation. Indeed, I use a frequency-space domain for the numerical method, which allows to consider all the frequency dependent terms. The spatial discretisation employs a discontinuous finite elements method (discontinuous Galerkin), which allows to take into account complex interfaces.The computation of the seismic attributes (velocities and attenuations) of complex porous media shows strong variations in respect with the frequency. Waveforms, computed without approximation, are strongly different if we take into account the full description of the medium or an homogenisation by averages. The last part of this work deals with the poroelastic parameters characterisation by inversion. For this, I develop a two-steps method: the first one consists in a classical inversion (tomography, full waveform inversion) of seismograms data to obtain macro-scale parameters (seismic attributes). The second step allows to recover, from the macroscale parameters, the poroelastic micro-scale properties. This downscaling step uses a semi-global optimisation method (neighbourhood algorithm), which allows the sampling of the full model space (thanks to the low numerical cost of the analytic direct model). With the a-priori knowledge of some parameters, a sensibility analysis shows that I can invert precisely skeleton parameters or the saturating fluid type, from the velocities only. Nevertheless, to recover the fluid saturation, it is preferable to use the attenuations. This two-steps procedure is tested on two realistic applications (reservoir monitoring and subsurface hydrogeophysics) and show that we can recover some constituve poroelastic parameters.

L’ambition de ce travail est la prédiction du champ d’onde incident réalistique, induit par des mouvement forts de sol, aux sites d’importance stratégique, comme des centrales nucléaires. À cette fin, un plateforme multi-outil est développé et exploité pour simuler les aspects différents d’un phénomène complexe et multi-échelle comme un tremblement de terre. Ce cadre computationnel fait face à la nature diversifiée d’un tremblement de terre par approche holistique local-régionale.Un cas d’étude complexe est choisie: le tremblement de terre MW6.6 Niigata-Ken Ch¯uetsu-Oki, qui a endommagé la centrale nucléaire de Kashiwazaki-Kariwa. Les effets de site non-linéaires observés sont à premier examinés et caractérisés. Dans la suite, le modèle 3D «de la faille au site» est construit et employé pour prédire le mouvement sismique dans une bande de fréquence de 0-7 Hz. L’effet de la structure géologique pliée au-dessous du site est quantifié en simulant deux chocs d’intensité modérée et en évaluant la variabilité spatiale des spectres de réponse aux différents endroits dans le site nucléaire. Le résultat numérique souligne le besoin d’une description plus détaillée du champ d’onde incident utilisé comme paramètre d’entrée dans la conception structurel antisismique de réacteurs nucléaires et des installations. Finalement, la bande de fréquences des signaux synthétiques obtenues comme résultat des simulations numériques est agrandie en exploitant la prédiction stochastique des ordonnées spectrales à courte période fournies par des Réseaux Artificiels de Neurones. / The ambition of this work is the prediction of a synthetic yet realistic broad-band incident wave-field, induced by strong ground motion earthquakes at sites of strategic importance, such as nuclear power plants. To this end, an multi-tool platform is developed and exploited to simulate the different aspects of the complex and multi-scale phenomenon an earthquake embodies. This multi-scale computational framework copes with the manifold nature of an earthquake by a holistic local-to-regional approach. A complex case study is chosen to this end: is the MW6.6 Niigata-Ken Ch¯uetsu-Oki earthquake, which damaged the Kashiwazaki-Kariwa nuclear power plant. The observed non-linear site-effects are at first investigated and characterized. In the following, the 3D source-to-site model is constructed and employed to provide reliable input ground motion, for a frequency band of 0-7 Hz. The effect of the folded geological structure underneath the site is quantified by simulating two aftershocks of moderate intensity and by estimating the spatial variability of the response spectra at different locations within the nuclear site. The numerical outcome stresses the need for a more detailed description of the incident wave-field used as input parameter in the antiseismic structural design of nuclear reactors and facilities. Finally, the frequency band of the time-histories obtained as outcome of the numerical simulations is enlarged by exploiting the stochastic prediction of short-period response ordinates provided by Artificial Neural Networks.

Nous présentons une nouvelle méthode pour la restauration d'images sismiques. Quand on l'observe, une image sismique est le résultat d'un système de dépôt initial qui a été transformé par un ensemble de déformations géologiques successives (flexions, glissement de la faille, etc) qui se sont produites sur une grande période de temps. L'objectif de la restauration sismique consiste à inverser les déformations pour fournir une image résultante qui représente le système de dépôt géologique tel qu'il était dans un état antérieur. Classiquement, ce procédé permet de tester la cohérence des hypothèses d'interprétations formulées par les géophysiciens sur les images initiales. Dans notre contribution, nous fournissons un outil qui permet de générer rapidement des images restaurées et qui aide donc les géophysiciens à reconnaître et identifier les caractéristiques géologiques qui peuvent être très fortement modifiées et donc difficilement identifiables dans l'image observée d'origine. Cette application permet alors d'assister ces géophysiciens pour la formulation d'hypothèses d'interprétation des images sismiques. L'approche que nous introduisons est basée sur un processus de minimisation qui exprime les déformations géologiques en termes de contraintes géométriques. Nous utilisons une approche itérative de Gauss-Newton qui converge rapidement pour résoudre le système. Dans une deuxième partie de notre travail nous montrons différents résultats obtenus dans des cas concrets afin d'illustrer le processus de restauration d'image sismique sur des données réelles et de montrer comment la version restaurée peut être utilisée dans un cadre d'interprétation géologique. / We present a new method for seismic image restoration. When observed, a seismic image is the result of an initial deposit system that has been transformed by a set of successive geological deformations (folding, fault slip, etc) that occurred over a large period of time. The goal of seismic restoration consists in inverting the deformations to provide a resulting image that depicts the geological deposit system as it was in a previous state. With our contribution, providing a tool that quickly generates restored images helps the geophysicists to recognize geological features that may be too strongly altered in the observed image. The proposed approach is based on a minimization process that expresses geological deformations in terms of geometrical constraints. We use a quickly-converging Gauss-Newton approach to solve the system. We provide results to illustrate the seismic image restoration process on real data and present how the restored version can be used in a geological interpretation framework.

A partir d’un ensemble de données interprétées et issues d’une analyse préalable par un opérateur expert (horizons, failles), l’objectif de la thèse est de proposer une segmentation d’une image sismique sous-jacente en parfaite cohérence avec les lois de la géologie. L’originalité de la démarche consistera à développer des techniques de segmentation d’images sismiques, entre autres basées sur des approches de type contours actifs, contraintes par des données interprétées en supplément de propriétés intrinsèques calculées par des procédés automatiques à partir de la donnée traitée sans nécessiter une quelconque supervision contrairement aux travaux existants. Un deuxième axe consistera à ordonnancer automatiquement les horizons (surfaces) interprétés et analyser finement chaque intervalle (le lieu existant entre deux horizons), en prenant en compte son contenu (amplitude, orientation, etc.). Tout cela aboutissant à la reconstruction du pseudo-temps géologique. / The objective of the thesis is to propose a segmentation of an underlying seismic image in perfect coherence with the results of a preliminary analysis by an expert (horizons, faults). laws of geology. The originality of the approach will be to develop techniques for segmenting seismic images, among others based on active contour type approaches, constrained by data interpreted in addition to intrinsic properties calculated by automatic processes from the data processed without requiring any supervision in contrast to existing work. A second axis will be to automatically schedule the horizons (surfaces) interpreted and to analyze each interval (the place between two horizons) finely, taking into account its content (amplitude, orientation, etc.). All this resulted in the reconstruction of the geological pseudo-time.

Nous modelisons les processus de glissement instable, en tenant compte de la geometrie souvent complexe des systemes de failles et des lois de frottement observees en laboratoire. Nous assimilons les failles a des interfaces de discontinuite dans un milieu elastique parfait, la croute, et utilisons une loi de frottement dependant du glissement deduite d'experiences qui montrent que la resistance du materiau diminue du seuil statique au seuil dynamique, proportionnellement au glissement. <br /><br />Nous proposons deux methodes numeriques. La premiere a pour but de simuler l'evolution temporelle spontanee d'un reseau de faille, soumis initialement a un champ de contraintes donne et auquel nous appliquons une perturbation initiale en vitesse. Elle utilise un schema de type Newmark en temps, et une discretisation spatiale en elements finis avec decomposition de domaine. Elle se revele capable de capturer efficacement les instabilites du glissement, et en particulier la phase d'initiation, qui precede la propagation de la rupture dynamique et qui se caracterise par une forme auto-similaire et une croissance exponentielle du glissement au cours du temps. Des experiences numeriques montrent que l'interaction, sur des segments de faille presentant un recouvrement significatif, se manifeste par l'existence de « zones d'ombre » dans lesquelles les contraintes sont dechargees et le glissement inhibe. En cas de recouvrement important de deux segments de faille, on observe une dissymetrie des profils de glissement, correspondant a la disparition de la singularite de contraintes a l'une des pointes de faille.<br /><br />Le deuxieme schema numerique realise l'analyse spectrale non lineaire du probleme de l'initiation « pseudo-linearise » autour de la position d'équilibre dans laquelle le système de failles est initialement au seuil de resistance statique. La non-linearite du probleme provient de la prise en compte des zones d'ombre dont on ne connait pas la geometrie a priori. Cette analyse permet de trouver le mode qui porte la signature de l'initiation, c'est-a-dire la forme auto-similaire remarquee plus haut. La version statique de cette analyse modale fournit un critere de stabilite des reseaux de failles, c'est-a-dire la valeur limite du taux d'affaiblissement au-dela de laquelle un episode de glissement donnera lieu a un evenement sismique.<br /><br />Nous faisons l'hypothese que le mode non lineaire statique, qui caracterise un comportement en affaiblissement a la limite de la stabilite, peut etre utilise pour decrire le glissement cumule a l'echelle tectonique sur un reseau de failles normales particulier en Afar, dont nous connaissons les glissements mesures en surface. Nous montrons qu'un choix judicieux du profil d'affaiblissement « equivalent » a l'echelle tectonique permet un bon accord entre le glissement observe et le mode statique. Nous en tirons des conclusions en termes d'interaction, mais aussi de propagation et/ou branchement des segments de failles. Enfin, nous decrivons deux autres applications : <br /><br />1) l'influence de la fracturation secondaire (endommagement) sur la forme du glissement<br />2) les parametres geometriques favorables a l'apparition d'une zone de relai entre deux segments de faille se propageant l'un vers l'autre.

La simulation de la propagation d'ondes pour des configurations 3D est un domaine de recherche très actif. Le principal avantage de la BEM est de ne discrétiser que les frontières du domaine. Elle est ainsi bien adaptée aux domaines infinis. Cependant, la BEM classique conduit à des matrices pleines et donc à des coûts de calcul et mémoire importants.<br />La FMM a permis d'augmenter de manière significative les capacités de la BEM dans beaucoup de domaines d'application. <br />Dans ce travail, la FMM est étendue aux équations de l'élastodynamique 3D dans le domaine fréquentiel, pour des domaines homogènes puis, grâce à une stratégie de couplage BE-BE, aux problèmes multi-domaines. D'autres améliorations de la méthode sont aussi présentées: préconditionnement, réduction du nombre de moments, développement multipôle pour les fonctions de Green du demi-espace. Des applications en sismologie sont présentées pour des modèles canoniques ainsi qu'au modèle de la vallée de Grenoble.

Ce travail apporte de nouvelles contraintes sur les interactions entre déformation et processus d'hydratation et de percolation de magma ou de fluides dans le manteau à l'aplomb d'une zone de subduction et leurs implications sur les propriétés sismiques dans le coin mantellique. Il se base sur l'analyse de péridotites à spinelle provenant du massif de Ronda (Espagne) et deux séries de xénolites issues de zones de subduction actives (Kamchatka, Papouasie-Nouvelle-Guinée). L'étude structurale, pétrologique et géochimique de ces échantillons montrent qu'ils ont subit une percolation réactive de magma ou de fluide synchrone d'une déformation de haute température basse contrainte cohérente avec les condition PT de la base de la lithosphère ou de l'asthénosphère. Cette percolation réactive est responsable d'un enrichissement en pyroxènes, qui est souvent localisé dans des bandes parallèles à la foliation. Cet enrichissement est associé à une dispersion de l'orientation cristallographique de l'olivine. Les systèmes de glissement dominant dans l'olivine sont dans tous les cas {0kl}[100], ce qui implique que la direction de polarisation rapide des ondes S dans la partie supérieure du coin mantellique est parallèle à la direction d'écoulement dans le manteau. La décroissance de l'intensité des OPR de l'olivine associée à l'enrichissement en pyroxènes a pour conséquence une réduction non négligeable de l'anisotropie sismique qui peut induire jusqu'à 33% d'erreur sur l'interprétation de l'épaisseur de la couche anisotrope. Un enrichissement en orthopyroxène peut aussi entraîner une baisse du rapport Vp/Vs, mais ne peut expliquer les rapports Vp/Vs <1,7 cartés dans certains avant-arcs. Cependant de telles valeurs peuvent être expliquées si l'anisotropie des péridotites du coin mantellique est prise en compte. Les analyses par spectroscopie infra-rouge montrent que les olivines des deux séries de xénolites étudiées contiennent moins d'eau que la saturation théorique calculée pour leurs températures d'équilibre dans le domaine du spinelle. Ces faibles concentrations en eau sont similaires à celles observées dans les olivines des péridotites à spinelle des autres zones de subduction. Elles enregistrent probablement la faible solubilité de l'eau dans l'olivine à des pressions relativement basses et la déshydratation au cours de l'exhumation des xénolites. Ces teneurs en eau mesurées ainsi que les estimations de la saturation théorique ne sont pas suffisantes pour changer la direction de glissement dominante dans l'olivine de [100] à [001].

L'optimisation intervient dans de nombreuses applications IFPEN, notamment dans l'estimation de paramètres de modèles numériques à partir de données en géosciences ou en calibration des moteurs. Dans ces applications, on cherche à minimiser une fonction complexe, coûteuse à estimer, et dont les dérivées ne sont pas toujours disponibles. A ces difficultés s'ajoutent la prise en compte de contraintes non linéaires et parfois l'aspect multi-objectifs. Au cours de cette thèse, nous avons développé la méthode SQA (Sequential Quadradic Approximation), une extension de la méthode d'optimisation sans dérivées de M.J.D. Powell pour la prise en compte de contraintes à dérivées connues ou non. Cette méthode est basée sur la résolution de problèmes d'optimisation simplifiés basés sur des modèles quadratiques interpolant la fonction et les contraintes sans dérivées, construits à partir d'un nombre limité d'évaluations de celles-ci. Si la résolution de ce sous-problème ne permet pas une progression pour l'optimisation originale, de nouvelles simulations sont réalisées pour tenter d'améliorer les modèles. Les résultats de SQA sur différents benchmarks montrent son efficacité pour l'optimisation sans dérivées sous contraintes. Enfin, SQA a été appliqué avec succès à deux applications industrielles en ingénierie de réservoir et en calibration des moteurs. Une autre problématique majeure en optimisation étudiée dans cette thèse est la minimisation multi-objectifs sous contraintes. La méthode évolutionnaire Multi-Objective Covariance Matrix Adaptation, adaptée à la prise en compte des contraintes, s'est révélée très performante dans l'obtention de compromis pour la calibration des moteurs.

Le but de ce travail est une étude pluridisciplinaire de modélisation sismique 2D et d'inversion d'une coupe géologique de la bordure méridionale du Vercors. L'originalité de ce travail réside dans l'emploi de techniques qui mènent à la caractérisation sismique de réservoir à partir d'une coupe géologique. Le modèle géologique est construit à partir de l'étude sur le terrain des ameurements de la bordure de la plateforme urgonienne (Barrémien - Aptien inférieur). Les corrélations cartographiques et séquentielles de coupes lithologiques de détail permettent d'établir une coupe géologique de 16 kilomètres de long pour 800 mètres d'épaisseur environ. Les calcaires hémipélagiques du Barrémien inférieur et les calcaires bioclastiques du prisme de bas niveau et de l'intervalle transgressif de la séquence BA3 montrent une diminution d'épaisseur en direction du Nord. Cette géométrie résulte d'une érosion des couches du Barrémien inférieur et~Ùne disposition en biseau d'aggradation ("onlaps") des paraséquences du prisme de bas niveau et de l'intervalle transgressif de la séquence BA3 sur la limite SbB3. Ces calcaires bioclastiques ont pu être datés de la base du Barrémien supérieur. Les Calcaires urgoniens montrent une position rétrogradante par rapport aux calcaires bioclastiques de Glandasse, indiquant qu'ils appartiennent à une formation transgressive dont l'âge est Barrémien supérieur. L'échantillonnage des différents Iithofaciès identifiés permet de réaliser la mesure en laboratoire des propriétés pétrophysiques et acoustiques (Vp, Vs et densité) de ces roches. Des relations linéaires entre les vitesses et la porosité ou l'argilosité sont déterminées. Ces. relations sont utilisées pour construire les modèles d'impédance, à partir des valeurs de porosité d'un réservoir réel. La construction des modèles d'impédance est réalisée en accord avec les hypothèses de la stratigraphie séquentielle, toutes les limites des paraséquences sont des réflect,eurs potentiels. La modélisation sismique a été réalisée par deux approches différentes, tracé de rayons à déport nul et différences fmies "multioffset". Les signaux sources employés ont un contenu fréquentiel comparable à ceux généralement extraits des profils sismiques de prospection pétrolière. L'interprétation des profils sismiques est comparée à la coupe géologique. Cette comparaison montre qu'il est difficile de retrouver le découpage en séquences de dépôt "il partir de l'interprétation des relations géométriques et des terminaisons des biseaux des événements sismiques. La dernière partie de ce travail est consacrée à la mise en oeuye de deux techniques de caractérisation sismique de réservoirs sur les profils simulés. Dans un premier temps, une inversion stratigraphique a été réalisée sur les données modélisées par tracé de rayons. Les sections d'impédance acoustiques obtenues sont caractérisées par des valeurs d'impédance proches du modèle d'origine. Dans un second temps, les attributs sismiques après sommation sont employés pour réaliser une interprétation lithologique des données sismiques. L'emploi de l'impédance sismique comme attribut permet une première discrimination lithologique. Un second attribut (attribut instantané) permet de discriminer les classes lithologiques présentant la même impédance sismique. Les profils lithologiques obtenus après prédiction sont comparables au modèle géologique.