i am unable to get a plot .....as my plot fig is coming empty .....please tel me if anyone notice my mistake i paste the code
2 visualizaciones (últimos 30 días)
Mostrar comentarios más antiguos
x(1)=-4; %initial position of particle in x direction
y(1)=1; %initial position of particle in y direction
nsteps=100;
x(i-1)=zeros(1,nsteps);
y(i-1)=zeros(1,nsteps);
R(i-1)=zeros(1,nsteps);
R(1)=sqrt(x(1)^2+y(1)^2);
p=1780;
f=1.2;
dp=100e-6;
m = 15.11e-4; %viscosity of fluid
Ux(1)=-0.986; %initial velocity in horizontal direction
Uy(1)=-0.00503; %initial velocity in vertical direction.
Ux(i-1)=(-1)*(1-1./(R(i-1)).^3+3.*y(1).^2./(2.*(R(i-1)).^5));
Uy(i-1)=3.*x(i-1).*y(i-1)./(2.*R(i-1).^5);
Ux=zeros(1,nsteps+1);
Uy=zeros(1,nsteps+1);
c=p*dp/f;
i=2:1:101;
RX(i-1)=Ux(i-1)*dp/m; %reynolds number in x direction
RY(i-1)=Uy(i-1)*dp/m; %reynolds number in y direction
ts=p*dp^2/(18*m);
cdu(1)=24/RX(1);
cdv(1)=24/RY(1);
if(RX(i-1) >= 3.0)
cdu(i-1)=(24/RX(i-1))+(4/RX(i-1)^(1/3));
elseif(RX(i-1) >= 0.5)
cdu(i-1)=(24/RX(i-1))+(3.6/RX(i-1)^(0.313));
elseif(RX(i-1) >= 0.1)
cdu(i-1)=(24/RX(i-1))+4.5;
elseif(RX(i-1) >=0.0001)
cdu(i-1)=24/RX(i-1);
end
while i>2;
h1=(1/(cdu(i-1)*RX(i-1)))*(c)*(-1.33);
h2=(1/(cdu(i-1)*RX(i-1)))*c*(-1.33);
h3=(1/(cdu(i-1)*RX(i-1)))*c*(-1.33);
h4=(1/(cdu(i-1)*RX(i-1)))*c*(-1.33);
j=ts/10;
x(i-1)=j/3*(h1+(4*h2)+(2*h3)+h4);
end
if(RY(i-1) >= 3.0)
cdu(i-1)=(24/RY(i-1))+(4/RY(i-1)^(1/3));
elseif(RY(i-1) >= 0.5)
cdu(i-1)=(24/RY(i-1))+(3.6/RY(i-1)^(0.313));
elseif(RY(i-1) >= 0.1)
cdu(i-1)=(24/RY(i-1))+4.5;
elseif(RY(i-1) >=0.0001);
g=9.81; %acceleration due to gravity
while i>2;
k1=sqrt((dp*f*g)/(3*cdv(i-1)*p));
k2=sqrt((dp*f*g)/(3*cdv(i-1)*p));
k3=sqrt((dp*f*g)/(3*cdv(i-1)*p));
k4=sqrt((dp*f*g)/(3*cdv(i-1)*p));
j=ts/10;
y(i-1)=j/3*(k1+(4*k2)+(2*k3)+k4);
end
end
plot(x,y)
0 comentarios
Respuestas (1)
Walter Roberson
el 24 de Abr. de 2013
Editada: Walter Roberson
el 24 de Abr. de 2013
You have
i=2:1:101;
but in some places (such as the "if" statements) you treat "i" as if it was a scalar instead of a vector.
When you have a construct such as
if R(i) > 5
... do something ...
end
when "i" is a vector, it does not mean that the code should be done only for the values of "i" for which R(i) > 5. In MATLAB, it would instead being interpreted as doing the code only if all R(i) > 5 .
Use a "for" loop, or study logical indexing.
2 comentarios
Walter Roberson
el 25 de Abr. de 2013
You initialize x(1) and y(1) first, and then assign zeros to x and y, so x and y will be all 0. Then you calculate R(1) based on x(1) and y(1), so R(1) will be 0. You divide by R (all 0) so you get NaN. Everything goes downhill from there.
Ver también
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!