Gene Golub SIAM Summer School 2012 Numerical Methods for Wave - PowerPoint PPT Presentation

Gene Golub SIAM Summer School 2012 Numerical Methods for Wave Propagation Finite Volume Methods Lecture 3 Randall J. LeVeque Applied Mathematics University of Washington R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

Outline • Shallow water equations with topography • Approximate Riemann solvers • f-wave formulation of wave-propagation method. • Well-balanced methods to preserve ocean-at-rest. • Dry state Riemann solvers R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

Great Tohoku Tsunami, 11 March 2011 Modeling and Simulating Tsunamis with an Eye to Hazard Mitigation, RJL and J. Behrens, SIAM News, May, 2011 http://www.siam.org/news/news.php?id=1882 R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

Great Tohoku Tsunami, 11 March 2011 Modeling and Simulating Tsunamis with an Eye to Hazard Mitigation, RJL and J. Behrens, SIAM News, May, 2011 http://www.siam.org/news/news.php?id=1882 R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

Great Tohoku Tsunami, 11 March 2011 R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

Shallow water equations h ( x, t ) = depth u ( x, t ) = velocity (depth averaged, varies only with x ) Conservation of mass and momentum hu gives system of two equations. mass flux = hu , momentum flux = ( hu ) u + p where p = hydrostatic pressure h t + ( hu ) x = 0 � � hu 2 + 1 2 gh 2 ( hu ) t + = 0 x Jacobian matrix: � � � 0 1 f ′ ( q ) = , λ = u ± gh. gh − u 2 2 u R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

Shallow water equations Hydrostatic pressure: Pressure at depth z > 0 below the surface is gz from weight of water above. Depth-averaged pressure is � h p = gz dz 0 � h � = 1 � 2 gz 2 � 0 = 1 2 gh 2 . R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

The Riemann problem Dam break problem for shallow water equations h t + ( hu ) x = 0 � 2 gh 2 � hu 2 + 1 ( hu ) t + x = 0 R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

The Riemann problem Dam break problem for shallow water equations h t + ( hu ) x = 0 � 2 gh 2 � hu 2 + 1 ( hu ) t + x = 0 R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

The Riemann problem Dam break problem for shallow water equations h t + ( hu ) x = 0 � 2 gh 2 � hu 2 + 1 ( hu ) t + x = 0 R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

The Riemann problem Dam break problem for shallow water equations h t + ( hu ) x = 0 � 2 gh 2 � hu 2 + 1 ( hu ) t + x = 0 R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

The Riemann problem Dam break problem for shallow water equations h t + ( hu ) x = 0 � 2 gh 2 � hu 2 + 1 ( hu ) t + x = 0 R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

Riemann solution for the SW equations in x - t plane Similarity solution: Solution is constant on any ray: q ( x, t ) = Q ( x/t ) Riemann solution can be calculated for many problems. Linear: Eigenvector decomposition. Nonlinear: more difficult. In practice “approximate Riemann solvers” used numerically. R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

An isolated shock If an isolated shock with left and right states q l and q r is propagating at speed s then the Rankine-Hugoniot condition must be satisfied: f ( q r ) − f ( q l ) = s ( q r − q l ) R m this can only hold for certain pairs q l , q r : For a system q ∈ l For a linear system, f ( q r ) − f ( q l ) = Aq r − Aq l = A ( q r − q l ) . So q r − q l must be an eigenvector of f ′ ( q ) = A . R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

An isolated shock If an isolated shock with left and right states q l and q r is propagating at speed s then the Rankine-Hugoniot condition must be satisfied: f ( q r ) − f ( q l ) = s ( q r − q l ) R m this can only hold for certain pairs q l , q r : For a system q ∈ l For a linear system, f ( q r ) − f ( q l ) = Aq r − Aq l = A ( q r − q l ) . So q r − q l must be an eigenvector of f ′ ( q ) = A . R m × m = A ∈ l ⇒ there will be m rays through q l in state space in the eigen-directions, and q r must lie on one of these. R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

An isolated shock If an isolated shock with left and right states q l and q r is propagating at speed s then the Rankine-Hugoniot condition must be satisfied: f ( q r ) − f ( q l ) = s ( q r − q l ) R m this can only hold for certain pairs q l , q r : For a system q ∈ l For a linear system, f ( q r ) − f ( q l ) = Aq r − Aq l = A ( q r − q l ) . So q r − q l must be an eigenvector of f ′ ( q ) = A . R m × m = A ∈ l ⇒ there will be m rays through q l in state space in the eigen-directions, and q r must lie on one of these. For a nonlinear system, there will be m curves through q l called the Hugoniot loci. R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

Riemann solution for the SW equations in x - t plane Nonlinear Riemann solution: If we know the 1-wave is a rarefaction wave, 2-wave is a shock, Can find h ∗ by solving: � � 1 � � � g h ∗ + 1 gh ∗ ) = u r + ( h ∗ − h r ) u l + 2( gh l − . 2 h r Expensive to do at every cell interface! R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

HLL Solver Harten – Lax – van Leer (1983): Use only 2 waves with s 1 = minimum characteristic speed s 2 = maximum characteristic speed W 1 = Q ∗ − Q ℓ , W 2 = Q r − Q ∗ Conservation implies unique value for middle state Q ∗ : s 1 W 1 + s 2 W 2 = f ( Q r ) − f ( Q ℓ ) ⇒ Q ∗ = f ( Q r ) − f ( Q ℓ ) − s 2 Q r + s 1 Q ℓ = . s 1 − s 2 R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

HLL Solver Harten – Lax – van Leer (1983): Use only 2 waves with s 1 = minimum characteristic speed s 2 = maximum characteristic speed W 1 = Q ∗ − Q ℓ , W 2 = Q r − Q ∗ Conservation implies unique value for middle state Q ∗ : s 1 W 1 + s 2 W 2 = f ( Q r ) − f ( Q ℓ ) ⇒ Q ∗ = f ( Q r ) − f ( Q ℓ ) − s 2 Q r + s 1 Q ℓ = . s 1 − s 2 Choice of speeds: • Max and min of expected speeds over entire problem, • Max and min of eigenvalues of f ′ ( Q ℓ ) and f ′ ( Q r ) . R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

HLLE Solver Einfeldt: Choice of speeds for gas dynamics (or shallow water) that guarantees positivity. Based on characteristic speeds and Roe averages: s 1 p (min( λ p i , ˆ λ p i − 1 / 2 = min i − 1 / 2 )) , p (max( λ p i +1 , ˆ λ p s 2 i − 1 / 2 = max i − 1 / 2 )) . where λ p i is the p th eigenvalue of the Jacobian f ′ ( Q i ) , ˆ λ p i − 1 / 2 is the p th eigenvalue using Roe average f ′ ( ˆ Q i − 1 / 2 ) R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

Approximate Riemann Solvers Approximate true Riemann solution by set of waves consisting of finite jumps propagating at constant speeds. Local linearization: Replace q t + f ( q ) x = 0 by q t + ˆ Aq x = 0 , where ˆ A = ˆ A ( q l , q r ) ≈ f ′ ( q ave ) . Then decompose r 1 + · · · α m ˆ q r − q l = α 1 ˆ r m to obtain waves W p = α p ˆ r p with speeds s p = ˆ λ p . R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

Approximate Riemann Solvers How to use? One approach: determine Q ∗ = state along x/t = 0 , � Q ∗ = Q i − 1 + W p , F i − 1 / 2 = f ( Q ∗ ) , p : s p < 0 A + ∆ Q = f ( Q i ) − F i − 1 / 2 . A − ∆ Q = F i − 1 / 2 − f ( Q i − 1 ) , R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

Approximate Riemann Solvers How to use? One approach: determine Q ∗ = state along x/t = 0 , � Q ∗ = Q i − 1 + W p , F i − 1 / 2 = f ( Q ∗ ) , p : s p < 0 A + ∆ Q = f ( Q i ) − F i − 1 / 2 . A − ∆ Q = F i − 1 / 2 − f ( Q i − 1 ) , Wave-propagation algorithm uses: � � s p W p , A + ∆ Q = s p W p . A − ∆ Q = p : s p < 0 p : s p > 0 Conservative only if A − ∆ Q + A + ∆ Q = f ( Q i ) − f ( Q i − 1 ) . This holds for Roe solver. R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012

Roe Solver Solve q t + ˆ Aq x = 0 where ˆ A satisfies ˆ A ( q r − q l ) = f ( q r ) − f ( q l ) . Then: • Good approximation for weak waves (smooth flow) • Single shock captured exactly: ⇒ q r − q l is an eigenvector of ˆ f ( q r ) − f ( q l ) = s ( q r − q l ) = A • Wave-propagation algorithm is conservative since � � s p i − 1 / 2 W p W p A − ∆ Q i − 1 / 2 + A + ∆ Q i − 1 / 2 = i − 1 / 2 = A i − 1 / 2 . R.J. LeVeque, University of Washington Gene Golub SIAM Summer School 2012