Euler 's Method for falling velocity ode..

2 views (last 30 days)
Mass = 100/6*pi*0.04^3;
g=9.8;
A=0.25*pi*0.04^2;
U=0.0001;
dt = 0.01;
for t=0:dt:15-dt
Re = 1.25*0.04*U/0.000015;
CD=24./Re+(2.6*Re/5)./(1+(Re/5).^1.52)+(0.411*(Re/263000).^(-7.94))./(1+(Re/263000).^(-8))+(0.25*Re/10^6)./(1+Re/10^6);
dUdt = g - (0.5*1.25*U(end)^2*A*CD)/Mass
U(end+1)=U(end)+dUdt*dt
end
plot(0:dt:15, U)
error occurs in "U(end+1)=U(end)+dUdt*dt" sentence
I can't find any problem... Could you find problem in this code? I think 'Re' and 'CD' also looks doubtful..
Thanks for reading my quenstion

Accepted Answer

Sugar Daddy
Sugar Daddy on 1 Jun 2020
Mass = 100/6*pi*0.04^3;
g=9.8;
A=0.25*pi*0.04^2;
U=0.0001;
dt = 0.01;
for t=0:dt:15-dt
Error is in this line as the U size increase Re size increases so put U(end) here too
% Re = 1.25*0.04*U/0.000015; % Error
Re = 1.25*0.04*U(end)/0.000015; % Correct
CD=24./Re+(2.6*Re/5)./(1+(Re/5).^1.52)+(0.411*(Re/263000).^(-7.94))./(1+(Re/263000).^(-8))+(0.25*Re/10^6)./(1+Re/10^6);
dUdt = g - (0.5*1.25*U(end)^2*A*CD)/Mass
U(end+1)=U(end)+dUdt*dt
end
plot(0:dt:15, U)

More Answers (0)

Categories

Find more on Programming in Help Center and File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!