3.1(b) Terminal Particle Velocity - Successive Substitution
g=9.80665
rhop=7860
rho=1.23
Dp=0.5e-3
vis=1.79e-5
vt0=2
Re0=vt0*Dp*rho/vis
CD0=if(Re0<0.1)then(24/Re0)else(24*(1+0.14*Re0^0.7)/Re0)
vt1=sqrt(4*g*(rhop-rho)*Dp/(3*CD0*rho))
err0=abs(vt1-vt0)
Re1=vt1*Dp*rho/vis
CD1=if(Re1<0.1)then(24/Re1)else(24*(1+0.14*Re1^0.7)/Re1)
vt2=sqrt(4*g*(rhop-rho)*Dp/(3*CD1*rho))
err1=abs(vt2-vt1)