# Help ploting different particles in a box

22 views (last 30 days)
MFDOOM on 27 Mar 2020
Commented: MFDOOM on 27 Mar 2020
I want to create a code that simulates multiple particles coliding with the walls of a box but I dont want the particles to colide with each other.
I alredy made a script for one particle, but Im having trouble making one for multiple ones.
So far I got this:
Here's the main script:
a=20;
b=10;
R=0.5;
dt=0.1;
np=3; %Number of particles
rx=a*rand(1,np);
ry=b*rand(1,np);
r=[rx;ry];
v=5*(2*rand(2,np)-1);
tmax=100;
t=0;
while (t<tmax)
for i=1:np
deltat(i,1)=(R-r(1,i))/v(1,i); deltat(i,2)=(a-R-r(1,i))/v(1,i);
deltat(i,3)=(R-r(2,i))/v(2,i); deltat(i,4)=(b-R-r(2,i))/v(2,i);
if deltat(i,1)<=4*eps
deltat(i,1)=1e6;
end
if deltat(i,2)<=4*eps
deltat(i,2)=1e6;
end
if deltat(i,3)<=4*eps
deltat(i,3)=1e6;
end
if deltat(i,4)<=4*eps
deltat(i,4)=1e6;
end
end
deltat_wall=[deltat(:,1),deltat(:,2),deltat(:3),deltat(:,4)]; %must be a matrix of np particles and 4 columns
[deltat, column]=min(deltat_wall,[],2);
[deltat,particle]=min(deltat);
wall=column(particle);
animacion(r,v,a,b,R,dt,deltat) % the function is just to plot the box and the posicion of the particles
r=r+v*deltat;
ep=[[1;0],[-1;0],[0;1],[0;-1]];
v=v-2*dot(v,ep(:,wall))*ep(:,wall);
t=t+deltat;
end
Here´s the function:
function []= animacion(r,v,a,b,R,dt,deltat)
x1=0;
x2=a;
y1=0;
y2=b;
side_a=[x1,x2,x2,x1,x1];
side_b=[y1,y1,y2,y2,y1];
for t=0:dt:deltat
plot(side_a,side_b,'b-','Linewidth',2)
hold on
plot(r(1,:),r(2,:),'ko','MarkerFaceColor','k', 'MarkerSize',28*R)
xlim([-10,30])
ylim([-10,20])
drawnow
hold off
r=r+v*dt;
end
end

Geoff Hayes on 27 Mar 2020
Can you comment on what this block of code is doing:
[deltat, column]=min(deltat_wall,[],2);
[deltat,particle]=min(deltat);
wall=column(particle);
animacion(r,v,a,b,R,dt,deltat) % the function is just to plot the box and the posicion of the particles
r=r+v*deltat;
ep=[[1;0],[-1;0],[0;1],[0;-1]];
v=v-2*dot(v,ep(:,wall))*ep(:,wall);
When I run your code with 2 or more particles, I see the following error
Error using dot (line 33)
A and B must be same size.
Can you explain what the above lines are doing? What are you expecting to happen when you have more than one particle?
MFDOOM on 27 Mar 2020
The block of code is to determinate what wall the 1st particle will hit and then when the 2sd will hit and the 3rd and the 4th ...
given the position and velocity that are column vectors with the first row for x and the second row for y. Can be used r only without code.
t = 0 while (t <tmax)
Within the cycle while he has to do the following:
a) calculate deltat for shock. Here, construct array 2 by 4 points_details, set the minimum values ​​equal to a large value and determine the smallest value, or the index of the smallest wall element indicates a wall. A wall can be identified by a numerical index. For example, wall 1 is on the left, wall 2 is on the right, wall 3 is the bottom and wall 4 is the top.
[deltat, column] = min (deltat_wall, [], 2); % determines the minimum between the columns of the matrix for each row [deltat, particula] = min (deltat); % determines the particle with shorter shock time
wall = column (particle); % determines the wall where this shock occurs
b) makes an animation animation (r, v, a, b, R, dt, deltat); this function only plots the box and position of the particle that moves with constant speed v.
c) updates the particle position with the old speed r = r + v * deltat
d) calculate a new speed after the shock. It is explained how it does not spelled out and in the list of frequently asked questions above, but repeat again. To collide with a wall, you can arrange the ep vector in a matrix with 2 rows and 4 columns for each column corresponding to a wall. For example, ep = [[1; 0], [- 1; 0], [0; 1], [0; -1]) where the numbering of the walls is given above. Each column of ep is the unit vector perpendicular to the wall with the column index of the ep array. Thus, it is possible to determine a collision with a wall with a wall index to obtain a velocity after the collision causes v = v-2 * point (v, ep (:, wall)) * ep (: wall). v on the right side of the speed signal is the speed before the shock. e) update time t = t + deltat . Time increases a shock to another in deltat .

Geoff Hayes on 27 Mar 2020
If the last step is to just calculate a new speed after the shock then this would apply to a single particle only. So we would modify the velocity of just that one particle instead of trying to modify the velocities for each. Try changing
v=v-2*dot(v,ep(:,wall))*ep(:,wall);
to
v(:, particle) = v(:, particle) - 2*dot(v(:, particle),ep(:,wall))*ep(:,wall);
where the only difference is how we index into v.

#### 1 Comment

MFDOOM on 27 Mar 2020
OMG thank you so much dude. You are amazing , you dont know how much this helped me :)