La modélisation de la propagation des ondes dans les milieux complexes 3d, inhomogènes, isotropes et viscoélastiques (milieux I.I.V.), est un sujet intéressant et largement d’actualité. cette modélisation aide pour la connaître la géométrie et les propriétés physiques des milieux géologiques. La viscosité du milieu de propagation est modélisée par la prise en compte de vitesses de propagation complexes. Le phénomène de propagation est modélisé par la théorie des rais avec la prise en compte de l’hypothese de faible viscosité (hypothèse non restrictive dans le domaine de la géophysique). Lles expressions analytiques des coefficients de réflexion transmission entre deux milieux viscoélastiques ont été établies avec une étude numérique de l’effet de la faible viscosité sur ces coefficients a été menée. Le code anray est, a l’origine, un code de simulation de la propagation des ondes dans les milieux inhomogènes 3d purement élastiques. son développement repose sur la théorie des rais. Le code anray a ensuite été modifie (anrayq), afin de l’adapter a la simulation de la propagation des ondes dans les milieux I.I.V. les résultats numériques du code anrayq, ont été confrontés aux résultats numériques issus du code 5kb, code de simulation base sur la méthode de réflectivité. Enfin, des données sismiques 3d, acquises lors de la campagne hydratech (Ifremer 2002), dans une région caractérisée par la présence des hydrates de gaz et du gaz libre, ont été simulées a l’aide du code anrayq. le profil de l’atténuation viscoélastique a été aussi déterminé a partir de ces données sismiques, afin de localiser les hydrates de gaz et le gaz libre dans les couches sédimentaires.

La modélisation de la propagation d'ondes pour des géométries bidimensionnelles avec des structures complexes associée à un comportement intrinsèque non-élastique quelconque, constitue un préalable essentiel à la compréhension et à l'interprétation des données sismiques, quelle que soit leur mise en œuvre, dans les milieux naturels. La résolution des équations de l'élastodynamique par la méthode des différences finies offre une flexibilité appréciable pour la représentation de modèles complexes, mais nécessite des développements numériques très spécifiques pour 1) inclure un modèle de comportement anélastique quelconque, 2) définir des conditions absorbantes efficaces afin de simuler un domaine ouvert, 3) définir des conditions de surfaces optimales afin de permettre l'évaluation des champs en présence de forts contrastes d'impédance, particulièrement en présence d'une topographie à l'interface air-sol. Dans le cadre de cette thèse plusieurs approches algorithmiques ont été développées et testées afin de répondre aux trois impératifs précédents. La discrétisation des équations de l'élastodynamique par différences finies dans le domaine temporel a été effectuée au second ordre en temps et quatrième ordre en espace (O(2,4)). Les conditions absorbantes sont de type 'Conditional Perfectly Matched Layer' afin de maximiser leur efficacité pour un grand nombre de configurations (caractéristiques du milieu, angles d'incidence). Le maillage est proposé soit sous sa forme classique dite orthogonale, soit avec une rotation appliquée (rotated grid) pour une évaluation des champs de vitesses et de contraintes plus précise à l'interface air-sol en présence de topographie. Concernant le modèle de comportement visco-élastique, une implémentation classique de type récursive est proposée, ainsi qu'une version originale de la méthode 'Auxiliary Differential Equation' pour la généralisation à des modèles visqueux quelconques (linéaire ou non-linéaire). Les validations numériques de rigueur et les comparaisons dans divers cas canoniques pour une validation expérimentale sur des assemblages aluminium/polyméthacrylate, puis un milieu granulaire disposé dans une cuve ont permis 1) d'asseoir la validité des modèles numériques vis à vis des temps de propagation et 2) d'établir les sources d'incertitudes sur les amplitudes essentiellement dues au caractère bidimensionnel des simulations. Par ailleurs, une étude du bassin de Tunis (Tunisie) initiée par une campagne de sismique réfraction et la modélisation de formes schématiques de bassin a été abordée. Cette dernière a permis de mettre en évidence les effets de site prévisibles liés à la forme et au type de remplissage d'un bassin sédimentaire.

La connaissance des structures internes de la Terre, à différentes échelles, présente des enjeux majeurs d’ordres économiques, humains, environnementaux et scientifiques. Diverses méthodes d’imagerie ont été développées en utilisant les informations contenues dans les ondes sismiques. La méthode d’inversion des formes d’ondes construit des images quantitatives haute résolution des paramètres physiques du sous-sol, en exploitant le champ d’onde complet, sous la forme d’un problème d’optimisation. Dans ce travail de thèse, je présente l’application de l’inversion des formes d’ondes en domaine fréquentiel, pour imager les paramètres viscoélastiques dans des géométries à deux dimensions à grands offsets. Dans un premier temps, les développements méthodologiques et algorithmiques sont présentés. La modélisation de la propagation des ondes P-SV en domaine fréquentiel, le problème direct du processus d’imagerie, est assurée par une méthode d’éléments finis Galerkin discontinus, assurant une grande flexibilité dans le choix des ordres d’interpolation et dans l’utilisation de maillages triangulaires non-structurés. Le problème inverse est résolu sous une forme linéarisée, afin de limiter le nombre de simulations directes, et utilise l’algorithme quasi-Newton L-BFGS permettant de tirer bénéfice de l’estimation « économique » du Hessien. Le processus global d’imagerie est implémenté sous la forme d’un algorithme massivement parallèle destiné aux calculateurs modernes à mémoire distribuée. Dans un deuxième temps , les algorithmes développés sont appliqués à des cas d’étude. Des applications sont menées dans des modèles synthétiques réalistes représentatifs d’environnements terrestres et marins. Ces études montrent les difficultés associées à la reconstruction des paramètres élastiques à partir des données mettant en jeu des ^phénomènes de propagations complexes (ondes converties, multiples, ondes de surfaces...) Des solutions sont proposées sous forme de processus hiérarchiques multi-échelles, afin de limiter les effets des non-linéarités du problème inverse et ainsi d’améliorer la convergence du processus vers le minimum global. Enfin la sensibilité de différentes normes et critères de minimisation est analysée, à partir de données bruitées issues de modèles synthétiques réalistes, ainsi que sous l’approximation acoustique pour un jeu de données réelles pétrolière. Ces tests montrent certaines limites du formalisme classique basé sur la norme L2 dans l’espace des données, tandis que la norme L1 apparaît comme alternative robuste pour l’inversion de données décimées en domaine fréquentiel.

LA SIMULATION DE LA PROPAGATION DES ONDES ELASTIQUES DANS DES MILIEUX ANISOTROPES OU VISCOELASTIQUES PEUT ETRE EFFECTUEE PAR LES METHODES DE DIFFERENCE FINIE. EN UTILISANT DES OPERATEURS DE CONVOLUTION POUR LES DERIVEES SPATIALES, ET UNE SERIE DE TAYLOR POUR LES DERIVEES TEMPORELLES, ON ATTEINT UNE SIMULATION TRES PRECISE. CETTE PRECISION EST CONFIRMEE PAR UNE ANALYSE DE DISPERSION DES VITESSES DE PHASE DANS LE CAS ANISOTROPE ET PAR LA COMPARAISON AVEC DES SOLUTIONS ANALYTIQUES DANS LE CAS VISCOELASTIQUE. CES METHODES NUMERIQUES SONT APPLIQUEES AFIN DE SIMULER UNE EXPERIENCE DE SISMIQUE REFLECTION MARINE. CETTE EXPERIENCE MONTRE QUE L'ON NE PEUT EXPLIQUER LES AMPLITUDES ENREGISTRES LOIN DE LA SOURCE SANS RAJOUTER A L'IMPEDANCE P UN DEUXIEME PARAMETRE ELASTIQUE, NOTAMMENT LE COEFFICIENT DE POISSON. A L'AIDE D'UNE METHODE ITERATIVE QUI MINIMISE L'ECART ENTRE LES DONNEES REELLES ET LES SISMOGRAMMES SYNTHETIQUES ON OBTIENT UN MODELE DE L'IMPEDANCE P AINSI QU'UN MODELE DU COEFFICIENT DE POISSON QUI PERMETTENT UNE MEILLEURE CONTRAINTE DE LA LITHOLOGIE

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.