Struture 3D de la croûte et du manteau supérieur des Alpes par inversion des données télésismiques de l'expérience AlpArray: méthodes et implications géodynamiques – LiSAlps

La méthode d’inversion des formes d’ondes complètes baptisée FWI (Full Waveform Inversion) s’est imposée comme l’approche d’imagerie phare depuis une décennie dans de nombreux domains comme l’exploration géophysique, la sismologie, l’imagerie médicale, le contrôle non destructif, … L’objectif au sens large est de reconstruire les propriétés constitutives d’un milieu à partir de mesures indirectes fournies par des ondes élastiques déclenchées par des sources naturelles ou artificielles et enregistrées par un dispositif de capteurs. Le modèle mathématique décrivant la propagation des ondes est une équation aux dérivés partielles, l’équation d’onde. Les propriétés du milieu sont contenues dans les coefficients de cette équation et les enregistrements au cours du temps collectés durant une campagne d’acquisition représentent les solutions de cette équation aux positions de l’espace où sont placés les capteurs. Cette estimation de paramètres est formulée sous forme d’un problème inverse résolu avec des méthodes d’optimisation locales par minimisation d’une fonction mesurant une distance entre les données enregistrées et simulées numériquement. L’attrait suscité par la FWI réside dans le fait que toute l’information contenue dans les mesures est utilisée en traitant chaque échantillon discrétisant chaque sismogramme de l’expérience. La plus value de cette approche totale est une résolution spatiale de l’ordre de la longueur d’onde et la possibilité théorique d’estimer toutes les propriétés du milieu gouvernant la propagation des ondes (coefficients élastiques, densité, atténuation). Le prix à payer est la nonlinéarité du problème inverse engendrée par la nature ondulatoire des signaux, son caractère mal posé exacerbé par l’éclairage incomplet de l’intérieur de la Terre depuis la surface et le coût calculatoire d’approches fondées sur la résolution complète de l’équation d’onde. En imagerie télésismique, le sources (les séismes) sont éloignées de la cible d’étude située sous le dispositif de stations. Cette configuration nécessite des protocoles spécifiques de modélisation des ondes où une première simulation est effectuée dans un modèle de Terre globale simplifié pour générer un champ d’onde de référence, qui est injecté aux frontières de la cible d’étude puis propagé en son sein jusqu’aux stations. Le challenge du problème inverse est d’atteindre la résolution visée en extrayant l’information contenue dans le champ d’onde incident diffracté vers l’avant mais aussi le champ d’onde retro diffracté par les hétérogénéités du milieu après une première réflexion sur la topographie. Cela soulève de multiples challenges méthodologiques liés à la simulation numérique des ondes dans des maillages de grande dimension et à la conception d’algorithmes d’optimisation robustes équipés de régularisations performantes pour traiter des mesures parcellaires et bruitées à partir desquelles sont estimées plusieurs classes de paramètres à l’empreinte contrastée.

Avant de lancer les premières inversions, un contrôle qualité des données est primordiale pour éviter d’injecter des données erronées ou inutilisables dans l’algorithme. Ces données aberrantes peuvent avoir été générées par des disfonctionnements de certaines stations ou par un rayonnement peu favorable de la source par rapport à l’azimuth d’arrivée à la station générant un rapport signal sur bruit médiocre pour l’onde P. Un premier jeu de données regroupant une vingtaine d’évènements et les 600 stations de l’expérience AlpArray a été sélectionné à partir de critères fondés sur l’analyse de la phase, de l’amplitude et de la polarisation des arrivée enregistrées par les stations multi composantes. Ce jeu de données est construit en sélectionnant une fenêtre temporelle d’une longueur à définir à partir du temps d’arrivée de l’onde P. La bande passante conservée pour cette première inversion est comprise entre 10s et 30s et le longueur de la fenêtre est limitée à 15s. La simulation des ondes dans la cible d’étude est effectuée avec la méthode des éléments spectraux sur maillage hexahédrique tandis que le calcul des champs d’onde dans la terre globale entre les séismes et les frontières de la cicle est effectuée avec le code AxiSEM. Le maillage du domaine de calcul est composé de 126 900 éléments discrétisés avec des polynomes de Lagrange d’ordre 4. Le modèle initial est construit à partir du modèle global AK135 et inclut la topographie des Alpes. L’inversion est effectuée avec un algorithme de quasi-Newton l-BFGS avec une soixantaine d’itérations. Des signatures de source initiales sont tirées du catalogue GlobalCMT et sont mises à jour à la volée lors de l’inversion de manière alternée avec la mise à jour du mileu paramétré par les vitesses de propagation P et S et la densité. Des coupes horizontales à 10km de profondeur du modèle de perturbation des vitesses P permettent de délimiter les principales structures tectoniques connues dans la croûte supérieure comme les grands bassins sédimentaires (Po, le bassin molassique, le bassin Rhone-Bresse) et le corps d‘Ivrée. Par ailleurs, une première carte du Moho a été établie révélant l’épaississement crustal sous les Alpes et les Apennins, un léger amincissement crustal sous le graben du Rhin et de la vallée du Rhone et une croûte très amincie dans le bassin Ligure. Concernant l’imagerie dans le manteau supérieur, des coupes verticales révèlent les panneaux plongeants Alpin et Apennins dont l’empreinte peut être suivie entre 75km et 350km de profondeur. Ce travail se poursuit pour fournir la prochaine génération de modèles élastiques des Alpes en utilisant plus d’évènements, en augmentant la fénêtre temporelle et en élargissant la bande passante.

Les perspectives à court terme de ce travail sont d’utiliser plus d’évènements dans l’inversion, utiliser une fenêtre d’écoute plus importante en vérifiant que cela ne rajoutte pas d’artefacts dans l’imagerie et étendre l’inversion à plus haute fréquence (0.3-Hz) pour améliorer la résolution de l’imagerie tout en vérifiant l’empreinte de l’acquisition parcellaire sur l’imagerie haute résolution. Cela concluera la première phase de travail délivrant un nouveau modèle de réference élastique isotrope de la croûte et du manteau supérieur des Alpes. Parallèlement à cela, un jeu de données contenant les ondes S sera construit pour améliorer les résultats obtenus à partir des ondes P et de sa coda en particulier la reconstruction de Vs. Ce jeu de données sera aussi utile pour des premiers tests de faisabilité visant à étendre la FWI à la reconstruction de l’anisotropie.
Une deuxième phase applicative sera focalisée sur l’imagerie de la zone d’étude centrée sur le noeud Ligure pour laquelle les données de la campagne CIFALPS2 complèteront celles d’AlpArray. Les mêmes stratégies numériques que pour la cible d’étude globale pourront être utilisées en poussant la FWI jusqu’à une fréquence de 1Hz. D’un point de vue méthodologique, nous analyserons les techniques de régularisation permettant de concilier ume imagerie haute résolution de structures contrastées comme celles que l’on peut observer dans la croûte avec l’imagerie d’un milieu plus lisse comme celui composant le manteau supérieur. Cela peut s’envisager avec des régularisations hybrides combinant des méthodes de variation totale et des régularisations classiques de Tikhonov. Ces régularisations doivent être combinées mathématiquement de manière à pouvoir adapter automatiquement le type de régularisation utilisée aux propriétés satistiques locales du milieu. Des régularisation fondées sur la norme l1 comme la variation totale sont également utiles pour minimiser l’empreinte des acquisitions parcellaires en termes d’aliasing sur l’imagerie. Pour mieux mesurer ces effets et calibrer les régularisations permettant de les minimiser, la conception d’un géomodèle des Alpes serait utile por appliquer la FWI sur un jeu de données contrôlé. Enfin, de nouvelles formulations de la FWI ont été récemment proposées pour reconstruire des champs d’onde plus précis que ceux prédit par le modèle initial. Ces champs d’onde sont calculés dans le modèle initial avec un terme de source supplémentaire de rappel aux données. On parlera de champs d’onde construits par assimilation de données. Ces champs d’ondes permettent de reconstruire des noyaux de sensibilité de la FWI plus précis lors de la mise à jour des paramètres. Cette approche est plus spécifiquement intéressante lorsque l’inversion est fortement non linéaire ce qui n’est pas le cas en sismologie en raison des basses fréquences disponibles. Néanmoins, cela ne doit pas empécher d’en évaluer l’apport potentiel à l’échelle télésismique.

En cours

Les Alpes constituent un laboratoire naturel privilégié pour étudier l'orogenèse en raison de affleurements exceptionnels. La structure profonde des Alpes et les processus tectoniques associés restent néanmoins controversés en raison de la complexité de la chaîne mais aussi de la faible résolution avec laquelle ces niveaux profonds ont été imagés. Pour palier ce manque, le consortium Européen AlpArray a déployé un réseau dense de stations large bande (inter-station < 52km) sur l'ensemble de l'arc Alpin ainsi que 30 stations sismiques sous-marines dans le bassin Ligure. Le projet LisAlps propose d'appliquer la « Full Waveform Inversion » (FWI) aux données télésismiques enregistrées durant AlpArray pour construire (1) un nouveau modèle de référence multi-paramétrique de la lithosphère et de l'asthénosphère à l'échelle de l'arc alpin (1500kmx700km) jusqu'à une profondeur de 700km à partir de l'ensemble du réseau et d'un catalogue de 300 séismes, (2) un modèle lithosphérique haute-résolution des Alpes Occidentales centré sur la zone complexe du noeud Ligure. Cette imagerie doit permettre de (a) préciser la géométrie et la pétrologie de la subduction continentale entre les plaques Européennes et Adriatiques, (b) améliorer la localisation de la sismicité dans les Alpes Occidentales, (c) prolonger à grandes profondeurs l'imagerie de la subduction continentale pour tester la continuité latérale des panneaux plongeants entre les Alpes Occidentales et les Alpes centrales, et préciser la géométrie du panneau plongeant dans la zone de transition entre les Alpes et les Apennins, (d) imager l'interaction entre les panneaux de vergence opposée des Alpes et des Apennins et comprendre le rôle de cette interaction dans la compression des Apennins du Nord, (e) comprendre la reprise en compression suivie de l'inversion structurale de la marge nord Ligure et l'origine de la microsismicité et des séismes occasionnels de forte magnitude dans le domaine nord Ligure. Répondre à ces questions soulève trois défis méthodologiques. Contrairement aux autres méthodes d'imagerie en sismologie, la FWI exploite la phase et l'amplitude de toutes les arrivées enregistrées dans une fenêtre temporelle significative. Cela confère à la FWI une résolution de l'ordre de la longueur d'onde et la capacité de reconstruire plusieurs propriétés physiques (Vp,Vs, densité, atténuation, anisotropie). L'échantillonnage souvent grossier des sources et des capteurs des acquisitions passives exposent néanmoins la FWI à de l'« aliasing ». Nous proposons de développer des méthodes de “compress sensing” via des régularisations fondées sur la norme l1 pour réduire ces artefacts. Deuxièmement, l'atténuation et l'anisotropie sont des marqueurs indirects de plusieurs propriétés importantes de l'intérieur de la Terre telles (composition, température, états de contrainte et de déformation). Bien que reconstruire ces paramètres de “second ordre” est un problème inverse mal posé, leur signature dans les données télésismiques est suffisante pour les imager par FWI. Enfin, la FWI est un problème inverse fortement non linéaire en raison du caractère oscillant des signaux sismiques. Pour étendre le régime linéaire de la FWI, un nouveau paradigme de la FWI reconstruit les champs d'onde avec une assimilation des données via une relaxation de l'équation d'onde avant d'estimer les paramètres de la Terre en minimisant les résidus de source générés par la relaxation (un problème linéaire). Cette approche est efficace si les ondes s'étant propagées dans l'ensemble du milieu à reconstruire ont une signature forte dans les données. C'est le cas des données télésismiques où les sources sont des ondes planes pénétrant la cible par ses bords avant d'interagir avec la topographie qui renvoie les ondes dans le milieu pour se réfléchir sur les discontinuités à imager. La faisabilité de ces approches a été vérifiée sur des cas synthétiques. Nous souhaitons tester leur potentiel sur les données d'AlpArray.