Runge-Kutta Method Runge-Kutta Method x n+ = x n +hf(t n x n ) dx dt dx x(t+h) x(t) (t) = lim dt h h t = t n h > dx dt (t x(t n +h) x(t n ) n) = lim x(t n+) x(t n ) h h h x = f(tx) x(t n+ ) x(t n ) f(t n x n ) h x(t n+ ) x(t n )+hf(t n x(t n )) x n+ = x n +hf(t n x n ) x n (explicit method) x n x n+ Jun 6 Version: naito@mathnagoya-uacjp
t = t n t = t n+ dx dt (t x(t n +h) x(t n ) n+) = lim x(t n+) x(t n ) h h h x = f(tx) x(t n+ ) x(t n ) f(t n+ x n+ ) h x(t n+ ) x(t n )+hf(t n+ x(t n+ )) x n+ = x n +hf(t n+ x n+ ) (backward Euler method) x n x n+ x n (implicit method) x (t) = x(t) x n+ = x n +hx n+ x n x n+ x n+ x n+ = h x n x (t) = sin(x(t)) x n+ = x n +hsin(x n+ ) x n x n+ x n+ y = x n y k+ = x n +hsin(y k ) x n+ h Jun 6 Version: naito@mathnagoya-uacjp
x (t) = f(x(t)) x (t) = f(tx(t)) X(t) = X(t) x = f(x(t)) x x(t n+ ) = x(t n )+hx (t n )+ h x (t n )+O(h ) = x(t n )+hf(x(t n ))+ h f(x (t n ))f(x(t n ))+O(h ) x n+ = x n +hf(x n ) x(t n+ ) x n+ = x(t n ) x n +h(f(x(t n )) f(x n ))+ h f (x(t n ))f(x(t n ))+O(h ) f f(p) f(q) L p q f x(t n+ ) x n+ (+Lh) x(t n ) x n +Mh x(t n ) = x(t n+ ) hx (t n+ )+ h x (t n+ )+O(h ) = x(t n+ ) hf(x(t n+ ))+ h f(x (t n+ ))f(x(t n+ ))+O(h ) x n+ = x n +hf(x n+ ) x(t n+ ) x n+ = x(t n ) x n +h(f(x(t n+ )) f(x n+ ))+ h f (x(t n+ ))f(x(t n+ ))+O(h ) x(t n+ ) x n+ ( Lh) x(t n ) x n +Mh ǫ n+ = Aǫ n +Mh ǫ n = x(t n ) x n t = t n x n x(t n ) O(h ) Jun 6 Version: naito@mathnagoya-uacjp
4 t = T = Nh x(t) x N Che TL x(t) x N Ch h h h x(t) x N = O(h p ) p > x(t n+ ) = x(t n )+hx (t n )+ h x (t n )+O(h ) = x(t n )+hf(x(t n ))+ h f(x (t n ))f(x(t n ))+O(h ) h k k = x n+ = x n +hf(x n )+ h f (x n )f(x n ) O(h ) (h k O(h k+ ) k f x(t n+ ) = x(t n )+hx (t n )+ h x (t n )+ h 6 x() (t)+o(h 4 ) x(t n ) = x(t n ) hx (t n )+ h x (t n ) h 6 x() (t)+o(h 4 ) x(t n+ ) x(t n ) = hx (t n )+ h x() (t)+o(h 5 ) x n+ = x n +hf(x n ) O(h ) Jun 6 Version: naito@mathnagoya-uacjp
5 x (t) = x(t) x() = x n+ = x n hx n x = x = +h {x n } x n n x(t) = e t x(t) (t ) x(t) x n x n+ = x n hx n t +ht = α α α < < α < x n = A α n +A α n x n A = x = C x = Cα A x n a n+ = a n +a n Jun 6 Version: naito@mathnagoya-uacjp
6 h = x = x 5 5 5 4 6 8 x = -x 4 8 6 4 4 6 8 x (t) = x(t) x() = x (t) = x(t) x() = x = (-x^) 8 6 4-5 5 5 5 4 x = sin(x) 4 5 5 5 5 4 6 8 x (t) = ( x(t) ) x() = x (t) = sin(x(t)) x() = x = -x 5 5-5 - -5-4 6 8 4 6 x = -x 5 5-5 - -5 - - -5 - -5 5 5 x (t) = x(t) x() = x () = Jun 6 Version: naito@mathnagoya-uacjp
7 x = -sin(x) - - - 4 6 8 4 6 - - x = -sin(x) - - x (t) = sin(x(t)) x() = x () = x = (-x^)x - x - - - 5 5 5 x = (-x^)x - x - - - - - - x (t) = ( x(t) )x (t) x(t) x() = 5 x () = X = -X/ X ^ 5-5 - - -5 5 X (t) = X(t)/ X(t) X() = () X () = () Double Pendulum 5 5-5 - -5 - -5 4 6 8 Double Pendulum 5 5-5 - -5 - -5 - -5 - -5 5 5 x () = x () = x () = x () = Jun 6 Version: naito@mathnagoya-uacjp
8 h = x = x 5 5 5 4 6 8 x = -x 4 8 6 4 4 6 8 x (t) = x(t) x() = x (t) = x(t) x() = x = (-x^) 8 6 4-5 5 5 5 4 x = sin(x) 4 5 5 5 5 4 6 8 x (t) = ( x(t) ) x() = x (t) = sin(x(t)) x() = x = -x 5 5-5 - -5-4 6 8 4 6 x = -x 5 5-5 - -5 - - -5 - -5 5 5 x (t) = x(t) x() = x () = Jun 6 Version: naito@mathnagoya-uacjp
9 x = -sin(x) - - - 4 6 8 4 6 - - x = -sin(x) - - x (t) = sin(x(t)) x() = x () = x = (-x^)x - x - - - 5 5 5 x = (-x^)x - x - - - - - - x (t) = ( x(t) )x (t) x(t) x() = 5 x () = X = -X/ X ^ 5-5 - - -5 5 X (t) = X(t)/ X(t) X() = () X () = () Double Pendulum 5 5-5 - -5 - -5 4 6 8 Double Pendulum 5 5-5 - -5 - -5 - -5 - -5 5 5 x () = x () = x () = x () = Jun 6 Version: naito@mathnagoya-uacjp
h = x = x 5 5 5 4 6 8 x = -x 4 8 6 4 4 6 8 x (t) = x(t) x() = x (t) = x(t) x() = x = (-x^) 8 6 4-5 5 5 5 4 x = sin(x) 4 5 5 5 5 4 6 8 x (t) = ( x(t) ) x() = x (t) = sin(x(t)) x() = x = -x 5 5-5 - -5-4 6 8 4 6 x = -x 5 5-5 - -5 - - -5 - -5 5 5 x (t) = x(t) x() = x () = Jun 6 Version: naito@mathnagoya-uacjp
x = -sin(x) - - - 4 6 8 4 6 - - x = -sin(x) - - x (t) = sin(x(t)) x() = x () = x = (-x^)x - x - - - 5 5 5 x = (-x^)x - x - - - - - - x (t) = ( x(t) )x (t) x(t) x() = 5 x () = X = -X/ X ^ 5-5 - - -5 5 X (t) = X(t)/ X(t) X() = () X () = () Double Pendulum 5 5-5 - -5 - -5 4 6 8 Double Pendulum 5 5-5 - -5 - -5 - -5 - -5 5 5 x () = x () = x () = x () = Jun 6 Version: naito@mathnagoya-uacjp
x = ±x x() = t = t = x = x relative error at t = x = -x relative error at t = Forward Euler Backward Euler Improved Euler Heun Heun Kutta Runge-Kutta Runge-Kutta-Gill e-6 Forward Euler Backward Euler Improved Euler Heun Heun Kutta Runge-Kutta Runge-Kutta-Gill e-6 e-8 e-8 e- e- e- e- e-4 e-4 e-6 e+6 e+7 e+8 e+9 e-6 e+6 e+7 e+8 e+9 x = x relative error at t = x = -x relative error at t = Forward Euler Backward Euler Improved Euler Heun Heun Kutta Runge-Kutta Runge-Kutta-Gill e-6 Forward Euler Backward Euler Improved Euler Heun Heun Kutta Runge-Kutta Runge-Kutta-Gill e-6 e-8 e-8 e- e- e- e- e-4 e-4 e-6 e-6 e+6 e+7 e+8 e+9 e-8 e+6 e+7 e+8 e+9 x (t) = x(t) x() = x (t) = x(t) x() = h 9 Jun 6 Version: naito@mathnagoya-uacjp