Boundary Approximations for Semi-Lagrangian Schemes Applied to - PowerPoint PPT Presentation

Overview Boundary treatment Analysis Linear solvers Conclusions Boundary Approximations for Semi-Lagrangian Schemes Applied to Hamilton-Jacobi-Bellman Equations Christoph Reisinger ⋆ Joint work with Julen Rotaetxe Arto ⋆ ⋆ Mathematical Institute University of Oxford Linz, 23 November 2016

Overview Boundary treatment Analysis Linear solvers Conclusions M OTIVATION Background: ◮ Monotone schemes have “wide stencils” (Motzkin and Wasow, 1952). ◮ They “overstep” domain boundaries. ◮ Standard iterative (multilevel) schemes perform poorly. We propose a(n) ◮ boundary treatment with same local accuracy as in interior; ◮ unconditionally stable and monotone implicit scheme; ◮ analysis of a convergent aggregation-based multigrid scheme.

Overview Boundary treatment Analysis Linear solvers Conclusions HJB EQUATIONS Consider α ∈A { L α [ u ]( t , x ) + f α ( t , x ) } = 0 , u t − inf ( t , x ) ∈ ( 0 , T ] × Ω , x ∈ ¯ u ( 0 , x ) = g ( x ) , Ω , u ( t , x ) = ψ ( x ) , ( t , x ) ∈ ( 0 , T ] × ∂ Ω , where Ω is a bounded domain, ¯ Ω := Ω ∪ ∂ Ω ⊆ R d , A is a compact set, and L α [ u ]( t , x ) = tr [ a α ( t , x ) D 2 u ( t , x )] + b α ( t , x ) Du ( t , x ) , a α = 1 2 σ α σ α, T with σ α ∈ R d × P , b α , f α , g and ψ take values in R d , R , R , and R .

Overview Boundary treatment Analysis Linear solvers Conclusions S EMI -L AGRANGIAN SCHEMES Define the L inear I nterpolation S emi- L agrangian ( LISL ) scheme by ( I ∆ x φ )( t , x + y α, + ( t , x )) − 2 ( I ∆ x φ )( t , x )+( I ∆ x φ )( t , x + y α, − ( t , x )) ∆ x [ I ∆ x φ ]( t , x ) := � M L α p p , p = 1 2 ∆ x as in Debrabant and Jakobsen (2013), where I is bilinear interpolation and √ y α, ± y α, ± ∆ x σ α P + 1 = ∆ xb α , = ± p for p ≤ P , and M = P + 1. p We consider ◮ an explicit Euler scheme: V n + 1 − V n α ∈A L α ∆ x [ I ∆ x V n ] = 0 ; − inf ∆ t ◮ an implicit Euler scheme: V n + 1 − V n α ∈A L α ∆ x [ I ∆ x V n + 1 ] = 0 . − inf ∆ t In practice, we replace A by A with |A| = N α .

Overview Boundary treatment Analysis Linear solvers Conclusions W IDE STENCILS AND BOUNDARIES √ The stencil “oversteps” in a layer of width k = ∆ x around the boundary. Figure: Truncation and extrapolation of the stencil for an elliptical domain and a mesh made of square cells. ∂ Ω Ω There are two situations where interpolation at x + y α, ± ( t , x ) is not possible: p A. x + y α, ± ∈ ¯ ( t , x ) / Ω (bottom left in Fig.); p B. x + y α, ± ( t , x ) ∈ ¯ Ω , but “its” element has vertices outside ¯ Ω (top right). p

Overview Boundary treatment Analysis Linear solvers Conclusions D IFFERENT TREATMENTS ∂ Ω Ω ◮ Constant extrapolation in stencil direction. ◮ Linear extrapolation in stencil direction. ◮ Stencil truncation: Consider instead ˆ L α ∆ x [ I ∆ x φ ]( t , x ) := M A α y α, + ( t , x )) − ( A α p + B α p ) φ ( t , x ) + B α y α, − p ( I ∆ x φ )( t , x + ˆ p ( I ∆ x φ )( t , x + ˆ ( t , x )) p p � , 2 ∆ x p = 1 where y α, ± = µ α, ± y α, ± where µ α, ± µ ≥ 0 : x + µ y α, ± ˆ , ( t , x ) = min � ( t , x ) ∈ ∂ Ω � . p p p p p and A α p ≡ A α p ( t , x ) and B α p ≡ B α p ( t , x ) , such that ˆ L α ∆ x a consistent approximation as ∆ x → 0.

Overview Boundary treatment Analysis Linear solvers Conclusions T HE TEST PROBLEM Problem A (see Section 9.3 from Debrabant and Jacobsen (2013)). It has exact solution � 3 � u ( t , x 1 , x 2 ) = 2 − t sin x 1 sin x 2 , and coefficients and control set are given by � 1 � � 3 � �� f α = cos 2 x 1 sin 2 x 2 + sin 2 x 1 cos 2 x 2 + 2 − t sin x 1 sin x 2 + 2 − t � − 2 sin ( x 1 + x 2 ) cos ( x 1 + x 2 ) cos x 1 cos x 2 , √ � sin ( x 1 + x 2 ) � b α = α, σ α = A = { α ∈ R 2 : α 2 1 + α 2 2 , 2 = 1 } . cos ( x 1 + x 2 ) ◮ The solution is periodic, and D. & J. (2013) use periodic boundary conditions in space for ( t , x 1 , x 2 ) ∈ [ 0 , T ] × [ − π, π ] 2 with T = 1 2 . ◮ Here, we use Dirichlet conditions .

Overview Boundary treatment Analysis Linear solvers Conclusions T RUNCATED STENCILS Figure: Problem A: Left the stencil over a 11 × 11 grid and N α = 10 equally spaced points in A . Right the histograms of the displacement from the central node for 641 × 641 grid. Stencil 4 B = 1.71 3 A = 1.26 C = 1.00 C = 1.00 B = 2.14 2 A = 1.18 A = 1.00 1 B = 1.00 x2 0 C = 1.00 C = 1.00 A = 2.14 -1 B = 1.00 C = 1.00 B = 1.28 B = 1.26 A = 1.00 -2 C = 1.00 -3 A = 2.29 -4 -4 -3 -2 -1 0 1 2 3 4 x1 (b) The radius of the stencil in σ α is 14.27 (a) The finite difference weights are A ≡ A α 2 , 1 ( x ) , B ≡ B α 2 , 1 ( x ) and for this grid, given by � σ α � 2 � √ ∆ x = 640 /π . C ≡ ( µ α 2 , 2 ( x )) − 1 , N α = 10.

Overview Boundary treatment Analysis Linear solvers Conclusions C ONSTANT AND LINEAR EXTRAPOLATION Table: Explicit Euler LISL method with N α = 40 for Problem A. (a) Constant extrapolation: L ∞ error over [ − π, π ) 2 3 ∆ t = ∆ x ∆ t = ∆ x 2 ∆ t = T ∆ t = ∆ x 2 N x 4 error rate error rate error rate error rate 41 1.36e+00 - 3.68e-01 - 3.72e-01 - 3.65e-01 - 81 1.89e+00 -0.48 2.61e-01 0.49 2.62e-01 0.51 2.60e-01 0.49 161 2.67e+00 -0.49 1.80e-01 0.54 1.80e-01 0.54 1.80e-01 0.53 321 3.77e+00 -0.50 1.27e-01 0.51 1.27e-01 0.51 1.27e-01 0.51 641 5.34e+00 -0.50 9.18e-02 0.47 9.18e-02 0.47 9.18e-02 0.46 (b) Linear extrapolation: L ∞ error over [ − π, π ) 2 3 ∆ t = ∆ x ∆ t = ∆ x 2 ∆ t = T ∆ t = ∆ x 2 4 N x error rate error rate error rate error rate 41 1.59e-01 - 1.04e-01 - 1.05e-01 - 1.03e-01 - 81 8.15e-02 0.96 5.25e-02 0.99 5.26e-02 1.00 5.22e-02 0.98 161 4.28e-02 0.93 5.62e-01 -3.42 5.63e-01 -3.42 5.58e-01 -3.42 321 2.75e-02 0.64 4.41e+03 -12.94 6.00e+03 -13.38 8.00e+03 -13.81 641 1.85e-02 0.57 2.77e+20 -55.80 2.70e+20 -55.32 1.37e+21 -57.25

Overview Boundary treatment Analysis Linear solvers Conclusions S TENCIL TRUNCATION ON TWO DOMAINS Table: Explicit Euler LISL method for Problem A. (a) L ∞ error over [ − π, π ) 2 3 ∆ t = ∆ x ∆ t = ∆ x 2 ∆ t = T ∆ t = ∆ x 2 N x 4 error rate error rate error rate error rate 41 1.42e-01 - 4.39e-02 - 4.39e-02 - 4.36e-02 - 81 1.04e-01 0.45 2.12e-02 1.05 2.11e-02 1.06 2.11e-02 1.05 161 7.36e-02 0.50 1.10e-02 0.94 1.10e-02 0.94 1.10e-02 0.94 321 5.28e-02 0.48 1.34e+23 -83.33 5.77e-03 0.93 5.76e-03 0.93 641 3.77e-02 0.48 5.07e+89 -221.17 3.10e-03 0.90 3.10e-03 0.89 (b) L ∞ error over [ − π 8 , 15 π 8 ) 2 3 ∆ t = ∆ x ∆ t = ∆ x 2 ∆ t = T ∆ t = ∆ x 2 4 N x error rate error rate error rate error rate 41 1.55e-01 - 4.71e-02 - 4.76e-02 - 4.67e-02 - 81 1.12e-01 0.47 1.57e+05 -21.67 7.90e+05 -23.98 2.11e-02 1.15 161 8.04e-02 0.47 1.02e+33 -92.39 1.30e+35 -97.06 1.10e-02 0.94 321 5.80e-02 0.47 6.73e+103 -235.26 5.96e+138 -344.35 5.76e-03 0.93 641 4.22e-02 0.46 8.17e+276 -574.97 NaN NaN 3.10e-03 0.89

Overview Boundary treatment Analysis Linear solvers Conclusions E MPIRICAL FINDINGS AND CONSISTENCY ◮ Constant extrapolation loses accuracy near boundary. ◮ Linear extrapolation is unstable. ◮ Stencil truncation has stricter CFL condition than in interior. ◮ CFL condition depends on boundary. The local LISL truncation error is y α, + y α, − O (∆ x ) if neither ˆ nor ˆ overstep; ◮ p p O (∆ x 1 / 2 ) y α, + y α, − if either ˆ or ˆ overstep; ◮ p p y α, + y α, − O ( 1 ) if both ˆ and ˆ overstep (cf. Max Jensen’s talk). ◮ p p (NB: The error is actually O (∆ x ) in the last case if the exact Dirichlet data are used.)

Overview Boundary treatment Analysis Linear solvers Conclusions L OCAL BOUNDARY REFINEMENT Define (with N ( x ) the corners of the element x is in) √ Ω ( 1 ) := { x ∈ Ω ∆ x : N ( x ± σ ∆ x ) �⊂ Ω ∆ x } , ∆ x √ � � Ω ( 2 ) � Ω ∆ x \ Ω ( 1 ) := N ( x ± σ ∆ x ) ∩ , ∆ x ∆ x x ∈ Ω ( 1 ) ∆ x Ω ( 3 ) Ω ∆ x \ (Ω ( 1 ) ∆ x ∪ Ω ( 2 ) := ∆ x ) . ∆ x ∆ x with mesh size ∆ x ( i ) and step k ( i ) , Now refine Ω ( i ) ∆ x ( 1 ) = O (∆ x 3 / 2 ) , k ( 1 ) = O (∆ x ) , √ ∆ x ( 2 ) = O (∆ x 3 / 2 ) , k ( 2 ) = ∆ x , √ ∆ x ( 3 ) = ∆ x , k ( 3 ) = ∆ x . Then the local error is O (∆ x ) everywhere, and the total # of points O ( | Ω ∆ x | ) . Practically not necessary so not done in the following computations.

Overview Boundary treatment Analysis Linear solvers Conclusions S TABILITY Recall ˆ L α ∆ x [ I ∆ x φ ]( t , x ) := M y α, + y α, − A α ( t , x )) − ( A α p + B α p ) φ ( t , x ) + B α p ( I ∆ x φ )( t , x + ˆ p ( I ∆ x φ )( t , x + ˆ ( t , x )) p p � . 2 ∆ x p = 1 The θ -scheme ( θ = 0 explicit, θ = 1 implicit) is monotone if   M A α p + B α p �  ≤ 1 . ( 1 − θ )∆ t n  2 ∆ x p = 1 This implies that for 0 ≤ θ < 1 the scheme is monotone ( − → ℓ ∞ − stable) if y α, + y α, − ∆ t ≤ C ∆ x and neither ˆ nor ˆ overstep; ◮ p p ∆ t ≤ C ∆ x 3 / 2 y α, + y α, − and either ˆ or ˆ overstep; ◮ p p y α, + y α, − ∆ t ≤ C ∆ x 2 and both ˆ and ˆ overstep. ◮ p p