Transcript Presentation

Discontinuous Galerkin Methods for Solving Euler
Equations
Mid Year Presentation/Status Update
December 6th, 2012
Andrey Andreyev ([email protected])
Advisor: James Baeder ([email protected])
Objectives:
Develop a Discontinuous Galerkin Method to solve the Euler Equations in
one dimension that allows for up to 3rd spatial order discretization
Develop a Discontinuous Galerkin Method to Solve the Euler Equation in
two dimensions that allows for up to 3rd order spatial discretization
Validate both codes against known solutions. Shock tube problem for the
one dimensional case and a known airfoil section for the two dimensional
case
The Euler Equations
General Form
One Dimensional Form

U   j  F(U)  0
t
 
 
U   ui 
E 
 


U  F(U)  0
t
x

 
U   u 
E
 


u j


F(U)   ui u j   ij p 
 u ( E  p) 
 j

 u 
 2

F(U)   u  p 
u ( E  p )


Spatial Discretization
Structured Mesh
Unstructured Mesh
Picture from: http://www.cglerlangen.com/downloads/Manual/ch09s16s01.html
Picture from:
http://ta.twi.tudelft.nl/users/wesselin/projects/un
structured.html
More on spatial discretization and accuracy
Picture of stencil
1st and 2ndorder first derivative
approximation
u j
x
u j
x
Figure From
Computational Gasdynamics: Laney2


u j 1  u j
x j 1  x j
u j 1  u j 1
2( x j 1  x j 1 )
Discretization, Conservation and Flux
Capturing
Require scheme to capture shocks and
other discontinuities “automatically” and
not using “shock fitting methods”
Higher spatial order shock capturing
schemes (>2nd order) tend to be more
oscillatory around the discontinuities
because of the larger stencils required
thus more points are contributing
around areas with large gradients
Figure from
Computational Gasdynamics: Laney 2
Overview of Current Computational
Approaches
In general, methods in Computational Fluid
Dynamics can be divided into three approaches:
Finite Difference Methods
Advantages:
• Ease of Implementation
• Easy to make higher order
Disadvantages:
• Only applicable on structured
grids
Finite Element
Advantages:
• Can be any order of accuracy
• Based on variational methods
• Applicable on unstructured grids
Disadvantages:
• More complex
• Not conservative!
• Naturally implicit (can be explicit with
modifications)
Finite Volume
Advantages:
• Naturally Conservative (captures
discontinuities in the flow field)
• Many upwinding possibilities
• Applicable on unstructured grids
Disadvantages:
• Difficult to devise stable higher
order scheme
General Discontinuous Galerkin Setup
1. Start with the Euler Equation:
U
   F(U)  0
t
2. Discretize the spatial domain
and assume and assume an
approximate solution on a perelement basis
U  uh  u j (t )vj
3. Multiply by weight function and integrate by parts
 vk 


ui , j (t )vjd   vk  Fid   vkFiR  dS  0
t


Note the boundary term has a different flux term.
In normal finite element, the boundary terms need
to enforce connectivity with neighboring elements.
In Discontinuous Galerkin Methods the boundary
fluxes are calculated using the Riemann Fluxes.
This enforces connectivity and allows for
discontinuities in the solution!
One-Dimensional Discontinuous Galerkin1


U  F(U)  0
t
x

 
U   u 
E
 
 u 


F(U)   u 2  p 
u ( E  p )


Require an approximation to the solution in the form of:
k
u   al u (jl ) (t )vl( j ) ( x)dx
h
j
al 
l 1
Define the shape function as:
v1j  1 v2j  x  x j v2j  ( x  x j ) 2 
x lj
( j)
2
[
v
(
x
)]
dx
l

Ij
1
x 2j
12
Define the degrees of freedom as:
1
u(jl ) (t )  l  u( x, t )vl( j ) dx
x I j
Note:
1st DOF is the cell average of
the conservative variables
In Galerkin method the weight functions are taken to be the same as the shape functions.
h
Multiplying the Euler Equations by the weight functions and substituting u for U and
integrating by parts, we obtain the following form:
d
1
u (jl ) 
[vl( j ) ( x j 1 / 2 )f (u j 1 / 2 )  vl( j ) ( x j 1/ 2 )f (u j 1/ 2 )]
l
dt
x j

1
x lj
h
F
(
u
( x, t )

I
j
d
vl( j ) d x  0
dx
Shape Functions over each element
Exact solution to the Riemann Problem: Interface fluxes2
The second term of in the last equation has not been defined yet. How do we
get the fluxes at the cell interfaces? The Riemann Problem has an exact
solution!
Consider an Euler Equation with the initial of:
u
u ( x, t o )   L
u R
x  xo
x  xo
Expansion Fan:
Figure from
Computational Gasdynamics: Laney 2
Consider that every cell interface is a Riemann problem!
Exact solution to Riemann problem is very expensive and we are not interested in
in the solution at all x/t. Look for a suitable approximation for x/t=0 only via Roe
Averages


 1 


r1   u RL 
1 2 
 2 u RL 

r2 
1


 u RL  aRL 
h  a u 
RL RL 
 RL
 RL 
2aRL

r2  
All equations taken from Laney2


 u RL  aRL 
h  a u 
RL RL 
 RL
 
2aRL
1
Slope Limiting for Stability1
Around discontinuities DOFs representing the gradients are very large
causing oscillations and instabilities. To remedy this problem slope limiters
are introduced to insure stability
Runge-Kutta Time Explicit Time Marching
•Time integration of the equations will be carried out using a higher order
Runge-Kutta technique.
•The space discretization in the previous slide converted the PDEs into a system
of ODEs in time.
• Using Higher Order Runge-Kutta, we carry out the time integration on a perelement basis
t  CFL  min 


u j a j 
x j
Note:
Time Step is calculated based
on the largest Eigenvalue i.e.
fastest information transfer
DG Method
Test Problems
The method will first be implemented on the one-dimensional version of the Euler
Equations to test the methods accuracy. Sod’s Shock Tube Problem will be used as the test
case since it has an exact solution and will test the scheme’s shock capturing ability.
Implemented using Fortran 95
Exact Solution
Image: Author Generated
Image:
http://en.wikipedia.org/wiki/File:SodShockTubeTest_
Regions.png
1st Order Spatial Discretization Results
Density Evolution
Mach Number Evolution
2nd Order Spatial Discretization Results
Density Evolution
Mach Number Evolution
3rd Order Spatial Discretization Results
Density Evolution
Mach Number
DG Method
Test Problems
Two Dimensional Airfoil on a structured mesh. Mesh provided by Dr. Baeder.
Boundary Conditions:
Tangential Flow around the airfoil (Inviscid Wall)
Undisturbed flow at the domain edge (Farfield)
Run at a variety of mach numbers to create different
flow regimes (subsonic, trans-sonic, supersonic) to
steady state
Compare the results to experimental data (large
database for many airfoils)
Compare results to other established computational
tools
Implemented serially using Fortran 95, then
parallelized using MPI
http://www.salomeplatform.org/forum/forum_10/213329959/view
DG Method Implementation Original
Schedule
10/31/12-
12/15/1302/15/1303/15/13-
One dimensional version. Apply to one dimensional problem with a known
solution to test accuracy and shock capturing abilities. Sod shock tube problem.
Will validate the 1-D version (serial)
Two dimensional version. Apply to 2-D airfoil problem using
provided grids (serial)
Validation of the two dimensional version using experimental airfoil results as
well as the results published in literature
Parallelization of the two dimensional. Validate using results from the serial
version
Questions???
References:
1. Bernardo Cockburn, Chi-Wang Shu. TVB Runge-Kutta Local Projection
Discontinuous Galerkin Method for Conservation Laws II: General Frame
Work. Mathematics of Computation Volume 52, Issue 186, yr. 1989
2. Culbert B. Laney. Computational Gasdynamics. Cambridge University
Press. 1998.