Reducing Uncertainty when a Approximating Solutions of ODEs W.H. Enright Department of Computer Science University of Toronto IFIP Working Conference Uncertainty Quantification in Scientific Computing Boulder, Colorado August 2011 a This research is supported by the Natural Science and Engineering Research Council of Canada. Reducing Uncertainty when Solving ODEs – p.1/30
The Focus of this Talk Tools for Reducing Uncertainty when Investigating Mathematical Models described by Systems of ODEs “Investigating” – not only “Approximating the Solution”. 1. Sensitivity Analysis - (of solution wrt parameters defining the problem) 2. Global Error Estimation 3. Estimating the "Conditioning" of a problem “ODEs” includes IVPs, BVPs, DDEs, DAEs and VIDEs. Reducing Uncertainty when Solving ODEs – p.2/30
Acknowledgement This work is part of an ongoing project and has benefited from numerous discussions and collaborations with Paul Muir Ken Jackson Silvana Ilie Hossein Zivari Piran Kante Easley Mohammad Shakourifar Li Yan William Zhang Bo Wang Reducing Uncertainty when Solving ODEs – p.3/30
An Effective ODE Solver Minimum Requirements: An Accurate Discrete Approximation is not Enough 120 22 20 100 18 80 16 14 60 12 40 10 8 20 6 0 4 0 5 10 15 20 25 30 35 40 10 20 30 40 50 60 70 80 90 100 110 An Accurate Continuous Extension is Necessary (aka Dense Output) 120 22 20 100 18 16 80 14 60 12 10 40 8 6 20 4 0 2 0 5 10 15 20 25 30 35 40 10 20 30 40 50 60 70 80 90 100 110 Reducing Uncertainty when Solving ODEs – p.4/30
Outline of Talk Current Scientific Computing Paradigm and its implications: Acceptability of an approximate solution Continuous RK Methods provide dense output for ODEs Defect Error Control for CRK Methods Measuring or quantifying the Reliability of a CRK Method Classes of ODE problems that can be Investigated by CRK-based Methods (IVPs, BVPs, DDEs, DAEs, and VIDEs) Useful Software Tools for Investigating and quantifying Important Properties of the Mathematical Model and its Approximate Solution. (sensitivity analysis, global error estimation, parameter fitting and condition number estimation.) Some Numerical Examples What is Next? Reducing Uncertainty when Solving ODEs – p.5/30
Scientific Computing Paradigm Mathematical Modelling in a Problem Solving Environment: Formulate the mathematical model of the system being investigated. (The model may involve parameters.) Approximate the exact solution of this model relative to a specified accuracy parameter, TOL . Visualize the approximate solution. Is the mathematical model well-posed and is the approximate solution stable? (may involve sensitivity analysis). Reducing Uncertainty when Solving ODEs – p.6/30
Implications for ODE Methods What is an acceptable approximate solution? An approximate solution must be easy to represent, display and manipulate. The accuracy (or quality ) of an approximate solution must be easy to measure and interpret. What are the implications for an ODE method? Solver must be easy to invoke –(only need specify those parameters necessary to define the problem). A discrete solution is not sufficient (as it is difficult to visualize and its accuracy is difficult to interpret). It should have a generic calling sequence so it is easy to adopt in a PSE. Reducing Uncertainty when Solving ODEs – p.7/30
Continuous Runge-Kutta Methods Consider an IVP defined by the system y ′ = f ( x, y ) , y ( a ) = y 0 , on [ a, b ] . A numerical method will introduce a partitioning a = x 0 < x 1 < · · · < x N = b and corresponding discrete approximations y 0 , y 1 · · · y N . The y i ’s are usually determined sequentially. On step i let z i ( x ) be the solution of the local IVP: z ′ i = f ( x, z i ( x )) , z i ( x i − 1 ) = y i − 1 , on [ x i − 1 , x i ] . Reducing Uncertainty when Solving ODEs – p.8/30
CRK methods (cont) A classical p th -order, s-stage, discrete RK formula determines s � y i = y i − 1 + h i ω j k j , j =1 where h i = x i − x i − 1 and the j th stage is defined by, s � k j = f ( x i − 1 + h i c j , y i − 1 + h i a jr k r ) . r =1 A Continuous extension (CRK) is determined by introducing (˜ s − s ) additional stages to obtain an order p approximation for any x ∈ ( x i − 1 , x i ) ˜ s b j ( x − x i − 1 � u i ( x ) = y i − 1 + h i ) k j , h i j =1 where b j ( τ ) is a polynomial of degree at least p and τ = x − x i − 1 . h i Reducing Uncertainty when Solving ODEs – p.9/30
CRK methods (cont) We consider O ( h p ) extensions, satisfying: ˜ s b j ( τ ) k j = z i ( x ) + O ( h p +1 � u i ( x ) = y i − 1 + h i ) . i j =1 The [ u i ( x )] N i =1 define a piecewise polynomial U ( x ) for x ∈ [ x 0 , x F ] . This is the approximate solution generated by the CRK method. U ( x ) ∈ C 0 [ x 0 , x F ] and will interpolate the underlying discrete RK values, y i , if b j (1) = ω j for j = 1 , 2 · · · s and b s +1 (1) = b s +2 (1) = · · · b ˜ s (1) = 0 . d d τ ( b j ( τ )) , and requiring Similarly a simple set of constraints on the that k s +1 = f ( x i , y i ) , k 1 = f ( x i − 1 , y i − 1 ) , will ensure U ′ ( x ) interpolates f ( x i , y i ) , f ( x i − 1 , y i − 1 ) and therefore U ( x ) ∈ C 1 [ x 0 , x F ] . Reducing Uncertainty when Solving ODEs – p.10/30
Defect Error Control for CRKs U ( x ) , the approximate solution, has an associated defect or residual, δ ( x ) ≡ f ( x, U ( x )) − U ′ ( x ) ≡ f ( x, u i ( x )) − u ′ i ( x ) , for x ∈ [ x i − 1 , x i ] . It can be shown that, for sufficiently differentiable f , δ ( x ) = G ( τ ) h p i + O ( h p +1 ) , i G ( τ ) = ˜ q 1 ( τ ) F 1 + ˜ q 2 ( τ ) F 2 + · · · + ˜ q k ( τ ) F k , where the ˜ q j are polynomials in τ that depend only on the CRK formula while the F j are constants (elementary differentials) that depend only on the problem. Methods can be implemented to adjust h i in an attempt to ensure that the maximum magnitude of δ ( x ) is bounded by TOL . The quality of an approx- imate solution can then be described in terms of the max of || δ ( x ) || /TOL . Reducing Uncertainty when Solving ODEs – p.11/30
Defect Error Control (cont) δ ( x ) = G ( τ ) h p i + O ( h p +1 ) , i G ( τ ) = ˜ q 1 ( τ ) F 1 + ˜ q 2 ( τ ) F 2 + · · · + ˜ q k ( τ ) F k . As h i → 0 the defect will then look like a linear combination of the known polynomials ˜ q j ( τ ) over [ x i − 1 , x i ] . In the special case where k = 1 the shape of the defect will be the same (as h i → 0 ) for all problems and all steps. That is, the defect will almost always ’converge’ to a multiple of ˜ q 1 ( τ ) , in which case the maximum should occur (as h i → 0 ) at τ = τ ∗ where τ ∗ is the location of the local extremum of ˜ q 1 ( τ ) . In this case we will refer to the defect control strategy as Strict Defect Control (SDC) . Reducing Uncertainty when Solving ODEs – p.12/30
Typical Shape of SDC Defects 1 0.8 0.6 defect 0.4 0.2 0 −0.2 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 δ ( τ ) /δ ( τ ∗ ) vs τ ) for each Plot of scaled defect vs τ (ie. step required to solve a typical problem with SDC CRK6 and TOL = 10 − 6 . Reducing Uncertainty when Solving ODEs – p.13/30
Cost of Strict Defect Control s p th − order, explicit, discrete RK : y i = y i − 1 + h i � ω j k j , j =1 ˜ s b j ( τ ) k j = z i ( x ) + O ( h p +1 ˜ � SDC : u i ( x ) = y i − 1 + h i ˜ ) . i j =1 Formula p s s ˜ CRK4 4 4 8 CRK5 5 6 12 CRK6 6 7 15 CRK7 7 9 20 CRK8 8 13 27 Table 1: Cost per step of some SDC-CRK formulas (Note that ˜ s ≈ 2 s .) Reducing Uncertainty when Solving ODEs – p.14/30
Strict Defect Control SDC CRKs are not unique (for a given discrete RK formula). Each SDC-CRK satisfies, q 1 ( τ ) ˆ q 2 ( τ ) ˆ k ( τ ) ˆ k ) h p +1 + O ( h p +2 q 1 ( τ ) F 1 h p δ ( x ) = ˜ i + (ˆ F 1 + ˆ F 2 + · · · · · · ˆ q ˆ F ˆ ) i i Potential Difficulties: q 1 ( τ ) may have a large maximum ( ˜ ˜ q 1 (0) = ˜ q 1 (1) = 0 and its ‘average’ value must be one). The ˆ q j ( τ ) may be large in magnitude relative to ˜ q 1 ( τ ) (and therefore h i would have to be small before the estimate is justified). (That is, before | h i ˆ q j ( τ ) | << | ˜ q 1 ( τ ) | .) | F 1 | may be zero (or very small) on isolated steps. For each p we have identified a particular SDC-CRK that minimizes these difficulties. Reducing Uncertainty when Solving ODEs – p.15/30
Optimal SDC CRK6 2 1 0 0.0 0.25 0.5 0.75 1.0 x −1 Figure 1: Plots of ˜ q 1 and ˆ q 7 for SDC CRK6. ˜ q 1 is represented q 2 · · · ˆ by the solid line and has the highest magnitude. Reducing Uncertainty when Solving ODEs – p.16/30
Optimal SDC CRK8 3 2 1 0 0.0 0.25 0.5 0.75 1.0 x −1 −2 Figure 2: Plots of ˜ q 1 and ˆ q 9 for SDC CRK8. ˜ q 1 is represented q 2 · · · ˆ by the solid line and has the highest magnitude. Reducing Uncertainty when Solving ODEs – p.17/30
Recommend
More recommend
