%Exploring stiffness in ODE

%Example: ODE for the evolution of the ball of flame radius on the match 

%delta is initial radius (critical parameter for the solution)

%delta not very small - problem is not very stiff
delta=0.01;
F=inline('y^2 - y^3','t','y');
opts=odeset('RelTol',1.e-4);
ode45(F,[0 2/delta],delta,opts);

%delta decreased by two orders of magnitude - stiffness!
delta=0.0001;
F=inline('y^2 - y^3','t','y');
opts=odeset('RelTol',1.e-4);
ode45(F,[0 2/delta],delta,opts);


delta=0.0001;
F=inline('y^2 - y^3','t','y');
opts=odeset('RelTol',1.e-4);
ode23s(F,[0 2/delta],delta,opts);

y=dsolve('Dy=y^2-y^3','y(0)=1/100');
y=simplify(y);
pretty(y)
ezplot(y,0,200);