常微分方程式をRで解く際の精度

$\frac{dy}{dx} = -\frac{\sin(x)}{y}, y(0) = 1$
という常微分方程式の解は
$y(x) = \sqrt{2\cos(x) -1}$
で、これをRで数値的に確かめたい。
ついでに、数値積分の手法（ここではEuler法とRunge-Kutta4次）によって、どのくらい厳密解とずれるのかも調べる。

Rで常微分方程式を解く際にはdeSolveパッケージを用いるのが楽。結果を算出するためのコードは以下。

library(deSolve)
#厳密解
exact <- function(t){sqrt(2*cos(t)-1)}
#deSolveパッケージのode関数に食わせるための常微分方程式
f <- function(t, y, p){list(-sin(t)/y)}
#始状態から順方向で＋終状態から逆方向に＋厳密解の結果を同時に計算する関数
solve <- function(method)
{
  t1 <- seq(0, 1, by = 0.05)
  t2 <- seq(1, 0, by = -0.05)
  rbind(
    data.frame(ode(c(y=exact(0)), t1, f, c(), method), type="F"),
    data.frame(ode(c(y=exact(1)), t2, f, c(), method), type="B"),
    data.frame(time=t1, y=exact(t1),type="E")
  )
}

結果をggplot2でEuler法 v.s. Runge-Kutta 4次で対比できるように描画すると・・・

library(ggplot2)
ggplot(solve("euler"), aes(x=time, y=y)) +
  geom_path(size=2, aes(color=type)) +
  ggtitle("Differential Equation Solve(Euler)")

ggplot(solve("rk4"), aes(x=time, y=y)) +
  geom_path(size=2, aes(color=type)) +
  ggtitle("Differential Equation Solve(Runge-Kutta 4th order)")