Off[General::spell1];
Clear[RungeKutta2D,f0,a,b,X0,n];
RungeKutta2D[f0_,a_,b_,X0_,n_] :=
Module[{f,h,j,F1,F2,F3,F4,T,X,Y},
  Set@@{f[t_,{x_,y_}],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<=n,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 + 2 F2 + 2 F3 + F4)/6 //N;
    T[[j+1]] = a + h j //N];
 Return[Transpose[{T,X}]] ];
On[General::spell1];

Clear[f,f1,f2,t,x];   
f1[t_,{x_,y_}] = x+ 4 y - Exp[t];
f2[t_,{x_,y_}] = x+ y + 2 Exp[t];
f[t_,{x_,y_}] ={f1[t,{x,y}],f2[t,{x,y}]}

RungeKutta2D[f[t,{x,y}], 0,1,{4,5/4},10]

Clear[RK,RKX,RKY];
 RK= RungeKutta2D[f[t,{x,y}], 0,1,{4,5/4},100];  
 RKX=Table[{RK[[i,1]],RK[[i,2,1]]}, {i,1,101}];
 RKY=Table[{RK[[i,1]],RK[[i,2,2]]} ,{i,1,101}];

Clear[graph1, graph2];
 graph1=ListPlot[RKX, PlotJoined->True,AspectRatio -> 0.5,
             DisplayFunction->Identity]; 
 graph2=ListPlot[RKY,PlotJoined->True,AspectRatio -> 0.5, 
             DisplayFunction->Identity]; 
 Show[graph1,graph2,DisplayFunction->$DisplayFunction];

ss=DSolve[{x'[t]==x[t]+4 y[t]-Exp[t], 
           y'[t]== y[t]+x[t]+ 2 Exp[t],
           x[0]==4, y[0]== 5/4}, {x[t], y[t]},t] 
xx[t_]= x[t] /. First[ss];  
yy[t_]= y[t] /. Last[ss] ;

Error1 = Table[RKX[[i,2]]-xx[(i-1)/100] //N, {i,1,101}];
Max[Abs[%]]

Error2= Table[RKY[[i,2]]-yy[ (i-1)/100]//N, {i,1,101}];
Max[Abs[%]]

Clear[f,f1,f2,t,x];   
f1[t_,{x_,y_}] = x+y^2 -t^3;
f2[t_,{x_,y_}] = y+x^3 + Cos[t];
f[t_,{x_,y_}] ={f1[t,{x,y}],f2[t,{x,y}]}

RK= RungeKutta2D[f[t,{x,y}], 1,-2,{3,1},10]

Clear[RK,RKX,RKY];
 RK= RungeKutta2D[f[t,{x,y}], 1,-2,{3,1},150];  
 RKX=Table[{RK[[i,1]],RK[[i,2,1]]}, {i,1,151}];
 RKY=Table[{RK[[i,1]],RK[[i,2,2]]} ,{i,1,151}];

Clear[graph1, graph2];
 graph1=ListPlot[RKX, PlotJoined->True,AspectRatio -> 0.5,
             DisplayFunction->Identity]; 
 graph2=ListPlot[RKY,PlotJoined->True,AspectRatio -> 0.5, 
             DisplayFunction->Identity]; 
 Show[graph1,graph2,DisplayFunction->$DisplayFu nction];

ss=NDSolve[{x'[t]==x[t]+y[t]^2 -t^3, y'[t]== y[t]+x[t]^3 + Cos[t], 
    x[1]==3, y[1]== 1}, {x[t], y[t]},{t,-2,1}];
xx[t_]= x[t] /. First[ss];  
yy[t_]= y[t] /. Last[ss] ;

Clear[solu1, solu2];
 solu1= Plot[xx[t], {t,-1.95,1}, AspectRatio -> 0.5, 
             DisplayFunction->Identity]; 
 solu2= Plot[yy[t], {t,-1.95,1}, AspectRatio -> 0.5, 
             DisplayFunction->Identity]; 
 Show[solu1,solu2,DisplayFunction->$DisplayFunction]

Clear[f,f1,f2,t,x];   
f1[t_,{x_,y_}] = y; 
f2[t_,{x_,y_}] = (-5 y -5 x)/2 ;
f[t_,{x_,y_}] ={f1[t,{x,y}],f2[t,{x,y}]}

RungeKutta2D[f[t,{x,y}], 0,Pi,{0,1/2},10]

Clear[RK,RKX,RKY];
 RK= RungeKutta2D[f[t,{x,y}], 0,Pi,{0,1/2},100];  
 RKX=Table[{RK[[i,1]],RK[[i,2,1]]}, {i,1,101}];

Clear[graph1, graph2];
 ListPlot[RKX, PlotJoined->True,AspectRatio ->0.25]

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

Plot[xx[t], {t,0, Pi}, AspectRatio-> 0.25 ]

Error1 = Table[RKX[[i,2]]-xx[Pi(i-1)/100] //N, {i,1,101}];
Max[Abs[%]]