function yp = duffing(t, y)

k = 0.25;
A = 0.4;
omega = 1.0;

yp = [y(2);
      -k * y(2) + y(1) * (1 - y(1)^2) + A * cos(omega * t)];
  