Modeling Physics with Differential-Algebraic Equations Lecture 3 Numerical Integration of DAEs COMASIC (M2) January 7th, 2019 Khalil Ghorbal k halil.ghorbal@inria.fr K. Ghorbal (INRIA) 1 COMASIC M2 1 / 22
Summary of lecture 2 Index Reduction • Given a DAE F ( x , ˙ x , t ), we have seen how to perform a structural analysis to numerically compute ˙ x function of x . The structural nonsingularity ensures that, generically, one can perform the computation following the block order suggested by the BLT decomposition. Thus, one is able to compute the numerical values of the derivatives given a consistent state of the system and carry on with a standard numerical integration . • Note: the structural index is not always equal to the differentiation index . K. Ghorbal (INRIA) 2 COMASIC M2 2 / 22
Summary of lecture 2 Index Reduction • Given a DAE F ( x , ˙ x , t ), we have seen how to perform a structural analysis to numerically compute ˙ x function of x . The structural nonsingularity ensures that, generically, one can perform the computation following the block order suggested by the BLT decomposition. Thus, one is able to compute the numerical values of the derivatives given a consistent state of the system and carry on with a standard numerical integration . • Note: the structural index is not always equal to the differentiation index . K. Ghorbal (INRIA) 2 COMASIC M2 2 / 22
Structural Analysis (Complement) Consider the following DAE z − ˙ xy − x ˙ y + 2 x + y − 3 = 0 ˙ z − xy = 0 x + y − 2 = 0 • Pantelides: structural index 1 • Differentiation index is 0 (simple linear system) • (Hidden) Cancellation problems are undecidable in general K. Ghorbal (INRIA) 3 COMASIC M2 3 / 22
Structural Analysis (Complement) Consider the following DAE z − ˙ xy − x ˙ y + 2 x + y − 3 = 0 ˙ z − xy = 0 x + y − 2 = 0 • Pantelides: structural index 1 • Differentiation index is 0 (simple linear system) • (Hidden) Cancellation problems are undecidable in general K. Ghorbal (INRIA) 3 COMASIC M2 3 / 22
Structural Analysis (Complement) Consider the following DAE z − ˙ xy − x ˙ y + 2 x + y − 3 = 0 ˙ z − xy = 0 x + y − 2 = 0 • Pantelides: structural index 1 • Differentiation index is 0 (simple linear system) • (Hidden) Cancellation problems are undecidable in general K. Ghorbal (INRIA) 3 COMASIC M2 3 / 22
Structural Analysis (Complement) Consider the following DAE z − ˙ xy − x ˙ y + 2 x + y − 3 = 0 ˙ z − xy = 0 x + y − 2 = 0 • Pantelides: structural index 1 • Differentiation index is 0 (simple linear system) • (Hidden) Cancellation problems are undecidable in general K. Ghorbal (INRIA) 3 COMASIC M2 3 / 22
Semi-Explicit DAE Index reduction transforms a fully implicit DAE to a semi-explicit DAE ˙ = f ( x , y , t ) x 0 = g ( x , y , t ) Integration Schemes • Backward Differentiation Formula (BDF) • Orthogonal Collocation K. Ghorbal (INRIA) 4 COMASIC M2 4 / 22
Outline 1 BDF Method 2 Collocation (Overview) 3 Brief Introduction to Modelica K. Ghorbal (INRIA) 4 COMASIC M2 4 / 22
Backward Differentiation Formula (BDF) Fixed Step Size • Implicit : next value not explicitly given. • Linear multistep : the next value is linearly related to the immediate previous (backward) values (eventually more than one). • For a fixed step size δ > 0, let t n = t 0 + n δ , and x s the approximate of the exact x ( t s ). • The general Backward Differentiation Formula: x n = linear combination of x n , x n − 1 , . . . , x 0 ˙ • For a Cauchy problem ˙ x = f ( x , t ) , x ( t 0 ) = x 0 , one obtains an implicit equation for x n : q � x n = a n x 1 − n + δ b q f ( x n , t n ) n =1 where the a n and b q depends only on the order q . • E.g., q = 1 (BDF1): x n = x n − 1 + δ f ( x n , t n ) (a.k.a. Backward Euler) K. Ghorbal (INRIA) 5 COMASIC M2 5 / 22
Backward Differentiation Formula (BDF) Fixed Step Size • Implicit : next value not explicitly given. • Linear multistep : the next value is linearly related to the immediate previous (backward) values (eventually more than one). • For a fixed step size δ > 0, let t n = t 0 + n δ , and x s the approximate of the exact x ( t s ). • The general Backward Differentiation Formula: x n = linear combination of x n , x n − 1 , . . . , x 0 ˙ • For a Cauchy problem ˙ x = f ( x , t ) , x ( t 0 ) = x 0 , one obtains an implicit equation for x n : q � x n = a n x 1 − n + δ b q f ( x n , t n ) n =1 where the a n and b q depends only on the order q . • E.g., q = 1 (BDF1): x n = x n − 1 + δ f ( x n , t n ) (a.k.a. Backward Euler) K. Ghorbal (INRIA) 5 COMASIC M2 5 / 22
Backward Differentiation Formula (BDF) Fixed Step Size • Implicit : next value not explicitly given. • Linear multistep : the next value is linearly related to the immediate previous (backward) values (eventually more than one). • For a fixed step size δ > 0, let t n = t 0 + n δ , and x s the approximate of the exact x ( t s ). • The general Backward Differentiation Formula: x n = linear combination of x n , x n − 1 , . . . , x 0 ˙ • For a Cauchy problem ˙ x = f ( x , t ) , x ( t 0 ) = x 0 , one obtains an implicit equation for x n : q � x n = a n x 1 − n + δ b q f ( x n , t n ) n =1 where the a n and b q depends only on the order q . • E.g., q = 1 (BDF1): x n = x n − 1 + δ f ( x n , t n ) (a.k.a. Backward Euler) K. Ghorbal (INRIA) 5 COMASIC M2 5 / 22
Backward Differentiation Formula (BDF) Fixed Step Size • Implicit : next value not explicitly given. • Linear multistep : the next value is linearly related to the immediate previous (backward) values (eventually more than one). • For a fixed step size δ > 0, let t n = t 0 + n δ , and x s the approximate of the exact x ( t s ). • The general Backward Differentiation Formula: x n = linear combination of x n , x n − 1 , . . . , x 0 ˙ • For a Cauchy problem ˙ x = f ( x , t ) , x ( t 0 ) = x 0 , one obtains an implicit equation for x n : q � x n = a n x 1 − n + δ b q f ( x n , t n ) n =1 where the a n and b q depends only on the order q . • E.g., q = 1 (BDF1): x n = x n − 1 + δ f ( x n , t n ) (a.k.a. Backward Euler) K. Ghorbal (INRIA) 5 COMASIC M2 5 / 22
Backward Differentiation Formula (BDF) Fixed Step Size • Implicit : next value not explicitly given. • Linear multistep : the next value is linearly related to the immediate previous (backward) values (eventually more than one). • For a fixed step size δ > 0, let t n = t 0 + n δ , and x s the approximate of the exact x ( t s ). • The general Backward Differentiation Formula: x n = linear combination of x n , x n − 1 , . . . , x 0 ˙ • For a Cauchy problem ˙ x = f ( x , t ) , x ( t 0 ) = x 0 , one obtains an implicit equation for x n : q � x n = a n x 1 − n + δ b q f ( x n , t n ) n =1 where the a n and b q depends only on the order q . • E.g., q = 1 (BDF1): x n = x n − 1 + δ f ( x n , t n ) (a.k.a. Backward Euler) K. Ghorbal (INRIA) 5 COMASIC M2 5 / 22
Coefficients of BDFq Taylor Series (Sundials Implementation) I x n = x n (Identity) N x n = x n +1 (Forward Shift) N − 1 x n = x n − 1 (Backward Shift) D x n = ˙ x n (Differential) = I − N ∆ (Backward Operator) Observe that N = ( I − ∆) − 1 (Operator Algebra) N x n = x n +1 = x n + δ D x n + δ 2 2! D 2 x n + δ 3 3! D 3 x n + · · · = ( I + δ D + δ 2 2! D 2 + δ 3 3! D 3 x n + · · · ) x n = e δ D x n Thus: N = e δ D K. Ghorbal (INRIA) 6 COMASIC M2 6 / 22
= ∆ + 1 2∆ 2 + 1 3∆ 3 + · · · ( I − ∆) − 1 � � δ D = ln( N ) = ln x n = δ D x n = ∆ x n + 1 2∆ 2 x n + 1 3∆ 3 x n + · · · δ ˙ BDFq: truncate at order q , for instance for q = 2 x n = ∆ x n + 1 2∆ 2 x n = ( x n − x n − 1 ) + 1 2( x n − 2 x n − 1 + x n − 2 ) δ ˙ Thus x n = 3 2 x n − 2 x n − 1 + 1 δ ˙ 2 x n − 2 x n = 4 3 x n − 1 − 1 3 x n − 2 + 2 3 δ f ( x n , t n ) K. Ghorbal (INRIA) 7 COMASIC M2 7 / 22
BDF for Semi-Explicit DAE of index 1 x = f ( x , y , t ) ˙ 0 = g ( x , y ) Numerical integration using the BDF2 Scheme Suppose x n − 1 and x n − 2 are known, x n and y n are computed by numerically solving the following system (e.g. Newton’s methods) at each iteration: x n = 4 3 x n − 1 − 1 3 x n − 2 + 2 3 δ f ( x n , y n , t n ) 0 = g ( x n , y n ) • BDF converges if m ≤ 6: x i − x ( t i ) = y i − y i ( t i ) = o ( δ m ) • BDF requires a consistent initial condition K. Ghorbal (INRIA) 8 COMASIC M2 8 / 22
BDF for Semi-Explicit DAE of index 1 x = f ( x , y , t ) ˙ 0 = g ( x , y ) Numerical integration using the BDF2 Scheme Suppose x n − 1 and x n − 2 are known, x n and y n are computed by numerically solving the following system (e.g. Newton’s methods) at each iteration: x n = 4 3 x n − 1 − 1 3 x n − 2 + 2 3 δ f ( x n , y n , t n ) 0 = g ( x n , y n ) • BDF converges if m ≤ 6: x i − x ( t i ) = y i − y i ( t i ) = o ( δ m ) • BDF requires a consistent initial condition K. Ghorbal (INRIA) 8 COMASIC M2 8 / 22
Newton’s Methods: Generalization Could be used to numerically approximate the roots of F : R k → R k ( F i continuously differentiable functions): • Multiply by the inverse of the Jacobian x n +1 = x n − J − 1 F ( x n ) F ( x n ) • Or solve the system of linear equations J F ( x n )( x n +1 − x n ) = − F ( x n ) Under some assumptions, the method converges quadratically towards a root of F . K. Ghorbal (INRIA) 9 COMASIC M2 9 / 22
Recommend
More recommend
