Asked by Neel Singh
on 15 Jun 2015

clear all clc rho = 2000; Cd = 0.44; %d = input('input diameter'); d = 0.01; %d = 9.99*10^-4; A = (d/2)^2*pi; %m = input('Enter multiple of mass: ')*2.8*10^-10; %m = 2*10^-10; m = rho * (4/3)*pi* (d/2)^3; dt = 0.01; N = 100; t = 0:dt:N*dt; T = dt * N; Sy = 0 ;

Vx(1) = 0; Ux = 0; Uy = 1; Sx(1) = 0;

for i = 1:N

Vy(1) = 0;

dSy(1) = 0;

if Sy <= 0.5

%For y direction

Vy(i+1) = (Cd*A*rho*(Uy-Vy(i))^2*dt)/(2*m) + Vy(i);

Ay(i) = (Vy(i+1)-Vy(i))/dt;

dSy(i) = ((Vy(i+1))^2 - (Vy(i))^2 )/ (2*Ay(i));

%disp(sum(Sy))

Sy = cumsum(dSy);

%For x direction

Vx(i+1) = Cd*A*rho/(2*m)*(Ux-Vx(i))^2*dt+Vx(i);

Ax(i) = (Vx(i+1)-Vx(i))/dt;

dSx(i) = (Vx(i)*dt);

Sx = cumsum(dSx);

else

Uy = 0;

Ux = 1;

Deacc = Cd*A*rho/(2*m)*(Uy-Vy(i))^2*dt

%For y direction

Vy(i+1) = Vy(i) - Deacc ;

%Vy(i+1) = Vy(i)- Vy(i-1);

Ay(i) = (Vy(i+1)-Vy(i))/dt;

dSy(i) = (Vy(i)*dt);

%disp(sum(Sy))

Sy = cumsum(dSy);

%For x direction

Vx(i+1) = Cd*A*rho/(2*m)*(Ux-Vx(i))^2*dt+Vx(i);

Ax(i) = (Vx(i+1)-Vx(i))/dt;

dSx(i) = ((Vx(i+1))^2 - (Vx(i))^2 )/ (2*Ax(i));

Sx = cumsum(dSx);

end

end

x=Sx;

y=Sy;

figure (1)

plot(x,y);

hold on;

h=plot(x(1),y(1),'r*');

This "d" value gives answers well, this is how its supposed to work. however changing it to 0.001 will give infinity to Vy, Vy, Sy, Sx and i cant figure out why. The only thing i found online is that it might be a round off from matlab at some point, driving everything divided by 0 to infinity. How can I work around this, so I can use small values of d ?

Thanks for any help,

Im new to this forum, so excuse my mistakes.

Answer by Jonathan Poirier
on 15 Jun 2015

First, never use i in a loop. Matlab use i for a complex number (Try to type i in command window >> i ans =

0.0000 + 1.0000i

Second, I never try to understand yous matematical equation, but -1.32e239 look inf, no?

If you really want to figure there big values, divide (by exemple) 1 e100 your value and in your y axis :

xlabel('data 1e100');

Neel Singh
on 15 Jun 2015

Hi thanks for the suggestions, however those are not the cause of the problem.

Sign in to comment.

Answer by Katalin
on 15 Jun 2015

I think it's not a rounding issue. In Vy(i+1) = (Cd*A*rho*(Uy-Vy(i))^2*dt)/(2*m) + Vy(i); The fraction (Cd*A*rho*dt)/(2*m) is either larger or smaller than 1 depending on the variable “d”.

It is equal to 1, when d = Cd*dt/(4/3) = 0.0033. When you are trying numbers below d = 0.0033 the system has a complete different behavior and it goes to extreme very fast.

Neel Singh
on 15 Jun 2015

Sign in to comment.

Opportunities for recent engineering grads.

Apply Today
## 0 Comments

Sign in to comment.