Off[General::spell1];
Clear[ABM,F0,a,b,x0,n];
ABM[F0_,a_,b_,x0_,n_] :=
Module[{f,f0,f1,f2,f3,f4,h,h2,j,k, X, T},
  Set@@{f[t_,x_],F0};
  h = (b - a)/n //N;
  T = Table[0,{n+1}];
  X = Table[0,{n+1}];
  T[[1]] = a //N;
  X[[1]] = x0 //N;
  For[j=1, j<=3,j++,
    F1 = h f[T[[j]],X[[j]]] //N;
    F2 = h f[T[[j]]+ h/2, X[[j]] + F1/2] //N;
    F3 = h f[T[[j]] + h/2, X[[j]]+ F2/2] //N;
    F4 = h f[T[[j]] + h,X[[j]] + F3] //N;
    X[[j+1]] = X[[j]] + (F1 + 2F2 + 2F3 + F4)/6 //N;
    T[[j+1]] = a + h j //N];
  f0 = f[T[[1]],X[[1]]];
  f1 = f[T[[2]],X[[2]]];
  f2 = f[T[[3]],X[[3]]];
  f3 = f[T[[4]],X[[4]]];
  h2 = h/24;
  For[k=4,k<=n,k++,
    P = X[[k]] + h2 (-9f0 + 37f1 - 59f2 + 55f3);
    T[[k+1]] = a + h k;
    f4 = f[T[[k+1]],P];
    X[[k+1]] = X[[k]] + h2 (f1 - 5f2 + 19f3 + 9f4);
    f0 = f1;
    f1 = f2;
    f2 = f3;
    f3 = f[T[[k+1]],X[[k+1]]] ];
  sol = Transpose[{T,X}];
  Return[sol] ];
On[General::spell1];

Clear[f,t,x];   
f[t_,x_] = Cos[t] - Sin[x] +t^2;

ABM[f[t,x],-1,1,3,10]

Clear[pts];
pts = ABM[f[t,x],-1,1,3,200];
ListPlot[pts,PlotJoined->True, AspectRatio -> 0.5];

Clear[ss];
ss=NDSolve[{x'[t]==f[t,x[t]],x[-1]==3},x[t],{t,-1,1}];
Plot[Evaluate[x[t]/.ss],{t,-1,1},AspectRatio -> 0.5];

Clear[xx];
xx[t_]= x[t] /. First[ss];
Error = Table[pts[[i,2]]-xx[-1+ (i-1) /100], {i,1,10}]

Error = Table[pts[[i,2]]-xx[-1+ (i-1) /100], {i,190,200}]

Clear[f,t,x];
f[t_,x_] = (t *x - x^2) / t^2;

Clear[pts];
pts = ABM[f[t,x],1,3,2,200 ];
ListPlot[pts,PlotJoined->True, AspectRatio -> 1];

Clear[ss,xx];
ss=DSolve[{x'[t]==f[t,x[t]],x[1]==2},x[t],t];
xx[t_]= x[t] /. First[ss]

Plot[xx[t], {t, 1,3},AspectRatio -> 1]

Error1 = Table[pts[[i,2]]-xx[1+ (i-1) /100], {i,1,10}]//N

Error2=Table[pts[[i,2]]-xx[1+ (i-1) /100], {i,90,100}]//N