Computational complexity of solving polynomial differential - PowerPoint PPT Presentation

Computational complexity of solving polynomial differential equations over unbounded domains Amaury Pouly ⋆, † Daniel Graça † , ‡ ⋆ Ecole Polytechnique, LIX, 91128 Palaiseau Cedex, France † CEDMES/FCT, Universidade do Algarve, C. Gambelas, 8005-139 Faro, Portugal ‡ SQIG /Instituto de Telecomunicações, Lisbon, Portugal July 8, 2013 A. Pouly, D. Graça (LIX, FCT) Complexity of solving PIVP July 8, 2013 −∞ / 18

Outline Introduction 1 Motivation Existing results Practice Theory Goal and result Complexity of solving PIVP 2 Crash course on numerical methods Euler method Taylor method Basic algorithm Enhanced algorithm Conclusion 3 A. Pouly, D. Graça (LIX, FCT) Complexity of solving PIVP July 8, 2013 −∞ / 18

Introduction Motivation Problem statement We want to solve: � y ′ = p ( y ) y ( t 0 )= y 0 where y : I ⊆ R → R n p : vector of polynomials Solve ? A. Pouly, D. Graça (LIX, FCT) Complexity of solving PIVP July 8, 2013 1 / 18

Introduction Motivation Problem statement We want to solve: � y ′ = p ( y ) y ( t 0 )= y 0 where y : I ⊆ R → R n p : vector of polynomials Solve ? ⊲ Compute y i ( t ) with arbitrary precision for any t ∈ I A. Pouly, D. Graça (LIX, FCT) Complexity of solving PIVP July 8, 2013 1 / 18

Introduction Motivation Problem statement We want to solve: � y ′ = p ( y ) y ( t 0 )= y 0 where y : I ⊆ R → R n p : vector of polynomials Solve ? ⊲ Compute y i ( t ) with arbitrary precision for any t ∈ I Example  c ′ ( t )= − s ( t )  c ( 0 )= 1   s ′ ( t )= c ( t ) s ( 0 )= 0 x ′ ( t )= 2 c ( t ) s ( t ) x ( t ) 2 x ( t )= 1   2 A. Pouly, D. Graça (LIX, FCT) Complexity of solving PIVP July 8, 2013 1 / 18

Introduction Motivation Problem statement We want to solve: � y ′ = p ( y ) y ( t 0 )= y 0 where y : I ⊆ R → R n p : vector of polynomials Solve ? ⊲ Compute y i ( t ) with arbitrary precision for any t ∈ I Example  c ′ ( t )= − s ( t )   c ( t )= cos ( t ) c ( 0 )= 1    s ′ ( t )= c ( t ) s ( t )= sin ( t ) s ( 0 )= 0 � 1 x ′ ( t )= 2 c ( t ) s ( t ) x ( t ) 2 x ( t )= 1 x ( t )=    1 + cos ( t ) 2 2 A. Pouly, D. Graça (LIX, FCT) Complexity of solving PIVP July 8, 2013 1 / 18

Introduction Motivation Motivation Theoretical complexity of solving differential equations Functions generated by the General Purpose Analog Computer (GPAC) Solve y ′ = f ( y ) where f is elementary (composition of polynomials, exponential,logarithms, (inverse) trigonometric functions, ...) A. Pouly, D. Graça (LIX, FCT) Complexity of solving PIVP July 8, 2013 2 / 18

Introduction Motivation Motivation Theoretical complexity of solving differential equations Functions generated by the General Purpose Analog Computer (GPAC) Solve y ′ = f ( y ) where f is elementary (composition of polynomials, exponential,logarithms, (inverse) trigonometric functions, ...) A. Pouly, D. Graça (LIX, FCT) Complexity of solving PIVP July 8, 2013 2 / 18

Introduction Motivation Motivation Theoretical complexity of solving differential equations Functions generated by the General Purpose Analog Computer (GPAC) Solve y ′ = f ( y ) where f is elementary (composition of polynomials, exponential,logarithms, (inverse) trigonometric functions, ...) A. Pouly, D. Graça (LIX, FCT) Complexity of solving PIVP July 8, 2013 2 / 18

Introduction Motivation Motivation Theoretical complexity of solving differential equations Functions generated by the General Purpose Analog Computer (GPAC) Solve y ′ = f ( y ) where f is elementary (composition of polynomials, exponential,logarithms, (inverse) trigonometric functions, ...) Example � y ′ = sin ( y ) y ′ = z   y ( 0 )= 1 z = sin ( y )   z ′ = u − − − − − − → z ( 0 )= sin ( 1 ) y ( 0 )= 1 u = cos ( y ) u ′ = − z u ( 0 )= cos ( 1 )   A. Pouly, D. Graça (LIX, FCT) Complexity of solving PIVP July 8, 2013 2 / 18

Introduction Existing results Practical Definition (Folklore) Numerical method: t i + 1 = t i + h and x i + 1 = f ( x 0 , . . . , x i ; h ) Local error: δ ih = � y ( t i ) − x i � ∞ Order: maximum ω such that δ h � h ω + 1 � n = O as h → 0 δ 3 y ( t ) x 3 δ 4 δ 2 x 2 x 4 x 1 x 0 t 0 t t 1 h t 2 t 3 t 4 A. Pouly, D. Graça (LIX, FCT) Complexity of solving PIVP July 8, 2013 3 / 18

Introduction Existing results Practical Definition (Folklore) Numerical method: t i + 1 = t i + h and x i + 1 = f ( x 0 , . . . , x i ; h ) Local error: δ ih = � y ( t i ) − x i � ∞ h ω + 1 � Order: maximum ω such that δ h � n = O as h → 0 Theorem (Folklore) Euler method has order 1 Runge-Kutta 4 (RK4) has order 4 ∀ ω , there exist methods of order ω (RK ω , Taylor) A. Pouly, D. Graça (LIX, FCT) Complexity of solving PIVP July 8, 2013 3 / 18

Introduction Existing results Practical Definition (Folklore) Numerical method: t i + 1 = t i + h and x i + 1 = f ( x 0 , . . . , x i ; h ) Local error: δ ih = � y ( t i ) − x i � ∞ h ω + 1 � Order: maximum ω such that δ h � n = O as h → 0 Theorem (Folklore) Euler method has order 1 Runge-Kutta 4 (RK4) has order 4 ∀ ω , there exist methods of order ω (RK ω , Taylor) A. Pouly, D. Graça (LIX, FCT) Complexity of solving PIVP July 8, 2013 3 / 18

Introduction Existing results Practical Definition (Folklore) Numerical method: t i + 1 = t i + h and x i + 1 = f ( x 0 , . . . , x i ; h ) Local error: δ ih = � y ( t i ) − x i � ∞ h ω + 1 � Order: maximum ω such that δ h � n = O as h → 0 Theorem (Folklore) Euler method has order 1 Runge-Kutta 4 (RK4) has order 4 ∀ ω , there exist methods of order ω (RK ω , Taylor) A. Pouly, D. Graça (LIX, FCT) Complexity of solving PIVP July 8, 2013 3 / 18

Introduction Existing results Practical Definition (Folklore) Numerical method: t i + 1 = t i + h and x i + 1 = f ( x 0 , . . . , x i ; h ) Local error: δ ih = � y ( t i ) − x i � ∞ Order: maximum ω such that δ h � h ω + 1 � n = O as h → 0 Theorem (Folklore) Euler method has order 1 Runge-Kutta 4 (RK4) has order 4 ∀ ω , there exist methods of order ω (RK ω , Taylor) Remark Difficult choice of h Quite efficient in practice A. Pouly, D. Graça (LIX, FCT) Complexity of solving PIVP July 8, 2013 3 / 18

Introduction Existing results Practical (Handwaving) Definition (Folklore) Adaptive method: t i + 1 = t i + h i and x i + 1 = f ( x 0 , . . . , x i ; h ) Local error: δ i = � y ( t i ) − x i � ∞ Error estimate: e i � δ i , → h i = g ( e i , x , t ) Idea Big steps when smooth and small error estimate Small steps when stiff and big error estimate Remark Unknown complexity Very efficient in practice A. Pouly, D. Graça (LIX, FCT) Complexity of solving PIVP July 8, 2013 3 / 18

Introduction Existing results And so ? Don’t we know everything ? A. Pouly, D. Graça (LIX, FCT) Complexity of solving PIVP July 8, 2013 4 / 18

Introduction Existing results And so ? Don’t we know everything ? Not quite! A. Pouly, D. Graça (LIX, FCT) Complexity of solving PIVP July 8, 2013 4 / 18

Introduction Existing results And so ? Don’t we know everything ? Not quite! � y ′ = p ( y ) y : I → R n where y ( t 0 )= y 0 p : vector of polynomials A. Pouly, D. Graça (LIX, FCT) Complexity of solving PIVP July 8, 2013 4 / 18

Introduction Existing results And so ? Don’t we know everything ? Not quite! � y ′ = p ( y ) y : I → R n where y ( t 0 )= y 0 p : vector of polynomials Issue #1: order ω , step size h � h ω + 1 � local error = O A. Pouly, D. Graça (LIX, FCT) Complexity of solving PIVP July 8, 2013 4 / 18

Introduction Existing results And so ? Don’t we know everything ? Not quite! � y ′ = p ( y ) y : I → R n where y ( t 0 )= y 0 p : vector of polynomials Issue #1: order ω , step size h local error � Kh ω + 1 A. Pouly, D. Graça (LIX, FCT) Complexity of solving PIVP July 8, 2013 4 / 18

Introduction Existing results And so ? Don’t we know everything ? Not quite! � y ′ = p ( y ) y : I → R n where y ( t 0 )= y 0 p : vector of polynomials Issue #1: order ω , step size h local error � Kh ω + 1 K depends on y and I !! A. Pouly, D. Graça (LIX, FCT) Complexity of solving PIVP July 8, 2013 4 / 18

Introduction Existing results And so ? Don’t we know everything ? Not quite! � y ′ = p ( y ) y : I → R n where y ( t 0 )= y 0 p : vector of polynomials Issue #1: order ω , step size h local error � Kh ω + 1 K depends on y and I !! Example: Euler method (Simplified) local error at step i � 1 2 h 2 � � p ′ ( y i ) � � ∞ A. Pouly, D. Graça (LIX, FCT) Complexity of solving PIVP July 8, 2013 4 / 18

Introduction Existing results And so ? Don’t we know everything ? Not quite! � y ′ = p ( y ) y : I → R n where y ( t 0 )= y 0 p : vector of polynomials Issue #1: order ω , step size h local error � Kh ω + 1 K depends on y and I !! Example: Euler method (Simplified) local error � 1 2 h 2 � � p ′ ( y i ) � p ′ ( y ( t )) � � � ∞ ⇒ O ( 1 ) = max ∞ ? � � t ∈ I A. Pouly, D. Graça (LIX, FCT) Complexity of solving PIVP July 8, 2013 4 / 18