Numerical Methods for Differential Equations 1.- Numerical Methods - PowerPoint PPT Presentation

Numerical Methods for Differential Equations 1.- Numerical Methods for DDEs Luis M. Abia, J. C. López Marcos, O. Angulo abia@mac.uva.es University of Valladolid Valladolid, (Spain) Euro Summer School Lipari (Sicilia, Italy) 13-26 September 2009– p. 1/38

Delay Differential Equations y ′ ( t ) = f ( t, y ( t ) , y ( t − τ ( t ))) , t 0 ≤ t ≤ T, t ∗ ≤ t ≤ t 0 , y ( t ) = φ ( t ) , f ( t, u , v ) , continuous and (locally) Lipschitz with respect its arguments u , v . t ∗ = τ, or 1. τ > 0 , constant , t ∗ = ´ 2. τ ≡ τ ( t ) ≥ 0 , ınf t ≥ t 0 ( t − τ ( t )) , and τ ( t ) continuous in [ t 0 , T ] or t ∗ = ´ 3. τ ≡ τ ( t, y ( t )) ≥ 0 , ınf t ≥ t 0 ( t − τ ( t, y ( t ))) , and τ ( t, u ) continuous and (locally) Lipschitz with respect to u . (State Dependent lag function case) y ′ ( t ) = f ( t, y ( t ) , y ( t − τ ( t )) , y ′ ( t − σ ( t ))) , t 0 ≤ t ≤ T, ( Neutral Differential Equations (NDE) ) We will denote α ( t ) = t − τ ( t, y ( t )) ≤ t, β ( t ) = t − σ ( t, y ( t )) ≤ t delayed arguments Euro Summer School Lipari (Sicilia, Italy) 13-26 September 2009– p. 2/38

Existence and Uniqueness of Solutions t ∗ ≤ t ≤ t 0 . y ′ ( t ) = f ( t, y ( t ) , y ( t − τ ( t, y ( t )))) , t ≥ t 0 , y ( t ) = φ ( t ) , Under ( H 1 ) There exists a τ 0 > 0 such that α ( t ) ≤ t − τ 0 , for t ∈ [ t 0 , T ] , or ( H ∗ R d . 1 ) There exists a τ 0 > 0 such that τ ( t, z ) ≥ τ 0 , for t ∈ [ t 0 , T ] , and z ∈ I (local) existence and uniqueness are derived easily from the existence and uniqueness theory for ODEs using the method of steps y ′ 1 ( t ) = f ( t, y 1 ( t ) , φ ( t − τ )) , t 0 ≤ t ≤ t 1 = t 0 + τ 0 , y ′ 2 ( t ) = f ( t, y 2 ( t ) , y 1 ( t − τ )) , t 1 ≤ t ≤ t 2 = t 1 + τ 0 , . . . y ′ m ( t ) = f ( t, y m ( t ) , y m − 1 ( t − τ )) , t m − 1 ≤ t ≤ t m = t m − 1 + τ 0 , Euro Summer School Lipari (Sicilia, Italy) 13-26 September 2009– p. 3/38

Existence and Uniqueness of Solutions t ∗ ≤ t ≤ t 0 . y ′ ( t ) = f ( t, y ( t ) , y ( t − τ ( t, y ( t )))) , t ≥ t 0 , y ( t ) = φ ( t ) , Theorem 1 (local existence) Let U and V be neighborhoods of Ψ( t 0 ) and Ψ( t 0 − τ ( t 0 , Ψ( t 0 ))) respectively, and assume that f ( t, u, v ) is continuous with respect to t and Lipschitz continuous with respect to u, v in [ t 0 , t 0 + h ] × U × V , for some h > 0 . Assume that the initial function Ψ( t ) is Lipschitz continuous for t ≤ t 0 and that the delay function τ ( t, y ) ≥ 0 is continuous with respect to t and Lipschitz continuous with respect to y in [ t 0 , t 0 + h ] × U . Then there exists a unique solution in [ t 0 , t 0 + δ ] for some δ > 0 and this solution depends continuously on the initial data R. D. Driver (1963), Hale (1986), Elsgolts and Norkin (1973) Euro Summer School Lipari (Sicilia, Italy) 13-26 September 2009– p. 4/38

Propagation of Discontinuities and Smoothing Assuming a discontinuity of the first derivative of the solution at t = t 0 = ξ 0 , 1 , and α ′ ( ξ 1 ,i ) � = 0 α ( ξ 1 ,i ) = t 0 , simple root f t + f y y ′ ( ξ 1 ,i ) + f x y ′ ( t 0 ) + α ′ ( ξ 1 ,i ) � = f t + f y y ′ ( ξ 1 ,i ) + f x Ψ ′ ( t 0 ) − α ′ ( ξ 1 ,i ) (first level primary discontinuity of y ′′ ) In general, solutions with odd multiplicity of α ( ξ k,j ) = ξ k − 1 ,i , for some ξ k − 1 ,i ( k -level primary discontinuity of y ( k +1) ) . Other discontinuities in the derivatives with respecto to t in the functions α ( t ) , f , and φ are called secondary discontinuities. Euro Summer School Lipari (Sicilia, Italy) 13-26 September 2009– p. 5/38

Propagation of Discontinuities and Smoothing Case of Constant Delay Propagation of Discontinuities (case of constant delay) 4 3 2 1 α (t) = t−2 0 −1 −2 0 0.5 1 1.5 2 2.5 3 3.5 4 ξ 0,1 = t 0 ξ 1,1 ξ 2,1 ξ 0 , 1 < ξ 1 , 1 < · · · < ξ k, 1 < · · · This is also the case when α ( t ) is an strictly increasing function. Euro Summer School Lipari (Sicilia, Italy) 13-26 September 2009– p. 6/38

Propagation of Discontinuities and Smoothing Case of vanishing discontinuity 10 5 α (t) 0 −5 −10 −15 −20 −25 ξ k+2,i ξ k,i+1 ξ ∗ ξ k+1,i+1 0 1 2 3 ξ k,i 4 ξ k+1,i 5 6 7 8 9 10 ξ k−1,i For general DDE, we will assume ( H 1 ) There exists a τ 0 > 0 such that α ( t ) ≤ t − τ 0 , for t ∈ [ t 0 , T ] , We replace [ t i , t i +1 ] with [ ξ i , ξ i +1 ] , i = 0 , 1 , 2 , . . . , where ξ 0 = t 0 , and for i ≥ 0 , ξ i +1 is the minimum root with odd multiplicity of α ( t ) = ξ i . (set of principal discontinuity points). Euro Summer School Lipari (Sicilia, Italy) 13 -26 September 2009– p. 7/38

Propagation of Discontinuities The extension of these ideas to general systems was developped by Willé and Baker (1992). In a system of DDEs discontinuities tracking can be complicated by discontinuities being propagated between solution components. This fact gives rise to the concept of strong and weak coupling and network dependency graphs. Strong coupling describes the propagation of discontinuities between different solution components by an ODE term. Weak coupling describes the propagation of discontinuities within the same solution component and between different solution components by a DDE term. For example, for the system y ′ 1 ( t ) = f 1 ( y 2 ( t ) , y 3 ( t − 1)) , y ′ t ≥ 0 2 ( t ) = f 2 ( y 3 ( t )) , y ′ 3 ( t ) = f 3 ( y 1 ( t ) , y 2 ( t − 1)) y 2 is strongly connected with y 1 , y 3 is strongly connected with y 2 , and y 1 is strongly connected with y 3 . However, y 3 is weakly connected with y 1 , and y 2 is weakly connected with y 3 . Euro Summer School Lipari (Sicilia, Italy) 13 -26 September 2009– p. 8/38

The Numerical Approach t ∗ ≤ t ≤ t 0 . y ′ ( t ) = f ( t, y ( t ) , y ( t − τ ( t ))) , t ≥ t 0 , y ( t ) = φ ( t ) , Variable step -size codes (with dense output) to solve for the y in the method of steps. . In the Method of Steps (constant lag scalar DDE) y ′ 1 ( t ) = f ( t, y 1 ( t ) , Ψ( t − τ )) , t 0 ≤ t ≤ t 1 = t 0 + τ, → ˜ y 1 ( t ) y ′ 2 ( t ) = f ( t, y 2 ( t ) , ˜ y 1 ( t − τ )) , t 1 ≤ t ≤ t 2 = t 1 + τ, → ˜ y 2 ( t ) . . . y ′ m ( t ) = f ( t, y m ( t ) , ˜ y m − 1 ( t − τ )) , t m − 1 ≤ t ≤ t m = t m − 1 + τ, → ˜ y m ( t ) Issues to consider for DDE 1. Dense Output to evaluate delayed solution values 2. The tracking of discontinuities in the solution 3. Vanishing lag delays and overlapping. Euro Summer School Lipari (Sicilia, Italy) 13-26 September 2009– p. 9/38 4. Linear stability properties

Continuous ODE Method y ′ ( t ) = g ( t, y ( t )) , t ≥ t 0 , y ( t 0 ) = y 0 . ODE method y n +1 = α n, 1 y n + · · · + α n,k y n − k +1 + h n +1 Φ( y n , . . . , y n − k +1 ; g, ∆ n ) , n ≥ k − 1 , y 0 , . . . , y k − 1 given , with ∆ n = { t n − k +1 , . . . , t n , t n +1 } Interpolation in [ t n , t n +1 ] , η ( t n + θh n +1 ) = β n, 1 ( θ ) y n + j n + · · · + β n,j n + i n +1 ( θ ) y n − i n + h n +1 Ψ( y n + j n , . . . , y n − i n ; θ, g, ∆ ′ n ) , 0 ≤ θ ≤ 1 , with ∆ ′ n = { t n − i n , . . . , t n + j n , t n + j n +1 } . η ( t n ) = y n , η ( t n +1 ) = y n +1 (continuity conditions). There existe an Ω > 0 , such that Ω − 1 h n +1+ i ≤ h n +1 ≤ Ω h n +1+ i , i = − k + 1 , . . . , 1 , Ω − 1 h n +1+ i ≤ h n +1 ≤ Ω h n +1+ i , i = − i n + 1 , . . . , j n , Euro Summer School Lipari (Sicilia, Italy) 13 -26 September 2009– p. 10/38

Continuous ODE Method We assume The increment functions Φ and Ψ in the continuous ODE method are Lipschitz continuous with respect their y arguments, and � Φ( y n , . . . , y n − k +1 ; ˜ g, ∆ n ) − Φ( y n , . . . , y n − k +1 ; g, ∆ n ) � g ∈ C 0 ≤ γ g sup � ˜ g ( t, y ) − g ( t, y ) � , ∀ ˜ t n − k +1 ≤ t ≤ t n +1 ,y g, ∆ ′ n ) − Ψ( y n + j n , . . . , y n − i n ; θ, g, ∆ ′ � Ψ( y n + j n , . . . , y n − i n ; θ, ˜ n ) � ≤ γ ′ g ∈ C 0 sup � ˜ g ( t, y ) − g ( t, y ) � , ∀ ˜ g t n − in ≤ t ≤ t n + jn +1 ,y   Let α n, 1 α n, 2 · · · α n,k − 1 α n,k   1 0 · · · 0 0       0 1 · · · 0 0 C n =     . . . .   . . . .   . . . .     0 0 · · · 1 0 be the companion matrix of the polynomial p n ( λ ) = λ k − α n, 1 λ k − 1 − · · · α n, 0 . Euro Summer School Lipari (Sicilia, Italy) 13 -26 September 2009– p. 11/38

Continuous ODE Method y ′ ( t ) = g ( t, y ( t )) , t ≥ t 0 , y ( t 0 ) = y 0 . The ODE method is consistent of order p , when p is the largest integer such that for all C q -continuous right-hand-side functions g and for all mesh points, we have y n +1 � = O ( h p +1 � u n +1 ( t n +1 ) − ˜ n +1 ) , ( h → 0) uniformly with respect to y ∗ n varying in bounded subsets, for n = 0 , . . . , N − 1 , where u n +1 ( t ) is the local solution given by the solution of  u ′ n +1 ( t ) = g ( t, u n +1 ( t )) , t n ≤ t ≤ t n +1 ,  u n +1 ( t n ) = y ∗  n and y n +1 ˜ = α n, 1 u n +1 ( t n ) + · · · + α n,k u n +1 ( t n − k +1 ) + h n +1 Φ( u n +1 ( t n ) , . . . , u n +1 ( t n − k +1 ); g, ∆ n ) , Euro Summer School Lipari (Sicilia, Italy) 13-26 September 2009– p. 12/38