Stability theory for finite-difference schemes using modified equations ´ Firas Dhaouadi 1 Emilie Duval 2 Sergey Tkachenko 3 Supervisors: Jean-Paul Vila 4 emy Baraille 5 R´ 1 Universit´ 2 Universit´ 3 Aix-Marseille Universit´ e Paul Sabatier e Grenoble Alpes e 4 Institut de Mathmatiques de Toulouse, INSA Toulouse 5 Service Hydrographique et Ocanographique de la Marine 22 August 2019 F. Dhaouadi, ´ E. Duval, S. Tkachenko (UPS, UGA, AMU) Project TOLOSA 22 August 2019 1 / 23
Introduction Modified equation Introduction to modified equations For instance, let us consider the scalar transport equation : ∂ u ∂ t + c ∂ u ∂ x = 0 where c > 0, and let us consider the following numerical scheme, for example : u n +1 + c u n i − u n − u n i − 1 i i = 0 (1) ∆ t ∆ x where ∆ t and ∆ x are respectively the discrete time step and the mesh size. After we use Taylor expansions in the vicinity of ( x i , t n ) ∂ t + ∆ t 2 ∂ 2 u n i + ∆ t ∂ u n = u ( x i , t n +1 ) = u ( x i , t n + ∆ t ) = u n u n +1 i ∂ t 2 + O (∆ t 3 ) i i 2 i − ∆ x ∂ u n ∂ x + ∆ x 2 ∂ 2 u n ∂ x 2 + O (∆ x 3 ) u n i − 1 = u ( x i − 1 , t n ) = u ( x i − ∆ x , t n ) = u n i i 2 and replace in the scheme (1) in order to get the scheme truncation error : ∂ 2 u ∂ 2 u ∂ u ∂ t = − c ∂ u ∂ x − ∆ t ∂ t 2 + c ∆ x ∂ x 2 + O (∆ t 2 , ∆ x 2 ) (2) 2 2 Now, for physical interpretation, we would like to have only space derivatives in the right hand side. F. Dhaouadi, ´ E. Duval, S. Tkachenko (UPS, UGA, AMU) Project TOLOSA 22 August 2019 2 / 23
Introduction Modified equation Introduction to modified equations For instance, let us consider the scalar transport equation : ∂ u ∂ t + c ∂ u ∂ x = 0 where c > 0, and let us consider the following numerical scheme, for example : u n +1 + c u n i − u n − u n i − 1 i i = 0 (1) ∆ t ∆ x where ∆ t and ∆ x are respectively the discrete time step and the mesh size. After we use Taylor expansions in the vicinity of ( x i , t n ) and replace in the scheme (1) in order to get the scheme truncation error : ∂ 2 u ∂ 2 u ∂ u ∂ t = − c ∂ u ∂ x − ∆ t ∂ t 2 + c ∆ x ∂ x 2 + O (∆ t 2 , ∆ x 2 ) (2) 2 2 F. Dhaouadi, ´ E. Duval, S. Tkachenko (UPS, UGA, AMU) Project TOLOSA 22 August 2019 2 / 23
Introduction Modified equation Introduction to modified equations For instance, let us consider the scalar transport equation : ∂ u ∂ t + c ∂ u ∂ x = 0 where c > 0, and let us consider the following numerical scheme, for example : u n +1 − u n + c u n i − u n i − 1 i i = 0 (1) ∆ t ∆ x where ∆ t and ∆ x are respectively the discrete time step and the mesh size. After we use Taylor expansions in the vicinity of ( x i , t n ) and replace in the scheme (1) in order to get the scheme truncation error : ∂ 2 u ∂ 2 u ∂ u ∂ t = − c ∂ u ∂ x − ∆ t ∂ t 2 + c ∆ x ∂ x 2 + O (∆ t 2 , ∆ x 2 ) (2) 2 2 Replace ∂ 2 u n ∂ t 2 by ∂ i ∂ t (2) in (2), then : ∂ 2 u ∂ 2 u ∂ u ∂ t + c ∂ u ∂ x = − c 2 ∆ t ∂ x 2 + c ∆ x ∂ x 2 + O (∆ t 2 , ∆ x 2 ) 2 2 � ∂ 2 u = c ∆ x � 1 − c ∆ t ∂ x 2 + O (∆ t 2 , ∆ x 2 ) 2 ∆ x F. Dhaouadi, ´ E. Duval, S. Tkachenko (UPS, UGA, AMU) Project TOLOSA 22 August 2019 2 / 23
Introduction Modified equation Introduction to modified equations For instance, let us consider the scalar transport equation : ∂ u ∂ t + c ∂ u ∂ x = 0 where c > 0, and let us consider the following numerical scheme, for example : u n +1 − u n + c u n i − u n i − 1 i i = 0 (1) ∆ t ∆ x where ∆ t and ∆ x are respectively the discrete time step and the mesh size. After we use Taylor expansions in the vicinity of ( x i , t n ) and replace in the scheme (1) in order to get the scheme truncation error : ∂ 2 u ∂ 2 u ∂ u ∂ t = − c ∂ u ∂ x − ∆ t ∂ t 2 + c ∆ x ∂ x 2 + O (∆ t 2 , ∆ x 2 ) (2) 2 2 Then, the modified equation is : � ∂ 2 u ∂ u ∂ t + c ∂ u ∂ x = c ∆ x � 1 − c ∆ t ∂ x 2 + O (∆ t 2 , ∆ x 2 ) 2 ∆ x The scheme is stable only if 1 − c ∆ t ∆ x ≥ 0 F. Dhaouadi, ´ E. Duval, S. Tkachenko (UPS, UGA, AMU) Project TOLOSA 22 August 2019 2 / 23
Introduction Heuristic stability theory Heuristic stability theory : heat equation Scheme PDE ∂ t − Q ∂ 2 u ∂ u u n +1 − u n − Q u n i +1 − 2 u n i + u n ∂ x 2 = 0 , Q > 0 i − 1 i i = 0 ∆ x 2 ∆ t This scheme is stable under the condition : Q ∆ t ∆ x 2 ≤ 1 2 Modified equation ∂ t = Q ∂ 2 u 6 Q ∆ t − ∆ x 2 � ∂ 4 u 360(∆ x 4 + 30 Q ∆ t ( − ∆ x 2 + 4 Q ∆ t )) ∂ 6 u ∂ u ∂ x 2 − Q ∂ x 4 + Q ∂ x 6 + O (∆ t 2 , ∆ x 4 ) � 12 F. Dhaouadi, ´ E. Duval, S. Tkachenko (UPS, UGA, AMU) Project TOLOSA 22 August 2019 3 / 23
Introduction Heuristic stability theory Heuristic stability theory : heat equation Scheme PDE ∂ t − Q ∂ 2 u ∂ u u n +1 − u n − Q u n i +1 − 2 u n i + u n ∂ x 2 = 0 , Q > 0 i − 1 i i = 0 ∆ x 2 ∆ t This scheme is stable under the condition : Q ∆ t ∆ x 2 ≤ 1 2 Modified equation ∂ t = Q ∂ 2 u 6 Q ∆ t − ∆ x 2 � ∂ 4 u 360(∆ x 4 + 30 Q ∆ t ( − ∆ x 2 + 4 Q ∆ t )) ∂ 6 u ∂ u ∂ x 2 − Q ∂ x 4 + Q ∂ x 6 + O (∆ t 2 , ∆ x 4 ) � 12 We look at the sign of the even order coefficients Q ≥ 0 − Q 12(6 Q ∆ t − ∆ x 2 ) ≤ 0 ⇔ Q ∆ t ∆ x 2 ≥ 1 6 ... F. Dhaouadi, ´ E. Duval, S. Tkachenko (UPS, UGA, AMU) Project TOLOSA 22 August 2019 3 / 23
Introduction Heuristic stability theory Heuristic stability theory : heat equation Scheme PDE ∂ t − Q ∂ 2 u ∂ u u n +1 − u n − Q u n i +1 − 2 u n i + u n ∂ x 2 = 0 , Q > 0 i − 1 i i = 0 ∆ x 2 ∆ t This scheme is stable under the condition : Q ∆ t ∆ x 2 ≤ 1 2 Modified equation ∂ t = Q ∂ 2 u 6 Q ∆ t − ∆ x 2 � ∂ 4 u 360(∆ x 4 + 30 Q ∆ t ( − ∆ x 2 + 4 Q ∆ t )) ∂ 6 u ∂ u ∂ x 2 − Q ∂ x 4 + Q ∂ x 6 + O (∆ t 2 , ∆ x 4 ) � 12 We look at the sign of the even order coefficients Q ≥ 0 − Q 12(6 Q ∆ t − ∆ x 2 ) ≤ 0 ⇔ Q ∆ t ∆ x 2 ≥ 1 6 ... This shows the limitations of the heuristic approach. F. Dhaouadi, ´ E. Duval, S. Tkachenko (UPS, UGA, AMU) Project TOLOSA 22 August 2019 3 / 23
Introduction Project goal and subject Project goal and subject Project goal To establish a clear link between the stability condition of a numerical scheme for a given PDE system and an associated modified equation. Next step : Once this is completed, retrieve the stability conditions presented in (P. Noble and J. P. Vila. “Stability theory for difference approximations of euler korteweg equations and application to thin film flows”. In: (2014), pp. 1–22. arXiv: arXiv:1304.3805v2 ) F. Dhaouadi, ´ E. Duval, S. Tkachenko (UPS, UGA, AMU) Project TOLOSA 22 August 2019 4 / 23
Introduction State of the art State of the art Presented the heuristic stability theory and also tried for non linear pdes. M. C. Hirt. “Heuristic Stability Theory for Finite-Difference Equations”. In: Journal of Computational Physics 2 (1968), pp. 339–355 The connection between the modified equation and the von Neumann (Fourier) method is established. R. F. Warming and B. J. Hyett. “The Modified Equation Approach to the Stability and Accuracy Analysis of Finite-Difference Methods”. In: Journal of Computational Physics 179 (1974), pp. 159–179 An alternative approach on how to derive the modified equation for linear problems Romuald Carpentier, Armel de la Bourdonnaye, and Larrouturou Bernard. “On the derivation of the modified equation for the analysis of linear numerical methods”. In: 31.4 (1997), pp. 459–470 F. Dhaouadi, ´ E. Duval, S. Tkachenko (UPS, UGA, AMU) Project TOLOSA 22 August 2019 5 / 23
Introduction Outline Linear scalar case 1 Von Neumann analysis Link with modified equations Algorithm Linear system case 2 Possible extensions Entropy stability Conclusions 3 F. Dhaouadi, ´ E. Duval, S. Tkachenko (UPS, UGA, AMU) Project TOLOSA 22 August 2019 6 / 23
Linear scalar case Von Neumann analysis Von Neumann stability analysis Given a linear PDE : u t + L x ( u ) = 0 assume we have a consistent one-step linear scheme given in general as : m r n r � � B q u n +1 A p u n j + q = j + p q = − m l p = − n l We replace every u n j by v ( k ) n exp ( ikj ∆ x ) and define the amplification factor as: � n r g ( k ) = v ( k ) n +1 p = − n l A p exp ( ipk ∆ x ) = � m r v ( k ) n q = − m l B q exp ( iqk ∆ x ) A necessary and sufficient stability condition is : | g ( k ) | ≤ 1 F. Dhaouadi, ´ E. Duval, S. Tkachenko (UPS, UGA, AMU) Project TOLOSA 22 August 2019 7 / 23
Linear scalar case Von Neumann analysis Von Neumann stability analysis In order to get a more practical formulation of this condition, we can show that the square of the modulus can be given by : | g ( k ) | 2 = 1 − 4 z r S ( z ) P ( z ) where : z = sin 2 ( k ∆ x / 2). r ≥ 1 is an integer. It is the maximum power of z that can be put as a common factor in the numerator. s α i z i is a polynomial function of z such that S (0) � = 0. � S ( z ) = i =0 d β i z i > 0 is a polynomial function of z such that P (0) = 1 � P ( z ) = i =0 Therefore, the stability condition becomes : | g ( k ) | 2 ≤ 1 ⇔ S ( z ) ≥ 0 ∀ z ∈ [0 , 1] F. Dhaouadi, ´ E. Duval, S. Tkachenko (UPS, UGA, AMU) Project TOLOSA 22 August 2019 8 / 23
Recommend
More recommend