# <Simulating a Bullet-Bullet Collision Using Python>

Written on

When I come across something fascinating, I often create a model to appreciate it fully. This time, I was captivated by the impressive slow-motion footage of a bullet-bullet collision featured on Smarter Every Day. You should check it out.

I am particularly intrigued by the circular debris that spreads out post-collision. My goal is to replicate this phenomenon using basic physics principles and Python. I will employ Web VPython, which offers convenient 3D objects.

## Basic Two-Ball Collision

Let's start with a straightforward scenario: two point particles approaching each other. How do we simulate their interaction upon collision? Naturally, we will apply the principles of physics.

When these two particles collide, there exists an interaction force between them (I'll refer to them as particles A and B for simplicity). The following diagram illustrates this interaction.

Note that the force exerted by B on A (F_B-A) is equal and opposite to the force from A on B (F_A-B), in accordance with Newton's third law. When a net force acts on an object, its momentum changes (as per the momentum principle).

Given that the forces are equal and opposite for the same duration, the momentum change of A is equal and opposite to that of B, implying that the total momentum remains constant (the principle of conservation of momentum).

Now, how do we model this scenario? How can we compute the forces acting on A and B? One approach is to utilize springs. Imagine if the two balls slightly overlap; we can apply a force proportional to the overlap distance. If there’s no overlap, the force is zero, akin to "springy" balls.

The direction of the force depends on the relative positions of A and B, making this model applicable for both glancing and direct collisions. To calculate the force, let’s denote the radii of the balls as R, with their vector positions as r_A and r_B. We also need an effective spring constant k to determine the force's magnitude.

Here’s the force equation (assuming overlap occurs). Remember, R signifies the radius of the balls.

The plan is to implement a numerical calculation by dividing the motion into small time intervals. During each interval, the following actions will take place:

- Calculate the vector distance from A to B. If this distance's magnitude is less than 2R, the balls overlap, and we need to compute the force.
- Calculate the forces acting on A and B (since they are equal and opposite, we only need to perform this step once).
- Use the forces to update the momentum of both balls (applying the momentum principle).
- Update the positions of both balls based on their momentum (assuming constant velocity).

Let’s translate this into Python. Below is the complete code, but I’ll highlight the key sections.

ballA = sphere(pos=vector(-0.03, 0, 0), radius=R, color=color.yellow) ballB = sphere(pos=vector(0.03, 0, 0), radius=R, color=color.cyan)

This code utilizes Web VPython, allowing us to create objects like sphere(). You must specify the position vector and radius for both balls (now modeled as bullets).

Next, we need to establish the initial conditions for momentum and other parameters.

v0 = 100 k = 50000000 ballA.m = m ballB.m = m ballA.p = m * vector(v0, 0, 0) ballB.p = m * vector(-v0, 0, 0)

Now, let’s focus on the loop that performs the calculations for each time interval.

while t < .6e-3:

rate(100)

ballA.F = vector(0, 0, 0)

ballB.F = vector(0, 0, 0)

r = ballB.pos - ballA.pos

if mag(r) < 2 * R:

ballA.F = -k * (2 * R - mag(r)) * norm(r)

ballB.F = -ballA.F

ballA.p = ballA.p + ballA.F * dt

ballB.p = ballB.p + ballB.F * dt

ballA.pos = ballA.pos + ballA.p * dt / ballA.m

ballB.pos = ballB.pos + ballB.p * dt / ballB.m

t = t + dt

A few notes:

- The rate(100) function controls the animation speed, ensuring no more than 100 iterations occur per second.
- The net force for each ball is set as a property, which will be useful when we introduce more balls.
- The magnitude of a vector is calculated using mag(r), while norm(r) gives the unit vector, both of which are built into VPython.

Here’s the visual output when the simulation runs.

While this simulation shows the balls colliding and bouncing off each other, it doesn’t accurately reflect the rapid bullet collisions.

## Modeling Bullets as Collections of Balls

In a real bullet-bullet collision, the bullets shatter into numerous smaller fragments. To simulate this in Python, we can represent each bullet as a cluster of tiny balls. When the two "bullets" collide, all of these small balls can interact, ideally producing a result similar to the video from Smarter Every Day.

Here’s a rough outline of this new model:

- Break each bullet into N individual balls.
- Randomly distribute these smaller balls within a spherical volume, mimicking the structure of a single bullet.
- As the two bullets (now made up of many tiny balls) approach each other, they all interact with one another.

To avoid redundantly coding for multiple balls (like 40), we will use lists. If you’re not familiar with lists, here’s a brief example.

Next, we must determine how to randomly distribute these tiny balls in a spherical formation. I’ll provide the code to create one bullet-ball and point out the interesting aspects.

# some constants R = (11.5e-3) / 2 # m - radius of bullet-ball m = 0.012 # mass of bullet N = 20 # number of balls dm = m / N # mass of each ball

n = 0 # this is a counter bs = [] # this is the empty list for the balls

rleft = vector(-0.03, 0.005, 0) # starting point for center of bullet-ball

v0 = 100 # initial velocity k = 500000 # spring constant rpiece = .8 * R / (N) ** (1 / 3) # radius of piece based on bullet size and N

The function overlap checks if a new position overlaps with any existing balls in the list:

def overlap(rt, radi, rlist):

# checks if point rt with radius radi overlaps with any object in rlist

for rtemp in rlist:

rdist = rtemp.pos - rt

if mag(rdist) < 2 * radi:

return Truereturn False

The loop below creates N balls:

while n < N:

# create a random vector

rt = R * vector(2 * random() - 1, 2 * random() - 1, 2 * random() - 1)

# check if within sphere and no overlap

if mag(rt) < R and not overlap(rt + rleft, rpiece, bs):

# add a ball to the list

bs.append(sphere(pos=rt + rleft, radius=rpiece, v=vector(v0, 0, 0)))

n += 1

### Remarks:

- An empty list is initialized with bs = [].
- If a random vector falls within the spherical volume and doesn’t overlap with existing balls, a new ball is created and added to the list.
- In Web VPython, the primary random number function is random(), which generates a number between 0 and 1. The transformation 2 * random() - 1 yields a range from -1 to 1.
- The overlap() function checks for overlaps, which is crucial to prevent premature "explosions" of the bullet before collision.

After creating the left bullet, we reset the counter and repeat the process for the second bullet. Below is what it looks like when we run this with 20 balls for each bullet.

Next, we’ll animate the two bullets moving toward each other and colliding. Here’s the second part of the code:

while t < .8e-3:

rate(200)

# calculate the forces

for b in bs:

# reset all net forces to zero

b.F = vector(0, 0, 0)

for b in bs:

for a in bs:

rab = b.pos - a.pos

if mag(rab) < 2 * rpiece:

a.F += -k * (2 * rpiece - mag(rab)) * norm(rab)# update momentum and position

for c in bs:

c.p += c.F * dt

c.pos += c.p * dt / c.m

# update time

t += dt

This loop iterates through all the short time intervals, checking which balls overlap to compute their forces. It is crucial to start with a net force of zero to avoid double-counting interactions.

Once the forces for each ball are determined, momentum and positions are updated accordingly.

Take a look at the results.

What if the bullets collide off-axis? I can shift the left bullet slightly upward and run the simulation again. Here’s the outcome.

This is closer to the effect I aim for, but there’s potential for further enhancement. Before proceeding, let’s verify the momentum conservation. If everything functions correctly, the total momentum of this system should remain constant. Below is a graph of the x-component of momentum (for all balls) against time.

While it may not appear constant, the changes in momentum are minimal compared to the total, which is acceptable.

## Partially Elastic Collisions

In the Smarter Every Day video, there's a scenario where a clump of bullets sticks together during the collision. In the animations above, every ball undergoes perfectly elastic collisions, meaning they never bond.

In a perfectly elastic collision, the total kinetic energy before and after remains constant. For two point particles, kinetic energy is the only form of energy, so it is inherently conserved. However, the collisions above involve springs, which temporarily store elastic potential energy while in contact, returning it as kinetic energy upon separation.

How can we introduce non-elastic behavior in this system? One approach is to model individual balls as collections of even smaller balls, but that would lead to an overly complex simulation. Instead, I propose a method using two distinct spring constants during collisions: a stronger constant for compression and a weaker one for separation.

Consider this uneven spring scenario in terms of energy. A typical spring performs negative work when compressed and positive work while pushing objects apart, leading to zero total work and conserved kinetic energy. However, if the outbound spring is weaker, the overall work becomes negative, resulting in a decrease in kinetic energy, which is precisely our goal.

To implement this in the program, we must determine whether the balls are approaching or moving apart. I’ll utilize the momentum vectors of both balls in conjunction with the vector connecting them.

Imagine two balls with positions and momentum vectors. I will calculate the momentum difference and then take the dot product with the relative position vector. This diagram illustrates the concept.

If the r vector and the momentum difference are aligned, the dot product will be positive, indicating compression. If zero, they are stationary; if negative, they are separating.

Returning to my colliding balls, I’ll incorporate this logic into the code.

r = ballB.pos - ballA.pos dp = ballB.p - ballA.p if mag(r) < 2 * R:

if dot(r, dp) > 0:

k = k2else:

k = k1

Here, k1 represents the compressive spring constant, while k2 is the weaker repelling constant. If k2 equals k1/2, the collision behaves as follows.

While it’s challenging to visually confirm kinetic energy conservation, I’ll print both momentum and kinetic energy before and after the collision.

Momentum remains conserved, but kinetic energy does not—precisely what I wanted. Now, I can apply this adjustment to my multi-body bullet simulation. Below is the updated section of the code (full code available here).

for b in bs:

for a in bs:

rab = b.pos - a.pos

dp = b.p - a.p

if mag(rab) < 2 * rpiece:

if dot(rab, dp) < 0:

k = k1else:

k = k2a.F += -k * (2 * rpiece - mag(rab)) * norm(rab)

Now for the moment you’ve been waiting for: a collision between two bullets, each composed of 50 small balls, with a rebound spring constant that’s half of the compressive constant. Note that the collision isn’t perfectly head-on for added interest.

This simulation is an improvement over the perfectly elastic collision scenario, but it still doesn’t replicate the video’s effect. I suspect achieving a realistic expanding ring of debris will require varying sizes of bullet fragments, with larger pieces remaining at the center. Overall, I'm pleased with the progress so far.