Plotting Several Data Sets in a For Loop
49 次查看(过去 30 天)
显示 更早的评论
I am working through a problem that requires me to calculate trajectories of an object based on the angle the object is initially thrown. The end goal of this portion of the problem is to have one graph that shows all of the trajectories on one plot. The angles I am trying to iterate through are 0 to 80 in steps of 10 degrees. However, when I run the script, I am getting one plot where all of the data is the same and all of the data (regardless of the angle) has the same output (as evidenced by the legend on the plot). If I change the angle variable to just one value (i.e. to 10, instead of 0:10:80), the program will successfully output a graph for that angle. Can someone help me sort out what my issue is?
clc
clear
% Constant Variables
m=0.175; %mass, kg
D=0.26; %Diameter,
A=(pi/4)*D^2; %Area, m^2
g=9.81; %m/s^2
rho_air=1.23; %Density of air, kg/m^3
AoA_rad=deg2rad(-4); %Angle of Attack for Drag Calc.
%Position Variables
x_prev=0;
y_prev=1;
x_new=0;
y_new=1;
%Time Variables
t_start=0;
t_step=.01;
t_end=7;
t_total=0;
time=t_start:t_step:t_end;
%Loop Variables/Matrices
Vel=14;
Coord=zeros(length(time),2);
count=1;
%Plot information
figure
hold on
title('Frisbee with Constant Velocity and Varying Angles')
xlabel ('Horizontal Distance (m)')
ylabel('Vertical Distance(m)')
grid on
legend
for j=0:10:80 %Angles to iterate through
for i=t_start:t_step:t_end
alpha_rad=deg2rad(j);
beta=deg2rad(90-j);
%Defining and Summing Forces
Cd=0.08+(2.72*alpha_rad)*(alpha_rad-AoA_rad)^2;
Cl=0.15+(1.4*alpha_rad);
Fd=(-Cd*rho_air*A*(Vel^2))/2;
Fl=(Cl*rho_air*A*(Vel^2))/2;
Sum_Fx=Fd*sin(beta)-Fl*sin(alpha_rad);
Sum_Fy=Fd*cos(beta)+Fl*cos(alpha_rad)-(m*g);
%Finding Accelerations and New Velocity
ax=Sum_Fx/m;
ay=Sum_Fy/m;
u0x=Vel*cos(alpha_rad);
v0y=Vel*sin(alpha_rad);
u0=u0x+i*ax;
v0=v0y+i*ay;
%Finding x and y coordinates
x_new=x_prev+u0x*i+0.5*ax*(i^2);
y_new=y_prev+v0y*i+0.5*-g*(i^2);
Coord(count,1)=x_new;
Coord(count,2)=y_new;
x_prev=x_new;
y_prev=y_new;
count=count+1;
Vel=sqrt((u0^2)+(v0^2));
end
plot(Coord(:,1),Coord(:,2))
ylim([0 10])
xlim([0 30])
end
0 个评论
回答(1 个)
William Rose
2023-3-9
You write:
x_new=x_prev+u0x*i+0.5*ax*(i^2);
which is not correct, given how you define variables. You define i as current time. xprev is updated on each inner loop pass. Therefore you should be using dt, the step size, rather than i=time, for this update step. Your equation for y_new has the same issue.
Inspect your inner loop code to identify calculations that can be pulled to the outer loop - that is, anything that is the same on each inner loop pass.
I would use NxM arrays for X, Y, Vel, etc., where N=length(time) and M=number of angles.
I would change the loop variables to integers and use them as indices for the arrays x(i,j) and y(i,j), which you will fill up as you go through the loops.
Here is a simple example for a ball launched at different angles. Your system is more complicated, so your equations for ax, ay will be different.
dt=.01; % time step (s)
time=0:dt:2; % time is not used in the loops, but you could use it to plot(t,x(:,j)), etc
angle=0:10:80;
N=length(time); M=length(angle);
g=-9.8; % y accel m/s^2
x=zeros(N,M); y=zeros(N,M); % allocate arrays for x, y
x(1,:)=zeros(1,M); % initial distance (m)
y(1,:)=ones(1,M); % initial height (m)
v0=14; % initial veloc (m/s)
for j=1:M
vx=v0*cos(angle(j)*pi/180);
vy=v0*sin(angle(j)*pi/180);
ax=0; % no x acceleration in my simple model
ay=g; % y acceleration
for i=2:N
x(i,j)=x(i-1,j)+dt*vx;
y(i,j)=max(y(i-1,j)+dt*vy,0); % don't allow y<0
vx=vx+dt*ax; % update vx
vy=vy+dt*ay; % update vy
end
plot(x(:,j),y(:,j))
hold on
legstr{j}=num2str(angle(j));
end
legend(legstr); grid on;
xlabel('Distance (m)'); ylabel('Height (m)');
That gives you some ideas.
7 个评论
William Rose
2023-3-9
@George, @Voss, my code is below. It is my attempt to implement the equaitons as I wrote them in an earlier comment. I think these are the equaitons that you are trying to solve. And I did it with a first order Euler method, which is more or less what you were doing - although your method is partly second order, because you incoporate acceleration in the update for x and y.
I don;t know why @Voss's modification of your code gives different answers than mine. If you really want to know, you could allocate arrays for Fx(i,j) and Fy(i,j) in each script and use them to save the forces. Then you can plot the forces in the two scripts for the different angles. This would likely provide insight.
William Rose
2023-3-9
Your script incldues the lines
u0x=Vel*cos(alpha_rad);
v0y=Vel*sin(alpha_rad);
After that, you use u0x and v0y to update x and y. In other words, the equations above are the x, y components of instantaneous velocity. I think the equations are incorrect, because they use alpha_rad, which is the initial launch angle. Alpha_rad does not get updated as the trajectory changes. It should. Or you can do what I did, which is to update vx and vy, based on the accelerations ax and ay.
vx=vx+dt*Fx/m; % update vx
vy=vy+dt*Fy/m; % update vy
Use vx and vy to update the positions:
x(i,j)=x(i-1,j)+dt*vx; % update x
y(i,j)=y(i-1,j)+dt*vy; % update y
That may explain the difference.
另请参阅
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!