Well-balanced DG scheme for Euler equations with gravity Praveen. C - PowerPoint PPT Presentation

Well-balanced DG scheme for Euler equations with gravity Praveen. C praveen@tifrbng.res.in Center for Applicable Mathematics Tata Institute of Fundamental Research Bangalore 560065 Oberseminar Mathematische Str¨ omungsmechanik Institut f¨ ur Mathematik-der Julius-Maximilians-Universit¨ at W¨ urzburg 16 April 2015 Supported by Airbus Foundation Chair at TIFR-CAM/ICTS 1 / 57

Euler equations with gravity Flow properties ρ = density , u = velocity p = pressure , E = total energy Gravitational potential Φ ; force per unit volume of fluid − ρ ∇ Φ System of conservation laws ∂ρ ∂t + ∂ ∂x ( ρu ) = 0 ∂t ( ρu ) + ∂ ∂ − ρ∂ Φ ∂x ( p + ρu 2 ) = ∂x ∂E ∂t + ∂ − ρu∂ Φ ∂x ( E + p ) u = ∂x 2 / 57

Euler equations with gravity Perfect gas assumption � � E − 1 γ = c p 2 ρu 2 p = ( γ − 1) , > 1 c v In compact notation   0 ∂ q ∂t + ∂ f  ∂ Φ ∂x = − ρ  ∂x ρu where     ρ ρu  , p + ρu 2 q = ρu f =    E ( E + p ) u 3 / 57

Hydrostatic solutions • Fluid at rest u e = 0 • Mass and energy equation satisfied • Momentum equation d p e dΦ d x = − ρ e (1) d x • Need additional assumptions to solve this equation • Assume ideal gas model p = ρRT, R = gas constant 4 / 57

Hydrostatic solutions • Isothermal hydrostatic state, i.e., T e ( x ) = T e = const , then � Φ( x ) � p e ( x ) exp = const (2) RT e Density ρ e ( x ) = p e ( x ) RT e • Polytropic hydrostatic state, then we have following relations ν 1 − − p e ρ − ν ν − 1 ν − 1 = const , p e T = const , ρ e T = const (3) e e e where ν > 1 is some constant. From (1) and (3), we obtain νRT e ( x ) + Φ( x ) = const (4) ν − 1 E.g., pressure is ν − 1 p e ( x ) = C 1 [ C 2 − Φ( x )] ν 5 / 57

Existing schemes • Finite volumes ◮ Isothermal case: Xing and Shu [3], well-balanced WENO scheme ◮ If ν = γ we are in isentropic case h ( x ) + Φ( x ) = const has been considered by Kappeli and Mishra [1]. ◮ Desveaux et al: Relaxation schemes, general hydrostatic states ◮ PC and Klingenberg: isothermal and polytropic, ideal gas • DG ◮ Yulong Xing: ideal gas, isothermal well-balanced 6 / 57

Why well-balanced scheme • Scheme is well-balanced if it exactly preserves hydrostatic solution. • General evolutionary PDE ∂ q ∂t = R ( q ) • Stationary solution q e R ( q e ) = 0 • We are interested in computing small perturbations q ( x, 0) = q e ( x ) + ε ˜ q ( x, 0) , ε ≪ 1 • Perturbations are governed by linear equation ∂ ˜ q ∂t = R ′ ( q e )˜ q 7 / 57

Why well-balanced scheme • Some numerical scheme ∂ q h ∂t = R h ( q h ) • q h,e = interpolation of q e onto the mesh • Scheme is well balanced if ∂ q h R h ( q h,e ) = 0 = ⇒ ∂t = 0 • Suppose scheme is not well-balanced R h ( q h,e ) � = 0 . Solution q h ( x, t ) = q h,e ( x ) + ε ˜ q h ( x, t ) 8 / 57

Why well-balanced scheme • Linearize the scheme around q h,e ∂ q h ) = R h ( q h,e ) + εR ′ ∂t ( q h,e + ε ˜ q h ) = R h ( q h,e + ε ˜ h ( q h,e )˜ q h or ∂ ˜ ∂t = 1 q h εR h ( q h,e ) + R ′ h ( q h,e )˜ q h • Scheme is consistent of order r : R h ( q h,e ) = Ch r � q h,e � ∂ ˜ ∂t = 1 q h εCh r � q h,e � + R ′ h ( q h,e )˜ q h • ε ≪ 1 then first term may dominate the second term; need h ≪ 1 • Well-balanced scheme: R h ( q h,e ) = 0 ∂ ˜ q h ∂t = R ′ h ( q h,e )˜ q h 9 / 57

Scope of present work • Nodal DG scheme ◮ Gauss-Lobatto-Legendre nodes ◮ Arbitrary quadrilateral cells in 2-D • Ideal gas model: well-balanced for isothermal hydrostatic solutions • Any consistent numerical flux 10 / 57

Finite element and well-balanced Consider conservation law with source term q t + f ( q ) x = s ( q ) which has a stationary solution q e f ( q e ) x = s ( q e ) Weak formulation : Find q ( t ) ∈ V such that d d t ( q, φ ) + a ( q, φ ) = ( s ( q ) , φ ) , ∀ φ ∈ V Finite element scheme : Find q h ( t ) ∈ V h such that d d t ( q h , φ h ) + a ( q h , φ h ) = ( s ( q h ) , φ h ) , ∀ φ h ∈ V h 11 / 57

Finite element and well-balanced This is not true in general since we need to use quadratures. FEM with quadrature : Find q h ( t ) ∈ V h such that d d t ( q h , φ h ) h + a h ( q h , φ h ) = ( s h ( q h ) , φ h ) h , ∀ φ h ∈ V h Let q e,h = Π h ( q e ) , Π h : V → V h , interpolation or projection For the above scheme to be well-balanced, we require that a h ( q e,h , φ h ) = ( s h ( q e,h ) , φ h ) h , ∀ φ h ∈ V h because if q h (0) = q e,h = ⇒ q h ( t ) = q e,h ∀ t 12 / 57

Mesh and basis functions • Partition domain into disjoint cells C i = ( x i − 1 2 , x i + 1 2 ) , ∆ x i = x i + 1 2 − x i − 1 2 • Approximate solution inside each cell by a polynomial of degree N • Map C i to a reference cell, say [0 , 1] x = ξ ∆ x i + x i − 1 (5) 2 • On this reference cell let ξ j , 0 ≤ j ≤ N be the Gauss-Lobatto-Legendre nodes, roots of (1 − ξ 2 ) P ′ N ( ξ ) • ℓ j ( ξ ) : nodal Lagrange basis functions using these GLL points, with the interpolation property ℓ j ( ξ k ) = δ jk , 0 ≤ j, k ≤ N 13 / 57

Mesh and basis functions • Basis functions in physical coordinates φ j ( x ) = ℓ j ( ξ ) , 0 ≤ j ≤ N • Derivatives of the shape functions φ j : apply the chain rule of differentiation d j ( ξ ) d ξ 1 d xφ j ( x ) = ℓ ′ ℓ ′ d x = j ( ξ ) ∆ x i • x j ∈ C i denote the physical locations of the GLL points x j = ξ j ∆ x i + x i − 1 2 , 0 ≤ j ≤ N 14 / 57

Discontinuous Galerkin Scheme Consider the single conservation law with source term ∂q ∂t + ∂f ∂x = s Solution inside cell C i is polynomial of degree N N � q h ( x, t ) = q j ( t ) φ j ( x ) , q j ( t ) = q h ( x j , t ) j =0 Also approximate flux N N � � f h ( x, t ) = f ( q h ( x j , t )) φ j ( x ) = f j ( t ) φ j ( x ) j =0 j =0 15 / 57

Discontinuous Galerkin Scheme Gauss-Lobatto-Legendre quadrature N � � φ ( x ) ψ ( x )d x ≈ ( φ, ψ ) h = ( φ, ψ ) N,C i = ∆ x i ω q φ ( x q ) ψ ( x q ) C i q =0 GLL weights ω q correspond to the reference interval [0 , 1] . Semi-discrete DG: For 0 ≤ j ≤ N d d t ( q h , φ j ) h + ( ∂ x f h , φ j ) h +[ ˆ 2 − f h ( x − 2 )] φ j ( x − f i + 1 2 ) i + 1 i + 1 (6) − [ ˆ 2 − f h ( x + 2 )] φ j ( x + f i − 1 2 ) = ( s h , φ j ) h i − 1 i − 1 where ˆ 2 = ˆ f ( q − 2 , q + f i + 1 2 ) is a numerical flux function. i + 1 i + 1 16 / 57

Numerical flux Numerical flux is consistent: ˆ f ( q , q ) = f ( q ) Def: Contact property The numerical flux ˆ f is said to satisfy contact property if for any two states q L = [ ρ L , 0 , p/ ( γ − 1)] q R = [ ρ R , 0 , p/ ( γ − 1)] and we have ˆ f ( q L , q R ) = [0 , p, 0] ⊤ • states q L , q R in the above definition correspond to a stationary contact discontinuity. • Contact Property = ⇒ numerical flux exactly support a stationary contact discontinuity. • Examples of such numerical flux: Roe, HLLC, etc. 17 / 57

Approximation of source term Let ¯ T i = temperature corresponding to the cell average value in cell C i Rewrite the source term in the momentum equation as (Xing & Shu) � Φ � ∂ � � s = − ρ∂ Φ − Φ ∂x = ρR ¯ T i exp ∂x exp R ¯ R ¯ T i T i Source term approximation: For x ∈ C i � ∂ � Φ( x ) N � � − Φ( x j ) � s h ( x ) = ρ h ( x ) R ¯ T i exp exp φ j ( x ) (7) R ¯ R ¯ T i ∂x T i j =0 Source term in the energy equation 1 ( ρu ) h s h ρ h 18 / 57

Well-balanced property Well-balanced property Let the initial condition to the DG scheme (6), (7) be obtained by interpolating the isothermal hydrostatic solution corresponding to a continuous gravitational potential Φ . Then the scheme (6), (7) preserves the initial condition under any time integration scheme. Proof : For continuous hydrostatic solution, by flux consistency ˆ ˆ 2 − f h ( x − 2 − f h ( x + f i + 1 2 ) = 0 , f i − 1 2 ) = 0 i + 1 i − 1 Above is true even if density is discontinuous, provided flux satisfies contact property. = ⇒ density and energy equations are well-balanced 19 / 57

Well-balanced property The flux f h has the form N � f h ( x, t ) = p j ( t ) φ j ( x ) , p j = pressure at the GLL point x j j =0 Isothermal initial condition, ¯ T i = T e = const . The source term evaluated 20 / 57

Well-balanced property at any GLL node x k is given by � ∂ � N � Φ( x k ) � − Φ( x j ) � s h ( x k ) = ρ h ( x k ) RT e exp exp ∂xφ j ( x k ) RT e RT e � �� � j =0 p k � ∂ � N � Φ( x k ) � − Φ( x j ) � = p k exp exp ∂xφ j ( x k ) RT e RT e j =0 � ∂ N � Φ( x k ) � � − Φ( x j ) � = p k exp exp ∂xφ j ( x k ) RT e RT e j =0 � ∂ N � Φ( x j ) � � − Φ( x j ) � = p j exp exp ∂xφ j ( x k ) RT e RT e j =0 N ∂xφ j ( x k ) = ∂ ∂ � = p j ∂xf h ( x k ) j =0 21 / 57

Well-balanced property Since ∂ x f h ( x k ) = s h ( x k ) at all the GLL nodes x k we can conclude that ( ∂ x f h , φ j ) h = ( s h , φ j ) h , 0 ≤ j ≤ N and hence the scheme is well-balanced for the momentum equation also. 22 / 57