Nous présentons une méthode Galerkin discontinue d'ordre élévé pour la simulation de la propagation d'ondes sismiques P-SV en milieu hétérogène et en deux dimensions d'espace. Le système vitesse-contraintes du premier ordre est obtenu en supposant un milieu linéaire, isotrope et viscoélastique, prenant ainsi en compte l'atténuation intrinsèque du milieu. La relation contraintes-déformations dans le domaine temporel étant une convolution, qui nécessiterait une approximation numérique très coûteuse, nous considérons la rhéologie d'un "generalized Maxwell body" (GMB) remplaçant la convolution par un jeu d'équations différentielles. Il en résulte un système vitesse-contraintes contenant des équations supplémentaires pour les fonctions anélastiques qui traduisent l'historique de déformation du matériau. Notre méthode numérique, applicable à des maillages triangulaires non structurés, est basée sur des flux centrés et un schéma saute-mouton en temps. L'extension en espace à l'ordre élevé est obtenue grâce à des polynômes de Lagrange, définis localement dans chaque élément. La méthode étant explicite, elle ne nécessite pas d'inversion de matrice de masse globale. La méthode est validée via des simulations numériques, notamment des comparaisons avec un schéma aux différences finies.