How to plot by using result data?

8 ビュー (過去 30 日間)
Thin Rupar Win
Thin Rupar Win 2021 年 11 月 10 日
コメント済み: Thin Rupar Win 2021 年 11 月 10 日
Dear Sir or Madam,
I am new to matlab programming I would like to know how to plot position_x , position_y , velocity x , velocity y by applying result data. Please show me how to perform the plot as I would like to see them.
Hopefully and Thank you very much for your help.
n_part=10;
kn=5;
kt=2/7*kn;
m=0.3;
g=9.81;
timestep=50;
v_init=0.5;
% initial of position, velocity, angular velocity
position_x(:,1:n_part)=rand(size(1:n_part));
position_y(:,1:n_part)=rand(size(1:n_part));
theta=2*pi*rand(size(position_x));
velocity_x=v_init*sin(theta);
velocity_y=v_init*cos(theta);
w_x=zeros(size(position_x));
w_y=zeros(size(position_y));
w_half_x=zeros(size(velocity_x));
w_half_y=zeros(size(velocity_y));
w_x_dot=zeros(size(velocity_x));
w_y_dot=zeros(size(velocity_y));
v_half=rand(size(velocity_x));
v_half_x=rand(size(velocity_x));
v_half_y=rand(size(velocity_y));
dt=0.001;
rad(1:n_part)=0.5;
acc_x=zeros(size(position_x));
acc_y=zeros(size(position_y));
acc_i=[0,0];
acc_j=[0,0];
w_i=[0,0];
w_j=[0,0];
Fn=[0,0];Fn_x=[0,0];Fn_y=[0,0];
for t=1:1:4
for i=1:n_part
v_half_x(t+1,i)=velocity_x(t,i)+0.5*dt*acc_x(t,i);
v_half_y(t+1,i)=velocity_y(t,i)+0.5*dt*acc_y(t,i);
position_x(t+1,i)=position_x(t,i)+v_half_x(t+1,i)*dt;
position_y(t+1,i)=position_y(t,i)+v_half_y(t+1,i)*dt;
% update of w_ part
w_half_x(t+1,i)=w_x(t,i)+dt*0.5*dt*w_x_dot(t,i);
w_half_y(t+1,i)=w_y(t,i)+dt*0.5*dt*w_y_dot(t,i);
j=1;
for j=i+1
if j<=n_part
lx=position_x(i)-position_x(j);
ly=position_y(i)-position_y(j);
vx=velocity_x(i)-velocity_x(j);
vy=velocity_y(i)-velocity_y(j);
root_xy=sqrt(lx.^2+ly.^2);
R_i=rad(i);
R_j=rad(j);
delta=(R_i+R_j)-root_xy;
if delta>0
N=[lx;ly];
n_vector=N/norm(N);
% relative vector to find tangent vector;
vij1=[vx;vy];
w_i=[w_x(i);w_y(i)];
w_j=[w_x(j);w_y(j)];
vij2=vij1+(R_i.*w_i+R_j.*w_j).*n_vector;
v_i_j=-n_vector.*(n_vector.*vij2);
% tangent vector
tan=v_i_j/norm(v_i_j);
Fn=kn.*delta;
Ft=kt.*v_i_j;
F_total=Fn.*n_vector+Ft.*tan+m*g;
Fn_x=Fn_x-F_total(1,:);
Fn_y=Fn_y+F_total(2,:);
%Momemtum of particles
M_i=R_i.*n_vector.*Ft;
M_j=R_j.*n_vector.*Ft;
% same radius for movement of inertia of circle
I=pi*R_i/4;
w_i=0;
w_j=0;
w_i=w_i+M_i'/I;
w_j=w_i+M_j'/I;
w_x_dot(t,i)=w_x_dot(t,i)+w_i(1);
w_x_dot(t,i)=w_x_dot(t,i)+w_j(1);
w_y_dot(t,i)=w_y_dot(t,i)+w_i(2);
w_y_dot(t,i)=w_y_dot(t,i)+w_j(2);
w_x_dot(t+1,:)=w_x_dot(t,:);
w_y_dot(t+1,:)=w_y_dot(t,:);
acc_i=0;
acc_i=acc_i+Fn_x./m;
acc_j=0;
acc_j=acc_j+Fn_y./m;
acc_x(t,i)=acc_x(t,i)+acc_i(1);
acc_y(t,i)=acc_y(t,i)+acc_i(2);
acc_x(t,j)=acc_x(t,j)+acc_j(1);
acc_y(t,j)=acc_y(t,j)+acc_j(2);
acc_x(t+1,:)=acc_x(t,:);
acc_y(t+1,:)=acc_y(t,:);
end
end
end
velocity_x(t+1,i)=v_half_x(t+1,i)+0.5*dt*acc_x(t+1,i);
velocity_y(t+1,i)=v_half_y(t+1,i)+0.5*dt*acc_y(t+1,i);
w_x(t+1,i)=w_half_x(t+1,i)+0.5*dt*w_x_dot(t+1,i);
w_y(t+1,i)=w_half_y(t+1,i)+0.5*dt*w_y_dot(t+1,i);
end
end
% my intension is using rad, position_x, position_y to create circle as
% particles and using velocity_ x and velocity _y to move them.
% position_x (t,i); position_y (t,i); velocity_x (t,i); velocity_y (t,i);
% rad=0.5;

回答 (1 件)

KSSV
KSSV 2021 年 11 月 10 日
How about using quiver?
quiver(position_x,position_y,velocity_x,velocity_y)
  1 件のコメント
Thin Rupar Win
Thin Rupar Win 2021 年 11 月 10 日
Sir, This is not working.

サインインしてコメントする。

カテゴリ

Help Center および File ExchangeVector Fields についてさらに検索

タグ

Community Treasure Hunt

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

Start Hunting!

Translated by