In Chapter 2, we consider the hydrostatic equilibrium figure of a rotating earth model with arbitrary radial density profile. We derive an exact non-linear partial differential equation describing the equilibrium figure. Perturbation theory is used to obtain approximate forms of this equation, and we show that the first-order theory is equivalent to Clairaut's equation. In Chapter 3, a method for parametrizing the possible equilibrium stress fields of a laterally heterogeneous earth model is described. In this method a solution of the equilibrium equations is first found that satisfies some desirable physical property. All other solutions can be written as the sum of this equilibrium stress field and a divergence-free stress tensor field whose boundary tractions vanish. In Chapter 4, we consider the minor vector method for the stable numerical solution of systems of linear ordinary differential equations. Results are presented for the application of the method to the calculation of seismic displacement fields in spherically symmetric, self-gravitating earth models. In Chapter 5, we present a new implementation of the direct solution method for calculating normal mode spectra in laterally heterogeneous earth models. Numerical tests are presented to demonstrate the validity and effectiveness of this method for performing large mode coupling calculations. In Chapter 6, we consider the theoretical basis for the viscoelastic normal mode method which is used in studies of seismic wave propagation, post-glacial rebound, and post-seismic deformation. We show how the time-domain solution to the viscoelastodynamic equation can be written as a normal mode sum in a rigorous manner.

Dans cette thèse, on s'intéresse à la propagation d'ondes en milieu fortement hétérogène modélisée par l'équation d'Helmholtz. Les méthodes numériques permettant de résoudre ce problème souffrent de dispersion numérique, en particulier à haute fréquence. Ce phénomène, appelé "effet de pollution", est largement analysé dans la littérature quand le milieu de propagation est homogène et l'utilisation de "méthodes d'ordre élevé" est souvent proposée pour minimiser ce problème. Dans ce travail, on s'intéresse à un milieu de propagation hétérogène, cas pour lequel on dispose de moins de connaissances. On propose d'adapter des méthodes éléments finis d'ordre élevé pour résoudre l'équation d'Helmholtz en milieu hétérogène, afin de réduire l'effet de pollution. Les méthodes d'ordre élevé étant généralement basées sur des maillages "larges", une stratégie multi-échelle originale est développée afin de prendre en compte des hétérogénéités de petite échelle. La convergence de la méthode est démontrée. En particulier, on montre que la méthode est robuste vis-a-vis de l'effet de pollution. D'autre part, on applique la méthode a plusieurs cas-tests numériques. On s'intéresse d'abord à des problèmes académiques, qui permettent de valider la théorie de convergence développée. On considère ensuite des cas-tests "industriels" appliqués à la Géophysique. Ces derniers nous permettent de conclure que la méthode multi-échelle proposée est plus performante que les éléments finis "classiques" et que des problèmes 3D réalistes peuvent être considérés. / The main objective of this work is the design of an efficient numerical strategy to solve the Helmholtz equation in highly heterogeneous media. We propose a methodology based on coarse meshes and high order polynomials together with a special quadrature scheme to take into account fine scale heterogeneities. The idea behind this choice is that high order polynomials are known to be robust with respect to the pollution effect and therefore, efficient to solve wave problems in homogeneous media. In this work, we are able to extend so-called "asymptotic error-estimate" derived for problems homogeneous media to the case of heterogeneous media. These results are of particular interest because they show that high order polynomials bring more robustness with respect to the pollution effect even if the solution is not regular, because of the fine scale heterogeneities. We propose special quadrature schemes to take int account fine scale heterogeneities. These schemes can also be seen as an approximation of the medium parameters. If we denote by h the finite-element mesh step and by e the approximation level of the medium parameters, we are able to show a convergence theorem which is explicit in terms of h, e and f, where f is the frequency. The main theoretical results are further validated through numerical experiments. 2D and 3D geophysica benchmarks have been considered. First, these experiments confirm that high-order finite-elements are more efficient to approximate the solution if they are coupled with our multiscale strategy. This is in agreement with our results about the pollution effect. Furthermore, we have carried out benchmarks in terms of computational time and memory requirements for 3D problems. We conclude that our multiscale methodology is able to greatly reduce the computational burden compared to the standard finite-element method

Les modèles géologiques numériques 2D et 3D permettent de comprendre l'organisation spatiale des roches du sous-sol. Ils sont également conçus pour réaliser des simulations numériques afin d’étudier ou de prédire le comportement physique du sous-sol. Pour résoudre les équations qui gouvernent les phénomènes physiques, les structures internes des modèles géologiques peuvent être discrétisées spatialement à l’aide de maillages. Cependant, la qualité des maillages peut être considérablement altérée à cause de l’inadéquation entre, d’une part, la géométrie et la connectivité des objets géologiques à représenter et, d’autre part, les contraintes requises sur le nombre, la forme et la taille des éléments des maillages. Dans ce cas, il est souhaitable de modifier un modèle géologique afin de pouvoir générer des maillages de bonne qualité permettant la réalisation de simulations physiques fidèles en un temps raisonnable. Dans cette thèse, j’ai développé des stratégies de réparation et de simplification de modèles géologiques 2D dans le but de faciliter la génération de maillages et la simulation de processus physiques sur ces modèles. Je propose des outils permettant de détecter les éléments des modèles qui ne respectent pas le niveau de détail et les prérequis de validité spécifiés. Je présente une méthode pour réparer et simplifier des coupes géologiques de manière locale, limitant ainsi l’extension des modifications. Cette méthode fait appel à des opérations d’édition de la géométrie et de la connectivité des entités constitutives des modèles géologiques. Deux stratégies sont ainsi explorées : modifications géométriques (élargissements locaux de l'épaisseur des couches) et modifications topologiques (suppressions de petites composantes et fusions locales de couches fines). Ces opérations d’édition produisent un modèle sur lequel il est possible de générer un maillage et de réaliser des simulations numériques plus rapidement. Cependant, la simplification des modèles géologiques conduit inévitablement à la modification des résultats des simulations numériques. Afin de comparer les avantages et les inconvénients des simplifications de modèles sur la réalisation de simulations physiques, je présente trois exemples d'application de cette méthode : (1) la simulation de la propagation d'ondes sismiques sur une coupe au sein du bassin houiller lorrain, (2) l’évaluation des effets de site liés à l'amplification des ondes sismiques dans le bassin de la basse vallée du Var, et (3) la simulation d'écoulements fluides dans un milieu poreux fracturé. Je montre ainsi (1) qu'il est possible d’utiliser les paramètres physiques des simulations, la résolution sismique par exemple, pour contraindre la magnitude des simplifications et limiter leur impact sur les simulations numériques, (2) que ma méthode de simplification de modèles permet de réduire drastiquement le temps de calcul de simulations numériques (jusqu’à un facteur 55 sur une coupe 2D dans le cas de l’étude des effets de site) tout en conservant des réponses physiques équivalentes, et (3) que les résultats de simulations numériques peuvent être modifiés en fonction de la stratégie de simplification employée (en particulier, la modification de la connectivité d’un réseau de fractures peut modifier les écoulements fluides et ainsi surestimer ou sous-estimer la quantité des ressources produites). / Numerical geological models help to understand the spatial organization of the subsurface. They are also designed to perform numerical simulations to study or predict the rocks physical behavior. The internal structures of geological models are commonly discretized using meshes to solve the physical governing equations. The quality of the meshes can be, however, considerably degraded due to the mismatch between, on the one hand, the geometry and the connectivity of the geological objects to be discretized and, on the other hand, the constraints imposed on number, shape and size of the mesh elements. As a consequence, it may be desirable to modify a geological model in order to generate good quality meshes that allow realization of reliable physical simulations in a reasonable amount of time. In this thesis, I developed strategies for repairing and simplifying 2D geological models, with the goal of easing mesh generation and simulation of physical processes on these models. I propose tools to detect model elements that do not meet the specified validity and level of detail requirements. I present a method to repair and simplify geological cross-sections locally, thus limiting the extension of modifications. This method uses operations to edit both the geometry and the connectivity of the geological model features. Two strategies are thus explored: geometric modifications (local enlargements of the layer thickness) and topological modifications (deletions of small components and local fusions of thin layers). These editing operations produce a model on which it is possible to generate a mesh and to realize numerical simulations more efficiently. But the simplifications of geological models inevitably lead to the modification of the numerical simulation results. To compare the advantages and disadvantages of model simplifications on the physical simulations, I present three applications of the method: (1) the simulation of seismic wave propagation on a cross-section within the Lorraine coal basin, (2) the site effects evaluation related to the seismic wave amplifications in the basin of the lower Var river valley, and (3) the simulation of fluid flows in a fractured porous medium. I show that (1) it is possible to use the physical simulation parameters, like the seismic resolution, to constrain the magnitude of the simplifications and to limit their impact on the numerical simulations, (2) my method of model simplification is able to drastically reduce the computation time of numerical simulations (up to a factor of 55 in the site effects case study) while preserving an equivalent physical response, and (3) the results of numerical simulations can be changed depending on the simplification strategy employed (in particular, changing the connectivity of a fracture network can lead to a modification of fluid flow paths and overestimation or underestimation of the quantity of produced resources).

Quantifying uncertainties in large-scale forward and inverse PDE simulations has emerged as a central challenge facing the field of computational science and engineering. The promise of modeling and simulation for prediction, design, and control cannot be fully realized unless uncertainties in models are rigorously quantified, since this uncertainty can potentially overwhelm the computed result. While statistical inverse problems can be solved today for smaller models with a handful of uncertain parameters, this task is computationally intractable using contemporary algorithms for complex systems characterized by large-scale simulations and high-dimensional parameter spaces. In this dissertation, I address issues regarding the theoretical formulation, numerical approximation, and algorithms for solution of infinite-dimensional Bayesian statistical inverse problems, and apply the entire framework to a problem in global seismic wave propagation. Classical (deterministic) approaches to solving inverse problems attempt to recover the “best-fit” parameters that match given observation data, as measured in a particular metric. In the statistical inverse problem, we go one step further to return not only a point estimate of the best medium properties, but also a complete statistical description of the uncertain parameters. The result is a posterior probability distribution that describes our state of knowledge after learning from the available data, and provides a complete description of parameter uncertainty. In this dissertation, a computational framework for such problems is described that wraps around the existing forward solvers, as long as they are appropriately equipped, for a given physical problem. Then a collection of tools, insights and numerical methods may be applied to solve the problem, and interrogate the resulting posterior distribution, which describes our final state of knowledge. We demonstrate the framework with numerical examples, including inference of a heterogeneous compressional wavespeed field for a problem in global seismic wave propagation with 10⁶ parameters.

La modélisation numérique de la propagation des ondes sismiques fait partie des études principales sur le calcul du mouvement sismique basées sur de différents schèmes numériques. La prise en compte du comportement nonlinéaire du sol est consideré désormais très important afin de pouvoir calculer la réponse du milieu cohéremment aux observations sous les sollications sismiques très fortes. En plus, le paramètre de pression interstitielle, qui pourrait emmener le sol aux phénomènes de liquéfaction, devient très important pour les sols saturés. Dans cette étude, dans un premier temps, la propagation des ondes sismiques est modelisée sur une composante (1C) dans les milieux linéaires et nonlinéaires en utilisant la méthode numérique des éléments spéctraux. Les rhéologies viscoélastique et nonlinéaire sont implementées par le méthode de technique des variables de mémoire et le modèle élastoplastique d’Iwan, respectivement. Ensuite, le modèle 1D - trois composantes (3C) est développé et une comparaison préalable sur l’effet de la considération des approches 1C et 3C est faite. L’effet de pression interstitielle est implementé dans le code 1D-3C et le site américain Wildlife Refuge Liquefaction Array (WRLA), qui a été frappé par le séisme de Superstition Hills en 1987 est étudié. Le changement de la réponse du sol sous les différents hypothèses de rhéologie du sol et de mouvement d’entrée est étudié. Le mouvement calculé est noté d’être amplifié pour les basses fréquences et atténué pour les hautes fréquences en raison de l’excès de pression interstitielle dans les sols liquéfiables. Par ailleurs, le sol devient plus nonlinéaire sous le chargement triaxial dans l’approche 3C et plus dilatant dû à la nonlinéarite élevée. En dépit de la similitude entre les accélérations et les vitesses en surface des approches 1C et 3C, une importante différence dans le déplacement en surface entre les deux approches est notée. Les analyses sont répétées pour deux sites japonais Kushiro Port et Onahama Port, qui ont été influencés par le séisme de Kushiro-Oki en 1993 et le séisme de la côte Pacifique de Tohoku en 2011, respectivement. Il a été montré que les changements apportés par la nonlinéarite ne sont pas identiques dans toute la gamme de fréquence concernée et l’influence du comportement des sols non-cohésives sur la propagation des ondes sismiques dépend fortement des propriétés du modèle et des conditions de chargement. Dernièrement, le code SEM est avancé en 2D en considerant les mêmes modèles implementés en 1D-3C pour la nonlinéarite du sol et les effets de pression interstitielle. Le code SEM 2D est mis en application dans un modèle de bassin sédimentaire dont la géometrie est assymmétrique et le profile du sol est composé des couches possédant différentes propriétés nonlinéaires. Le modèle est étudié par les analyses totale et effective pour les propagation des ondes P-SV et SH. La differentiation du mouvement calculé en surface est très importante sous les chargements avec les signaux d’entrée synthétique et réel. L’analyse effective résulte en plus de déformations dans les couches superficielles par rapport à l’analyse totale.De plus, la durée de propagation des ondes augmente à l’intérieur du bassin et les reflections aux frontières de bassin-rocher entraînent plus de nonlinéarite dans les coins du bassin. Cette thèse révèle la possibilité de la modélisation du comportement nonlinéaire du sol en prenant en compte l’effet de pression interstitielle dans les études de la propagation des ondes sismiques en couplant les modèles différents avec la méthode des éléments spéctreaux. Ces analyses contribuent à l’identification et la compréhension des phénomènes majeures qui se déroulent dans les couches superficielles en respectant les conditions locales et les mouvements d’entrées, ce quirend ce travail très important pour les études spécifiques de sites / Numerical modeling of seismic wave propagation has been a major topic on ground motion studies using a number of different numerical integration schemes. The consideration of soil nonlinearity holds an important place in order to achieve simulations consistent with real observations for strong seismic shaking. Additionally, in the presence of strong ground motion in saturated soils, pore pressure becomes an important parameter to take into account for related phenomena such as flow liquefaction and cyclic mobility. In this study, first, one component (1C) - seismic wave propagation is modeled in linear and nonlinear media in 1D based on the spectral element numerical method. Viscoelastic and nonlinear soil rheologies are implemented by use of the memory variables technique and Iwan’s elastoplastic model, respectively. Then, the same study is extended to a 1D - three component (3C) model and a preliminary comparison on the effect of using 1C and 3C approaches is made. Then, the influence of excess pore pressure development is included in the 1D-3C model and the developped numerical model is applied to realistic case on the site of Wildlife Refuge Liquefaction Array (USA) which is affected by the 1987 Superstition Hills event. The ground motion modification for different assumptions of the soil rheology in the media and different input motions is studied. The calculated motion is found to be amplified on low frequency and damped in high frequency range due to excess pore pressure development. Furthermore, the soil is found to be more nonlinear under triaxial loading in 3C approach and more dilative due to higher nonlinearity. Despite the similitude in surface acceleration and velocity results, significant differences in surface displacement results of 1C and 3C approaches are remarked. Similar analyses are performed on two Japanese sites Kushiro Port and Onahama Port, which are influenced by the 1993 Kushiro-Oki and the 2011 off the Pacific coast of Tohoku earthquakes, respectively. It has been shown that the nonlinearity-related changes are not homogeneous all over the concerned frequency band and the influence of cohesionless soil behavior on wave propagation is highly dependent on model properties and loadingconditions. Lastly, the 2D SEM code is developped by taking into account soil nonlinearity and pore pressure effects similary to 1D-3C SEM code. The developped 2D SEM code is applied to a 2D sedimentary basin site where the basin geometry is asymmetrical and soil profile consists of layers with different nonlinearity properties. Total and effective stress analyses are performed on the 2D basin for P-SV and SH zave propagation models. The calculated surface motion is shown to differ significantly under synthetic and realistic input motion loading conditions and the resultant deformation in superficial layers is found to be very high in effective stress analysis compared to total stress analysis. Also, wave propagation takes longer time inside basin media and the reflections on bedrock-basin boundaries lead the soil in basin edges to higher nonlinearity. This study shows the possibility of modeling nonlinear soil behavior including pore pressure effects in seismic wave propagation studies by coupling different models with spectral element method. These analyses help identifying and understanding dominant phenomena occurring in superficial layers, depending on local conditions and input motions. This is of great importance for site-specific studies