Simulation of Elastic Collision

[MatinSAR]

Homework Statement

Create a simulation to model elastic collisions among particles in a 2D square area. Ensure the particles reflect elastically off the boundaries.

Relevant Equations

$\vec P_i = \vec P_f$
$T_i = T_f$

Hello everyone,

Can we use the conversion of kinetic energy (T) and momentum ($\vec P$) to find the velocities after an elastic collision in a 2D problem? Can we treat objects like particles?

 

[kuruman]

Sure. Start with the equations giving the final velocities of the two particles in terms of the initial velocities in one dimension. (See http://hyperphysics.phy-astr.gsu.edu/hbase/elacol2.html). Note that this "one" dimension is along a line joining the centers of the masses and the ignored, "second", dimension is the direction perpendicular to that. In a two dimensional simulation and in the absence of friction, the component of the velocity of each mass in the perpendicular direction does not change. What changes is the parallel component of the velocity of each mass according to the equations provided in the reference.

Can you find the parallel component for an arbitrary direction of the velocity of each mass? Hint: Find a unit vector in the direction from mass 1 to mass 2.

 

[MatinSAR]

Can you find the parallel component for an arbitrary direction of the velocity of each mass? Hint: Find a unit vector in the direction from mass 1 to mass 2.

Given that we know the positions of each particle, we can use $\dfrac{\vec{r}_2 - \vec{r}_1}{|\vec{r}_2 - \vec{r}_1|}$ as the unit vector pointing from particle 1 to particle 2.

I believe I understand your concept. However, there's a problem. In my code, I can't predict which two particles will collide to use the line connecting them as our one-dimensional reference, where velocity changes only in that specific dimension.
Additionally, I don't need an analytical solution. As long as I have the correct equations that provide the magnitude and direction of each particle's speed after a collision, I can solve the system of equations numerically.

 

[MatinSAR]

This was my intended approach.

$$m_1u_1=m_1v_1\cos\theta_1 + m_2v_2\cos\theta_2$$$$0=m_1v_1\sin\theta_1 - m_2v_2\sin\theta_2$$$$m_1u^2_1=m_1v^2_1+m_2v^2_2$$
There are two issues here. First, we have four unknowns and only three equations. Second, is this really a 2D problem? I mean, in 2D, can we consider particles to only have head-on collisions even if they have orthogonal velocities?

 

[erobz]

This was my intended approach.

View attachment 355282​

$$m_1u_1=m_1v_1\cos\theta_1 + m_2v_2\cos\theta_2$$$$0=m_1v_1\sin\theta_1 - m_2v_2\sin\theta_2$$$$m_1u^2_1=m_1v^2_1+m_2v^2_2$$
There are two issues here. First, we have four unknowns and only three equations. Second, is this really a 2D problem? I mean, in 2D, can we consider particles to only have head-on collisions even if they have orthogonal velocities?

How are they flying off in 2 directions if they are particles without dimension. I could see this scenario as something that could happen with things colliding that have some physical geometry (where a glancing blow is possible), but point masses...I don't get it.

To me, there are three scenarios $m_1$ stops - $m_2$ moves directly to the right, $m_1$ moves to the left - $m_2$ to the right, or they both continue to the right.

 

[MatinSAR]

How are they flying off in 2 directions if they are particles without dimension. I could see this scenario as something that could happen with things colliding that have some physical geometry (where a glancing blow is possible), but point masses...I don't get it.

To me, there are three scenarios $m_1$ stops, $m_2$ moves directly to the right, $m_1$ moves to the left, $m_2$ to the left, and they both continue to the right.

Third problem. Perhaps I should visit my professor's office tomorrow to discuss these issues with him.

 

[Steve4Physics]

Homework Statement: Create a simulation to model elastic collisions among particles in a 2D square area. Ensure the particles reflect elastically off the boundaries.

The problem statement seems very vague Some questions that spring to mind are:
- are the particles to be modelled as points or as discs?
- will this be a real-time simulation?
- how many (roughly) particles are in the square?
- what are the rough relative sizes of the particles, box and average particle-separation?
Or are you allowed to decide these sorts of parameters for yourself?

Knowing the key parameters could help you to decide what simplifications/approximations (if any) you might want to make. This might have an impact on your approach.

The level of sophistication of the simulation will be limited by how many person-hours of effort are allocated for the work – and the level of your own programming skills! Is there a ball-park figure?

A good start point (IMO) would be to establish answers to the above (and any other relevant) questions before diving into the detail. Then formulate an overall strategy for how to tackle the job.

So maybe there are some additional things to discuss with your professor!

 

[erobz]

Third problem. Perhaps I should visit my professor's office tomorrow to discuss these issues with him.

I could be missing something; someone will straighten this out. Spend some time thinking about it today, and then go to your prof. I was just throwing in my two cents about specifically that scenario you proposed. An oblique collision could send particles (without physical geometry) off in different directions.

 

[kuruman]

How are they flying off in 2 directions if they are particles without dimension. I could see this scenario as something that could happen with things colliding that have some physical geometry (where a glancing blow is possible), but point masses...I don't get it.

To me, there are three scenarios $m_1$ stops - $m_2$ moves directly to the right, $m_1$ moves to the left - $m_2$ to the right, or they both continue to the right.

There is no mention in the original post about point masses. As you say, there can be no 2d collision between point masses. I would assume discs of radii $R_1$ and $R_2$. If $\mathbf r_1$ and $\mathbf r_2$ are their positions, the criterion for a collision to occur is $|\mathbf r_2-\mathbf r_1|\leq(R_1+R_2).$

To @MatinSAR:
I am in the process of answering your posts #3 and #4. Please be patient.

 

[MatinSAR]

Or are you allowed to decide these sorts of parameters for yourself?

Initially, all the decisions are mine. Once I complete my work based on my preferences, I will discuss it with him, and he will add more details. I recall he mentioned that he might ask me to apply it to random shapes like squares or pentagons, and I would need to consider rotation and rotational inertia as well. I really hope that won't be the case since it's complicated enough already.

 

[Steve4Physics]

... I recall he mentioned that he might ask me to apply it to random shapes like squares or pentagons, and I would need to consider rotation and rotational inertia as well.

That seems to imply that you are expected (at least initially) to model the 'particles' as discs rather than as points.

Squares or pentagons? - aagh!

really hope that won't be the case since it's complicated enough already.

Yes indeed!

 

[kuruman]

Initially, all the decisions are mine. Once I complete my work based on my preferences, I will discuss it with him, and he will add more details. I recall he mentioned that he might ask me to apply it to random shapes like squares or pentagons, and I would need to consider rotation and rotational inertia as well. I really hope that won't be the case since it's complicated enough already.

I think you should first simulate two identical colliding discs and no spin to see what the issues are. Then consider two discs of unequal radii, then consider many disks all of unequal radii. Then you will see that squares and polygons can become a nightmare.

The main point of view to adopt is that you don't need an analytical expression to describe the mass trajectories. Since you are doing a simulation, you can use the computer to calculate the next position in equal intervals $\Delta t$ using the straight line trajectory equation $$\mathbf r_{\text{new,before}} =\mathbf r_{\text{old,before}}+\mathbf v_{\text{before}}~\Delta t.\tag{1}$$ If the collision criterion is satisfied within such an interval, you assume that the position remains the same the over that interval but the velocities change to their after the collision values then keep calculating straight line motion trajectories using $$\mathbf r_{\text{new,after}} =\mathbf r_{\text{old,before}}+\mathbf v_{\text{after}}~\Delta t.\tag{2}$$

 

[MatinSAR]

Squares or pentagons? - aagh!

At first, I was about to accept it as my project, but I'm really glad that didn't happen.

Then you will see that squares and polygons can become a nightmare.

It's already clear enough at this point.

The main point of view to adopt is that you don't need an analytical expression to describe the mass trajectories. Since you are doing a simulation, you can use the computer to calculate the next position in equal intervals $\Delta t$ using the straight line trajectory equation $$\mathbf r_{\text{new,before}} =\mathbf r_{\text{old,before}}+\mathbf v_{\text{before}}~\Delta t.\tag{1}$$ If the collision criterion is satisfied within such an interval, you assume that the position remains the same the over that interval but the velocities change to their after the collision values then keep calculating straight line motion trajectories using $$\mathbf r_{\text{new,after}} =\mathbf r_{\text{old,before}}+\mathbf v_{\text{after}}~\Delta t.\tag{2}$$

That's exactly what I did for my first simulation with just two disks. It's an effective method, but it doesn't work for disks with high speed because the position intervals become too large, causing the code to miss collisions. Reducing the time step (dt) to a very small value introduces its own problems.

Thank you, @kuruman, @Steve4Physics, and @erobz, for your time. I'll discuss this with him tomorrow. I hope this issue will be resolved with my professor's assistance tomorrow. However, if any problems persist regarding the physics of the problem, I will reach out again.

 

[haruspex]

Depending on what you want to do with the result of the simulation, there may be an alternative to taking tiny steps of time. For each pair of particles and each particle-boundary pair you can calculate if and when they will next collide, then set your $\Delta t$ as the minimum. This can greatly speed up the simulation if you are not, e.g., creating a video. The result is also more accurate usually.

 

[MatinSAR]

if you are not, e.g., creating a video.

Thank you for your reply. I am creating a gif.

 

[kuruman]

You can do what you want because it's your simulation. If I were simulating this, I would use the following piece of flowchart that requires only the 1-d equations linking the final velocity components to the initial velocity components.
Start with initial positions $\mathbf r_{1i},~\mathbf r_{2i}$ and velocities $\mathbf v_{1i},~\mathbf v_{2i}.$
Calculate the tentative final positions using $\mathbf r_{1f}= \mathbf r_{1i}+\mathbf v_{1i}\Delta t$ and $\mathbf r_{2f}= \mathbf r_{1i}+\mathbf v_{2i}\Delta t$.
If the tentative final positions do not imply a collision, then they become the actual final positions in which case the become the initial positions and you go back to step 1.
If the tentative final positions imply a collision, then the initial positions are kept, and you go back to step 1 with $\mathbf v_{1f},~\mathbf v_{2f}$ becoming $\mathbf v_{1i},~\mathbf v_{2i}.$ I show below how to get these.

A unit vector parallel to the center-to-center direction from mass 1 towards mass 2 is, in terms of the Cartesian axes that you have adopted, $$ \mathbf{\hat n}_{\parallel} =\frac{\mathbf r_{2}-\mathbf r_{1}}{|\mathbf r_{2}-\mathbf r_{1}|} =\frac{ (x_2-x_1)~\mathbf {\hat x}+ (y_2-y_1)~\mathbf {\hat y} } {\left[(x_2-x_1)^2+ (y_2-y_1)^2 \right]^{1/2}}$$ Once you have the parallel unit vector, the perpendicular unit vector can be written as $$\mathbf{\hat n}_{\perp}=\frac{ (y_2-y_1)~\mathbf {\hat x}- (x_2-x_1)~\mathbf {\hat y} } {\left[(x_2-x_1)^2+ (y_2-y_1)^2 \right]^{1/2}}.$$ The velocities before the collision have components in the parallel and perpendicular directions given by, $$ v_{k,\parallel}^{\text{before}} = \mathbf{v}_k^{\text{before}}\cdot\mathbf{\hat n}_{\parallel}~; ~~
v_{k,\perp}^{\text{before}} = \mathbf{v}_k^{\text{before}}\cdot\mathbf{\hat n}_{\perp}~~~~(k=1,2).
$$After the collision, the perpendicular components are the same as before. The parallel components are related to the "before" components as given in the reference in post #2.

What we have done here is calculate the final velocities in a rotated frame where one has a 1d collision along the center-to-center distance which is the direction along which the momentum transfer takes place.

 

[MatinSAR]

You can do what you want because it's your simulation. If I were simulating this, I would use the following piece of flowchart that requires only the 1-d equations linking the final velocity components to the initial velocity components.
Start with initial positions $\mathbf r_{1i},~\mathbf r_{2i}$ and velocities $\mathbf v_{1i},~\mathbf v_{2i}.$
Calculate the tentative final positions using $\mathbf r_{1f}= \mathbf r_{1i}+\mathbf v_{1i}\Delta t$ and $\mathbf r_{2f}= \mathbf r_{1i}+\mathbf v_{2i}\Delta t$.
If the tentative final positions do not imply a collision, then they become the actual final positions in which case the become the initial positions and you go back to step 1.
If the tentative final positions imply a collision, then the initial positions are kept, and you go back to step 1 with $\mathbf v_{1f},~\mathbf v_{2f}$ becoming $\mathbf v_{1i},~\mathbf v_{2i}.$ I show below how to get these.

A unit vector parallel to the center-to-center direction from mass 1 towards mass 2 is, in terms of the Cartesian axes that you have adopted, $$ \mathbf{\hat n}_{\parallel} =\frac{\mathbf r_{2}-\mathbf r_{1}}{|\mathbf r_{2}-\mathbf r_{1}|} =\frac{ (x_2-x_1)~\mathbf {\hat x}+ (y_2-y_1)~\mathbf {\hat y} } {\left[(x_2-x_1)^2+ (y_2-y_1)^2 \right]^{1/2}}$$ Once you have the parallel unit vector, the perpendicular unit vector can be written as $$\mathbf{\hat n}_{\perp}=\frac{ (y_2-y_1)~\mathbf {\hat x}- (x_2-x_1)~\mathbf {\hat y} } {\left[(x_2-x_1)^2+ (y_2-y_1)^2 \right]^{1/2}}.$$ The velocities before the collision have components in the parallel and perpendicular directions given by, $$ v_{k,\parallel}^{\text{before}} = \mathbf{v}_k^{\text{before}}\cdot\mathbf{\hat n}_{\parallel}~; ~~
v_{k,\perp}^{\text{before}} = \mathbf{v}_k^{\text{before}}\cdot\mathbf{\hat n}_{\perp}~~~~(k=1,2).
$$After the collision, the perpendicular components are the same as before. The parallel components are related to the "before" components as given in the reference in post #2.

What we have done here is calculate the final velocities in a rotated frame where one has a 1d collision along the center-to-center distance which is the direction along which the momentum transfer takes place.

Impressive insights! Thank you, @kuruman. I will proceed with this method and see where I get with it.

 

[haruspex]

If the tentative final positions

What if the overlap of the discs is only for part of the way between the two positions? Are you assuming $\Delta t$ is so small we can ignore such glancing blows?

then they become the actual final positions

… except that there are all the other possible collisions these two particles might have before that. Having thought about it some more, I think the approach I described in post #14 avoids some pitfalls. Given that a video is to be generated, I would:

• choose an interframe time, $t_{frame}$
• for each frame {
  • $starttime=clock()$
  • initialise $t_{next}=t_{frame}$
  • while $t_{next}>0$ {
    • $t_{collide}=t_{next}$
    • go through all possible collisions to see if any happen before $t_{collide}$, reducing $t_{collide}$ accordingly
    • compute all positions and velocities after $t_{collide}$ elapses
    • subtract $t_{collide}$ from $t_{next}$
  • }
  • generate frame
  • $now=clock()$
  • if $(starttime+t_{frame}>now)$ sleep($starttime+t_{frame}-now$)
  • // watch out for timestamp wrapping
• }

 

[MatinSAR]

Thank you for your reply @haruspex ...

What if the overlap of the discs is only for part of the way between the two positions? Are you assuming Δt is so small we can ignore such glancing blows?

This problem existed in my initial code too, as mentioned in post #13. If I start with velocities smaller than 100, then vΔt for Δt=0.01 would be smaller than 1. So, I use disks with a radius greater than 1. I understand this is a limitation, but since this is only a university project, I will ask my professor if this velocity limit is a negative point or not.

except that there are all the other possible collisions these two particles might have before that.

I do not understand this part. Could you please explain it more?

Having thought about it some more, I think the approach I described in post #14 avoids some pitfalls

But isn't it more complicated to predict which two particles are going to collide?

 

[haruspex]

I do not understand this part. Could you please explain it more?

But isn't it more complicated to predict which two particles are going to collide?

You have to take the particles in pairs to check if they have collided. According to post #16, if the first pair you pick, A and B, have not collided in time $\Delta t$ then you update their positions on the basis of non collision. But then you find A would have collided with C, so the update is wrong. Likewise, if you find A and B have collided and update them you could then find one of those would have collided with C sooner.
The only way to avoid this is to calculate for every pair (including particle and wall pairs) if and when they would next collide and recompute the whole ensemble at the earliest collision time you find.
The rest of the logic in post #18 is to arrange that the frames are equally spaced in elapsed time.

 

[MatinSAR]

You have to take the particles in pairs to check if they have collided. According to post #16, if the first pair you pick, A and B, have not collided in time Δt then you update their positions on the basis of non collision. But then you find A would have collided with C, so the update is wrong.

Why check for just a pair? We can update the positions of all particles and then check them all in pairs. The issue arises when three or more particles collide. In such cases, I'll print an error message and break the code. Am I missing something?

The only way to avoid this is to calculate for every pair (including particle and wall pairs) if and when they would next collide and recompute the whole ensemble at the earliest collision time you find.

Isn't it difficult to predict when they will collide? Each particle has the chance to collide with every other particle.

 

[haruspex]

Why check for just a pair? We can update the positions of all particles and then check them all in pairs.

If I understand your algorithm, you update all positions as though no collisions occur then check for overlapping objects. Immediately there is the possibility that you overlook a collision because the overlap occurred only briefly during the interval.
Next, how do you update the positions? Do you leave them overlapping or backtrack to where they would have made contact?

The issue arises when three or more particles collide.

How do you detect that?

In such cases, I'll print an error message and break the code. Am I missing something?

So the run is aborted? You could find that happening most of the time, and the runs that succeed will not be a representative sample. Besides, that's poor design when you can so easily make it always work.

Isn't it difficult to predict when they will collide?

It is just a matter of solving a quadratic, no? If the quadratic has no solution then they do not collide.

Each particle has the chance to collide with every other particle.

Sure, but you can calculate if and when any given pair of objects (particles and walls) would collide in the absence of all other objects then determine which happens first.

 

[MatinSAR]

If I understand your algorithm, you update all positions as though no collisions occur then check for overlapping objects. Immediately there is the possibility that you overlook a collision because the overlap occurred only briefly during the interval.

Yes this problem with my algorithm was mentioned before:

This problem existed in my initial code too, as mentioned in post #13. If I start with velocities smaller than 100, then vΔt for Δt=0.01 would be smaller than 1. So, I use disks with a radius greater than 1. I understand this is a limitation, but since this is only a university project, I will ask my professor if this velocity limit is a negative point or not.

Next, how do you update the positions? Do you leave them overlapping or backtrack to where they would have made contact?

For particles that collide, I calculate their velocities after the collision using their past positions and take another step using these velocities. For particles that don’t collide, I use their current positions and velocities and take another step.

How do you detect that?

The code check if any two pairs collide, It might detect more than one, and we cannot choose one pair.

So the run is aborted? You could find that happening most of the time, and the runs that succeed will not be a representative sample. Besides, that's poor design when you can so easily make it always work.

Is that really easy?
What formulas you suggust using to calculate velocities after each collision?

It is just a matter of solving a quadratic, no? If the quadratic has no solution then they do not collide.

By difficult, I meant we have many particles... But I don't understand why quadratic? The velocities are constant.

Sure, but you can calculate if and when any given pair of objects (particles and walls) would collide in the absence of all other objects then determine which happens first.

Great!!! I did not think of it this way. Thank you! I will think on this more.

 

[MatinSAR]

@haruspex I think I understand your point now, or at least I understand a better way to do it.
As you suggusted we predict which particles (or a particle and wall) are going to collide and when. Then we use formulas mentioned by @kuruman for velocities after collisions. Right? Because we know which particles are going to collide we can use their initial positions to use the coordinate in which only in one direction velocity changes. Then after each collision using updated position and velocity of those two particles and we calcluate which particels are going to collide first again.

 

[haruspex]

@haruspex I think I understand your point now, or at least I understand a better way to do it.
As you suggusted we predict which particles (or a particle and wall) are going to collide and when. Then we use formulas mentioned by @kuruman for velocities after collisions. Right? Because we know which particles are going to collide we can use their initial positions to use the coordinate in which only in one direction velocity changes. Then after each collision using updated position and velocity of those two particles and we calcluate which particels are going to collide first again.

Yes.

By difficult, I meant we have many particles...

Whether it is three or 1000, the logic (software) is the same, just takes longer to execute. Ultimately, there would be a limit on how fast the particles could move in the video.

But I don't understand why quadratic? The velocities are constant.

You can write the x(t) and y(t) functions for each particle in the form $x_i(t)=x_{i0}+v_{xi0}t$. The equation $(x_1(t)-x_2(t))^2+(y_1(t)-y_2(t))^2=4R^2$ gives a quadratic in t.

 

[MatinSAR]

Yes.

Thank you for your time and valuable help @haruspex ...
So I'm going to try this algorithm and see where I get with it..

Whether it is three or 1000, the logic (software) is the same, just takes longer to execute. Ultimately, there would be a limit on how fast the particles could move in the video.

You can write the x(t) and y(t) functions for each particle in the form $x_i(t)=x_{i0}+v_{xi0}t$. The equation $(x_1(t)-x_2(t))^2+(y_1(t)-y_2(t))^2=4R^2$ gives a quadratic in t.

I get it.