function yd = f(y, t)
  yd(1) = y(2);
  yd(2) = -sin(y(1));
endfunction

x0 = [3, 0];

t = linspace(0, 20, 500);

x = lsode("f", x0, t);

plot(t,x);

