deff('y=f(t,x)','y=-150*x+30')
//
N=evstr(x_dialog('nombre de points N','50'));
eps=evstr(x_dialog('perturbation epsilon','1E-3'));
//
T=1;
dt=T/N;t=0:dt:T;
//////////// Euler explicite//////
//
timer();
y0=1/5+eps;
yeu=y0;
for i=1:N
  y0=y0+dt*f(t(i),y0);
  yeu=[yeu,y0]  
end
//
disp('temps euler ex:');disp(timer());
xset('window',0)
clf()
plot2d(t,yeu,rect=[0,0,1,0.5])
//
///////////////////////////////////
// Euler implicite
/////////////////////////////////
y0=1/5+eps;
yeuimp=y0;
timer();
for i=1:N
  deff('y=g(x)','y=x-y0-dt*f(t(i),x)');
  y0=fsolve(y0,g); 
  yeuimp=[yeuimp,y0]  
end
disp('temps euler imp:');disp(timer());
//
plot2d(t,yeuimp,2,rect=[0,0,1,0.5])
xtitle('schema euler explicite (noir) vs implicite (bleu)');
