function [y] = vdp(t, x) c = 0.4; y(1) = x(2); y(2) = c * (1 - x(1) * x(1)) * x(2) - x(1); endfunction