Numerical ODE Solving: BUILT-IN

Clear[n,A, f,y,TMax]
n=12;
A=RandomReal[{-1,1},{n,n}];
f[y_]:= A.y
y0=RandomReal[{-1,1},n];
TMax=3;

sol=NDSolve[{
y'[t]==f[y[t]],
y[0]==y0
},
y[t],
{t, 0, TMax}]

ySol[t_]= y[t]/.sol[[1]]
Plot[ySol[t],{t, 0, TMax}]

Numerical ODE Solving: Explicit Euler Fixed Step

(Subscript[y, n+1])=(Subscript[y, n])+h f[(Subscript[y, n])]

Clear[n,A, f,y,TMax]
n=12;
A=RandomReal[{-1,1},{n,n}];
f[y_]:= A.y
y0=RandomReal[{-1,1},n];
NMax=12;
h=0.1;
{t,y}={0, y0};
Table[
{t,y}={t,y}+h {1, f[y]},
{NMax}]

Numerical ODE Solving: Explicit RK Fixed Step

Clear[n,A, f,y,TMax]
n=12;
A=RandomReal[{-1,1},{n,n}];
f[y_]:= A.y
y0=RandomReal[{-1,1},n];
NMax=12;
h=0.1;
{t,y}={0, y0};
Table[
k1=f[y];
k2=f[y+((h )/2)k1];
k3=f[y+(h/2)k2];
k4=f[y+h k3];
{t,y}={t,y}+h {1, (k1+2 k2 + 2 k3 + k4)/6},
{NMax}]

Numerical ODE Solving: The problems with explicit solvers.

Stiffness

Numerical ODE Solving: Implicit/Backwards Euler

(Subscript[y, n+1])=(Subscript[y, n])+h f[(Subscript[y, n+1])]

Big problem.  We need to solve an equation for (Subscript[y, n+1])!  This almost always looks like a Newton type iteration with an initial guess.

What does one step look like in code.

Clear[n,A, f,y,TMax]
n=12;
A=RandomReal[{-1,1},{n,n}];
f[y_]:= A.y + 0.01 Sin[y]
y0=RandomReal[{-1,1},n]

h=0.9;
{t,y}={0, y0};
G[z_List]:= y+h f[z]-z
FindRoot[
G[z]==0,
{z, y}]

NMax=12;
Table[
{t,y}={t,y}+h {1, f[y]},
{NMax}]