Radioactive Decay -- Mid-Point solution
Mid-Point vs exact
restart;
# rhs function
rhsN := unapply( -lambda*N, N );
lambda := 2.25;
# time grid
dt := 0.1:
Nt := 30:
# initial value
NN[1] := 10000:
tt[1] := 0.0:
# loop over time
for it from 2 to Nt
do
NN_half := NN[it-1] + rhsN( NN[it-1] ) * dt / 2.0;
NN[it] := NN[it-1] + rhsN( NN_half ) * dt;
tt[it] := tt[it-1] + dt;
end do:
exact := unapply( NN[1] * exp( -lambda * t ), t );
plot( [ [ [ tt[ii], NN[ii] ] $ii=1..Nt ], exact(t) ],
t=tt[1]..tt[Nt], style=[ point, line ] );
`Check`;
Nt, tt[Nt], NN[Nt], exact( tt[Nt] );