function tracevpd(a, t0, T, Neuler, Nptmil, Nrk4, Node)
    // Approximation euler
    u = euler(a, t0, T, Neuler, vdp);
    
    // Approximation point milieu
    v = pointmilieu(a, t0, T, Nptmil, vdp);
    
    // Approximation runge kutta 4
    w = rk4(a, t0, T, Nrk4, vdp);
    
    // Approximation ode
    tode = linspace(t0, T, Node+1);
    y = ode(a, t0, tode, vdp)
    
    // Euler
    plot(u(1,:), u(2,:), '+');
    
    // Pt milieu
    plot(v(1,:), v(2,:), 'yo');
    
    // RK4
    plot(w(1,:), w(2,:), 'black.');
    
    // ODE
    plot(y(1,:), y(2,:), 'rx');
    
    title('VDP');
endfunction