SEM-DG Approximation for elasto-acoustics
Hélène Barucq, Henri Calandra, Aurélien Citrain, Julien Diaz, Christian Gout
To cite this version:
Hélène Barucq, Henri Calandra, Aurélien Citrain, Julien Diaz, Christian Gout. SEM-DG Approxi-
mation for elasto-acoustics. MATHIAS 2018 Computational Science Engineering & Data Science by
TOTAL, Oct 2018, Serris, France. �hal-01872812v2�
HAL Id: hal-01872812
[Link]
Submitted on 29 Oct 2018
HAL is a multi-disciplinary open access L’archive ouverte pluridisciplinaire HAL, est
archive for the deposit and dissemination of sci- destinée au dépôt et à la diffusion de documents
entific research documents, whether they are pub- scientifiques de niveau recherche, publiés ou non,
lished or not. The documents may come from émanant des établissements d’enseignement et de
teaching and research institutions in France or recherche français ou étrangers, des laboratoires
abroad, or from public or private research centers. publics ou privés.
SEM-DG Approximation for elasto-acoustics
Hélène Barucq1 , Henri Calandra2 , Aurélien Citrain3,1 , Julien Diaz1 and Christian Gout3
1 Team project Magique.3D, INRIA, E2S UPPA, CNRS, Pau, France.
2 TOTAL SA, CSTJF, Pau, France.
3 INSA Rouen-Normandie Université, LMI EA 3226, 76000, Rouen.
MATHIAS 2018 October 22-24
The authors thank the M2NUM project which is co-financed by the European Union with the European regional development fund (ERDF, HN0002137) and
by the Normandie Regional Council.
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 1 / 24
Why using hybrid meshes?
Useful when the use of unstructured grid is non-sense (e.g. medium with a layer of water).
Well suited for the coupling of numerical methods in order to reduce the computational cost
and improve the accuracy.
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 2 / 24
Why using hybrid meshes?
water
water sand
sandstone
salt
Useful when the use of unstructured grid is non-sense (e.g. medium with a layer of water).
Well suited for the coupling of numerical methods in order to reduce the computational cost
and improve the accuracy.
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 2 / 24
Elastodynamic system
x ∈ Ω ⊂ Rd , t ∈ [0, T ], T > 0 :
∂v
ρ(x) (x, t) = ∇ · σ(x, t),
∂t
∂σ (x, t)
= C (x)(v (x, t)).
∂t
With:
ρ(x) the density,
C (x) the elasticity tensor,
(x, t) the deformation tensor,
v (x, t), the wavespeed,
σ(x, t) the strain tensor.
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 3 / 24
Elasticus software
Software written in Fortran for wave propagation simulation in the time domain
Features
Simulation:
on various types of meshes (unstructured triangles and tetrahedra),
on heterogeneous media (acoustic, elastic and elasto-acoustic).
Discontinuous Galerkin (DG) based on unstructured triangles and unstructured tetrahedra,
with various time-schemes : Runge-Kutta (2 or 4), Leap-Frog,
with multi-order computation(p-adaptivity)...
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 4 / 24
Table of contents
1 Numerical Methods
2 Comparison DG/SEM on structured quadrangle mesh
3 DG/SEM coupling
4 3D extension
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 5 / 24
1 Numerical Methods
Discontinous Galerkin Method (DG)
Spectral Element Method (SEM)
Advantages of each method
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 6 / 24
Discontinuous Galerkin Method
Use discontinuous functions :
Degrees of freedom necessary on each cell :
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 7 / 24
Spectral Element Method
General principle
Finite Element Method (FEM) discretization + Gauss-Lobatto quadrature,
Gauss-Lobatto points as degrees of freedom (gives us exponential convergence on L2 -norm).
Z N+1
X
f (x)dx ≈ ωj f (ξj ),
j=1
ϕi (ξj ) = δij .
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 8 / 24
Advantages of each method
DG
Element per element computation ( hp-adaptivity).
Time discretization quasi explicit (block diagonal mass matrix).
Simple to parallelize.
Robust to brutal changes of physics and geometry
SEM
Couples the flexibility of FEM with the accuracy of the pseudo-spectral method.
Simplifies the mass and stiffness matrices (mass matrix diagonal).
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 9 / 24
2 Comparison DG/SEM on structured quadrangle mesh
Description of the test cases
Comparative tables
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 10 / 24
Description of the test cases
Physical parameters General context
Acoustic homogeneous medium.
Four different meshes : 10000
cells, 22500 cells, 90000 cells,
250000 cells.
CFL computed using power
iteration method.
Leap-Frog time scheme.
Eight threads parallel execution
with OpenMP.
P wavespeed 1000 m.s −1
Density 1 kg .m−3
Second order Ricker Source in Pwave
(fpeak = 10Hz)
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 11 / 24
Comparative tables
Error computed as the difference between an analytical and a numerical solution for each
method.
Three cases considered : DG without penalization terms, DG with penalization terms and
SEM.
CFL L2-error CPU-time Nb of time steps
DG(α = 0) 3.18e-3 2e-1 5.13 629
SEM 4.9e-3 5e-2 0.80 409
Figure: DG not penalized and SEM comparison on the 10000 cells case
CFL L2-error CPU-time(s) Nb of time steps
DG(α = 0) 2.12e-3 7e-1 18.11 943
SEM 3.26e-3 4e-2 3.54 613
Figure: DG not penalized and SEM comparison on the 20000 cells case
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 12 / 24
Comparative tables
Error computed as the difference between an analytical and a numerical solution for each
method.
Three cases considered : DG without penalization terms, DG with penalization terms and
SEM.
CFL L2-error CPU-time Nb of time steps
DG(α = 0) 3.18e-3 2e-1 5.13 629
SEM 4.9e-3 5e-2 0.80 409
Figure: DG not penalized and SEM comparison on the 10000 cells case
CFL L2-error CPU-time(s) Nb of time steps
DG(α = 0) 2.12e-3 7e-1 18.11 943
SEM 3.26e-3 4e-2 3.54 613
Figure: DG not penalized and SEM comparison on the 20000 cells case
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 12 / 24
Comparative tables
Error computed as the difference between an analytical and a numerical solution for each
method.
Three cases considered : DG without penalization terms, DG with penalization terms and
SEM.
CFL L2-error CPU-time Nb of time steps
DG(α = 0) 3.18e-3 2e-1 5.13 629
SEM 4.9e-3 5e-2 0.80 409
Figure: DG not penalized and SEM comparison on the 10000 cells case
CFL L2-error CPU-time(s) Nb of time steps
DG(α = 0) 2.12e-3 7e-1 18.11 943
SEM 3.26e-3 4e-2 3.54 613
Figure: DG not penalized and SEM comparison on the 20000 cells case
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 12 / 24
Comparative tables
Error computed as the difference between an analytical and a numerical solution for each
method.
Three cases considered : DG without penalization terms, DG with penalization terms and
SEM.
CFL L2-error CPU-time Nb of time steps
DG(α = 0.5) 2e-3 3e-2 7.93 1000
SEM 4.9e-3 5e-2 0.80 409
Figure: DG penalized and SEM comparison on the 10000 cells case
CFL L2-error CPU-time(s) Nb of time steps
DG(α = 0.5) 1.33e-3 2e-2 32.98 1502
SEM 3.26e-3 4e-2 3.54 613
Figure: DG penalized and SEM comparison on the 20000 cells case
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 13 / 24
Comparative tables
Error computed as the difference between an analytical and a numerical solution for each
method.
Three cases considered : DG without penalization terms, DG with penalization terms and
SEM.
CFL L2-error CPU-time Nb of time steps
DG(α = 0.5) 2e-3 3e-2 7.93 1000
SEM 4.9e-3 5e-2 0.80 409
Figure: DG penalized and SEM comparison on the 10000 cells case
CFL L2-error CPU-time(s) Nb of time steps
DG(α = 0.5) 1.33e-3 2e-2 32.98 1502
SEM 3.26e-3 4e-2 3.54 613
Figure: DG penalized and SEM comparison on the 20000 cells case
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 13 / 24
Comparative tables
Error computed as the difference between an analytical and a numerical solution for each
method.
Three cases considered : DG without penalization terms, DG with penalization terms and
SEM.
CFL L2-error CPU-time Nb of time steps
DG(α = 0.5) 2e-3 3e-2 7.93 1000
SEM 4.9e-3 5e-2 0.80 409
Figure: DG penalized and SEM comparison on the 10000 cells case
CFL L2-error CPU-time(s) Nb of time steps
DG(α = 0.5) 1.33e-3 2e-2 32.98 1502
SEM 3.26e-3 4e-2 3.54 613
Figure: DG penalized and SEM comparison on the 20000 cells case
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 13 / 24
Comparative tables
Error computed as the difference between an analytical and a numerical solution for each
method.
Three cases considered : DG without penalization terms, DG with penalization terms and
SEM.
CFL L2-error CPU-time Nb of time steps
DG(α = 0.5) 2e-3 3e-2 7.93 1000
SEM 4.9e-3 5e-2 0.80 409
Figure: DG penalized and SEM comparison on the 10000 cells case
CFL L2-error CPU-time(s) Nb of time steps
DG(α = 0.5) 1.33e-3 2e-2 32.98 1502
SEM 3.26e-3 4e-2 3.54 613
Figure: DG penalized and SEM comparison on the 20000 cells case
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 13 / 24
Comparative tables
Error computed as the difference between an analytical and a numerical solution for each
method.
Three cases considered : DG without penalization terms, DG with penalization terms and
SEM.
CFL L2-error CPU-time(s) Nb of time steps
DG(α = 0.5) 2e-3 3e-2 7.93 1000
SEM 2e-3 3e-2 2.12 1000
Figure: DG penalized and SEM comparison using the same CFL on a 10000 thousands cells mesh
CFL L2-error CPU-time(s) Nb of time steps
DG (α = 0.5) 1.33e-3 2e-2 32.98 1502
SEM 1.33e-3 2e-2 8.67 1502
Figure: DG penalized and SEM comparison using the same CFL on a 20000 thousands cells mesh
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 14 / 24
Comparative tables
Error computed as the difference between an analytical and a numerical solution for each
method.
Three cases considered : DG without penalization terms, DG with penalization terms and
SEM.
CFL L2-error CPU-time(s) Nb of time steps
DG(α = 0.5) 2e-3 3e-2 7.93 1000
SEM 2e-3 3e-2 2.12 1000
Figure: DG penalized and SEM comparison using the same CFL on a 10000 thousands cells mesh
CFL L2-error CPU-time(s) Nb of time steps
DG (α = 0.5) 1.33e-3 2e-2 32.98 1502
SEM 1.33e-3 2e-2 8.67 1502
Figure: DG penalized and SEM comparison using the same CFL on a 20000 thousands cells mesh
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 14 / 24
Comparative tables
Error computed as the difference between an analytical and a numerical solution for each
method.
Three cases considered : DG without penalization terms, DG with penalization terms and
SEM.
CFL L2-error CPU-time(s) Nb of time steps
DG(α = 0.5) 2e-3 3e-2 7.93 1000
SEM 2e-3 3e-2 2.12 1000
Figure: DG penalized and SEM comparison using the same CFL on a 10000 thousands cells mesh
CFL L2-error CPU-time(s) Nb of time steps
DG (α = 0.5) 1.33e-3 2e-2 32.98 1502
SEM 1.33e-3 2e-2 8.67 1502
Figure: DG penalized and SEM comparison using the same CFL on a 20000 thousands cells mesh
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 14 / 24
Advantages of each method
DG
Element per element computation ( hp-adaptivity).
Time discretization quasi explicit (block diagonal mass matrix).
Simple to parallelize.
Robust to brutal changes of physics and geometry
SEM
Couples the flexibility of FEM with the accuracy of the pseudo-spectral method.
Simplifies the mass and stiffness matrices (mass matrix diagonal).
Reduces the computational costs on structured quadrangle cells in comparison with DG
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 15 / 24
3 DG/SEM coupling
Hybrid meshes structures
Variational formulation
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 16 / 24
Hybrid meshes structures
Aim at coupling Pk and Qk structures.
Need to extend or split some structures (e.g. neighbour indexes).
Define new face matrices:
Z Z Z
MijK ,L = φK L K ,L
i φj , Mij = ψiK ψjL , MijK ,L = φK L
i ψj .
K ∩L K ∩L K ∩L
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 17 / 24
Variational formulation
Global context
Domain in two parts : Ωh,1 (structured quadrangles + SEM), Ωh,2 (unstructured triangles
+ DG).
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 18 / 24
Variational formulation
SEM variational formulation :
Z Z Z
ρ∂t v1 · w1 = − σ1 · ∇w1 + (σ1 n1 ) · w1 ,
Ωh,1
Ωh,1 Γout,1
Z Z Z
∂t σ1 : ξ1 = − (∇(C ξ1 )) · v1 + (C ξ1 n1 ) · v1 .
Ωh,1 Ωh,1 Γout,1
DG variational formulation :
Z Z Z Z
ρ∂t v2 · w2 = − σ2 · ∇w2 + (σ2 n2 ) · w2 + {{σ2 }}[[w2 ]] · n2 ,
Ωh,2
Ωh,2 Γout,2 Γint
Z Z Z Z
∂t σ2 : ξ2 =− (∇(C ξ2 )) · v2 + (C ξ2 n2 ) · v2 + {{v2 }}[[C ξ2 ]] · n2 .
Ωh,2 Ωh,2 Γout,2 Γint
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 18 / 24
Variational formulation
Z Z Z Z
ρ∂ t v1 · w 1 + ρ∂ t v2 · w 2 = − σ 1 · ∇w 1 − σ2 · ∇w2
Ωh,1 Ωh,2 Ωh,1 Ωh,2
Z Z Z
+ (σ1 n1 ) · w1 + (σ2 n2 ) · w2 + {{σ2 }}[[w2 ]] · n2
Γout,1
Γout,2 Γint
Z
[[σw ]] · n,
+
Γ1/2
Z
Z Z Z
∂ σ : ξ + ∂ σ : ξ = − (∇(C ξ )) · v − (∇(C ξ2 )) · v2
t 1 1 t 2 2 1 1
Ωh,1 Ωh,2 Ωh,1 Ωh,2
Z Z Z
(C ξ1 n1 ) · v1 + (C ξ2 n2 ) · v2 + {{σ2 }}[[w2 ]] · n2
+
Γout,1
Γout,2 Γint
Z
+
[[(C ξ)v ]] · n.
Γ1/2
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 19 / 24
Variational formulation
Z Z Z Z
ρ∂ t v1 · w 1 + ρ∂ t v2 · w 2 = − σ 1 · ∇w 1 − σ2 · ∇w2
Ωh,1 Ωh,2 Ωh,1 Ωh,2
Z Z Z
+ (σ1 n1 ) · w1 + (σ2 n2 ) · w2 + {{σ2 }}[[w2 ]] · n2
Γout,1
Γout,2 Γint
Z
{{σ}}[[w ]] · n + [[σ]]{{w }} · n,
+
Γ1/2
Z
Z Z Z
∂ σ : ξ + ∂ σ : ξ = − (∇(C ξ )) · v − (∇(C ξ2 )) · v2
t 1 1 t 2 2 1 1
Ωh,1 Ωh,2 Ωh,1 Ωh,2
Z Z Z
(C ξ1 n1 ) · v1 + (C ξ2 n2 ) · v2 + {{σ2 }}[[w2 ]] · n2
+
Γout,1
Γout,2 Γint
Z
+
[[C ξ]]{{v }} · n + {{C ξ}}[[v ]] · n.
Γ1/2
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 19 / 24
Variational formulation
Z Z Z Z
ρ∂ t v1 · w 1 + ρ∂ t v2 · w 2 = − σ 1 · ∇w 1 − σ2 · ∇w2
Ωh,1 Ωh,2 Ωh,1 Ωh,2
Z Z Z
+ (σ1 n1 ) · w1 + (σ2 n2 ) · w2 + {{σ2 }}[[w2 ]] · n2
Γout,1
Γout,2 Γint
Z
+
{{σ}}[[w ]] · n ( (((}}
+[[σ]]{{w (( · n,
Γ1/2
Z
Z Z Z
∂ σ : ξ + ∂ σ : ξ = − (∇(C ξ )) · v − (∇(C ξ2 )) · v2
t 1 1 t 2 2 1 1
Ωh,1 Ωh,2 Ωh,1 Ωh,2
Z Z Z
(C ξ1 n1 ) · v1 + (C ξ2 n2 ) · v2 + {{σ2 }}[[w2 ]] · n2
+
Γout,1
Γout,2 Γint
Z
((
+
[[C ξ]]{{v }} · n ( +{{C
(( ξ}}[[v
(( ]] · n.
Γ1/2
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 19 / 24
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 20 / 24
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 20 / 24
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 20 / 24
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 20 / 24
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 20 / 24
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 20 / 24
General settings
Figure: Hexa/Tet boundary configuration
Only deal with a simple case of 3D hybrid meshes : one hexahedron has only two tetrahedra
as neighbour.
Extend SEM in 3D (basis functions...).
Require introducing a new matrix which handles the rotation cases between two elements.
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 21 / 24
General settings
Figure: Hexa/Tet boundary configuration
Only deal with a simple case of 3D hybrid meshes : one hexahedron has only two tetrahedra
as neighbour.
Extend SEM in 3D (basis functions...).
Require introducing a new matrix which handles the rotation cases between two elements.
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 21 / 24
General settings
Figure: Hexa/Tet boundary configuration
Only deal with a simple case of 3D hybrid meshes : one hexahedron has only two tetrahedra
as neighbour.
Extend SEM in 3D (basis functions...).
Require introducing a new matrix which handles the rotation cases between two elements.
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 21 / 24
General settings
Figure: Hexa/Tet boundary configuration
Only deal with a simple case of 3D hybrid meshes : one hexahedron has only two tetrahedra
as neighbour.
Require introducing a new matrix which handles the rotation cases between two elements.
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 21 / 24
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 22 / 24
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 22 / 24
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 22 / 24
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 22 / 24
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 22 / 24
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 22 / 24
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 22 / 24
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 22 / 24
Conclusion and perspectives
Conclusion
1 SEM is more efficient on structured quadrangle mesh than DG
2 Build a variational formulation for DG/SEM coupling and find a CFL condition that ensures
stability
3 Show the utility of using hybrid meshes and method coupling (reduce computational cost,...)
Perspectives
Implement DG/SEM coupling on the code (2D) X
Develop DG/SEM coupling in 3D X
Develop PML in the hexahedral part
Add a local time-stepping scheme
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 23 / 24
Thank you for your attention !
Questions?
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 24 / 24
Error-order graphic
100 DG
SEM
Target
L2 -error
10−1
10−2
1 2 3 4 5 6
order
Figure: L2 -error comparison on a 10000 cells mesh
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 24 / 24
Error-order graphic
100 DG
SEM
Target
L2 -error
10−1
10−2
1 2 3 4 5 6
order
2
Figure: L -error comparison on a 10000 cells mesh
CFL L2-error CPU-time Nb of time steps
DG 2e-3 3e-2 7.93 1000
SEM(order five) 2.13e-3 3e-2 9.06 943
Figure: SEM and DG comparison with fixed error
Aurélien Citrain Coupling DG/SEM MATHIAS 2018 October 22-24 24 / 24