(* DIRECTION FIELDS *)
(* Generate a direction field. Instructions from:
http://webpages.math.luc.edu/~lauve/courses/132-sp2011/slopefields_howto.pdf *)

plot1 = VectorPlot[  {1, Cos[x y ] }, {x,-Pi, Pi}, {y, -2, 2}, Frame->False, Axes ->True, AxesStyle->Black]

(* Solve the equation numerically *)

Clear[y];
solutionN = NDSolve[ {y'[x] == Cos[x y[x]], y[0]==0},y,{x,-3,3}]

(* Plot the numerical solution *)

plot2=Plot[Evaluate[y[x]/.solutionN],{x,-3,3}, PlotStyle->Red]

(* Plot the numerical direction field and the numerical solution together *)

Show[plot1,plot2]

(* EULER'S METHOD.
Instructions from: http://mathematica.stackexchange.com/questions/11924/ndsolve-with-euler-method *)

Clear[y];
y=y/.First[NDSolve[{y'[x]== Cos[ x * y[x] ],y[0]==0},y,{x,0,1.5},StartingStepSize->0.5,Method->{"FixedStep",Method->"ExplicitEuler"}]];
grid=Table[{x,y[x]},{x,0,1.5,0.5}];

(* Show the list of points (x_m, y_m).
Plot the numerical solution: plot the points (x_m,y_m) and joint them by straight lines *)

grid
ListLinePlot[grid]

(* Redo with smaller step size *)

Clear[y];
y=y/.First[NDSolve[{y'[x]== Cos[ x * y[x] ],y[0]==0},y,{x,0,1.5},StartingStepSize->0.1,Method->{"FixedStep",Method->"ExplicitEuler"}]];
grid=Table[{x,y[x]},{x,0,1.5,0.1}];

grid
ListLinePlot[grid]

(* Show the list of points (x_m, y_m).
Plot the numerical solution: plot the points (x_m,y_m) and joint them by straight lines *)

Clear[y];
y=y/.First[NDSolve[{y'[x]== Cos[ x * y[x] ],y[0]==0},y,{x,0,1.5},StartingStepSize->0.01,Method->{"FixedStep",Method->"ExplicitEuler"}]];
grid=Table[{x,y[x]},{x,0,1.5,0.01}];

grid
ListLinePlot[grid]