Ball bounce physics, system seems to gain energy (bounce higher)

Question

I'm attempting to model the bouncing of balls in a pattern in python with pygame. Something is causing the physics to be incorrect - the balls GAIN energy, i.e. they bounce fractionally higher over time. (I have included a 'speed factor' which can be increase to see the effect I describe.)

Here is my code:

import pygame
import math

# Window dimensions
WIDTH = 800
HEIGHT = 600

# Circle properties
RADIUS = 5
NUM_BALLS = 20
GRAVITY = 9.81  # Gravitational acceleration in m/s²
SPEED_FACTOR = 10  # Speed multiplier for animation

# Circle class
class Circle:
    def __init__(self, x, y, vx, vy, color):
        self.x = x
        self.y = y
        self.vx = vx
        self.vy = vy
        self.color = color

    def update(self, dt):
        # Update positions
        self.x += self.vx * dt
        self.y += self.vy * dt

        # Apply gravity
        self.vy += GRAVITY * dt

        # Bounce off walls
        if self.x - RADIUS < 0 or self.x + RADIUS > WIDTH:
            self.vx *= -1
            self.x = max(RADIUS, min(WIDTH - RADIUS, self.x))  # Clamp x within bounds
        if self.y - RADIUS < 0 or self.y + RADIUS > HEIGHT:
            self.vy *= -1
            self.y = max(RADIUS, min(HEIGHT - RADIUS, self.y))  # Clamp y within bounds

    def draw(self, screen):
        pygame.draw.circle(screen, self.color, (int(self.x), int(self.y)), RADIUS)

# Initialize Pygame
pygame.init()
screen = pygame.display.set_mode((WIDTH, HEIGHT))

circles = []

# Calculate circle arrangement
circle_radius = RADIUS * 2  # Diameter of an individual ball
circle_diameter = NUM_BALLS * circle_radius  # Diameter of the circle arrangement
circle_center_x = WIDTH // 2
circle_center_y = HEIGHT // 2
angle_increment = 2 * math.pi / NUM_BALLS

for i in range(NUM_BALLS):
    angle = i * angle_increment
    x = circle_center_x + math.cos(angle) * circle_diameter / 2
    y = circle_center_y + math.sin(angle) * circle_diameter / 2
    vx = 0
    vy = 0
    hue = i * (360 // NUM_BALLS)  # Calculate hue value based on the number of balls
    color = pygame.Color(0)
    color.hsla = (hue, 100, 50, 100)  # Set color using HSL color space
    circles.append(Circle(x, y, vx, vy, color))

# Game loop
running = True
clock = pygame.time.Clock()
prev_time = pygame.time.get_ticks()  # Previous frame time

while running:
    for event in pygame.event.get():
        if event.type == pygame.QUIT:
            running = False

    current_time = pygame.time.get_ticks()
    dt = (current_time - prev_time) / 1000.0  # Time elapsed in seconds since the last frame
    dt *= SPEED_FACTOR  # Multiply dt by speed factor

    prev_time = current_time  # Update previous frame time

    screen.fill((0, 0, 0))  # Clear the screen

    for circle in circles:
        circle.update(dt)
        circle.draw(screen)

    pygame.display.flip()  # Update the screen

pygame.quit()

I don't fully understand how the gravity factor works, but assume it's just an acceleration. Where is the extra energy coming in to the system?

Answer 1

When you 'clamp' the ball to the position on the edge, you're moving it up slightly. Theoretically, a perfect ball bouncing like this would have the same magnitude of velocity going up and down for the same y-coordinate, but you change the position for the same velocity. On the way back up, the ball gets a head-start, and can go ever so slightly higher.

To solve this, you could do some sort of calculation to determine how much velocity would be lost to gravity in that small distance and diminish your vy accordingly, but there's probably a better way to do it

EDIT

another solution would be to leave the y coordinate as-is and simply change the shape of the ball into an appropriately sized ellipse so that it falls within the boundary. Then change it back into a circle when it no longer intersects a boundary. you would need to have sufficiently short time steps, and a limit on velocity or else your ellipses might become very weird.

Answer 2

Your problem I think is the fact that you are doing a perfectly elastic collision. What I mean is that your vertical velocity is just multiplied by a factor of -1 on hitting the edges. What this does is that ball never looses vertical velocity.

Try changing the vy code to this:

if self.y - RADIUS < 0 or self.y + RADIUS >= HEIGHT: # You don't want it to equal the height. It should be strictly less.
    self.vy *= -0.95 # Instead of a factor of 1. This will dampen the ball thereby reducing the max height it can go.
    self.y = max(RADIUS, min(HEIGHT - RADIUS, self.y))  # Clamp y within bounds

Answer 3

alex_danielssen's answer touches on the reason your simulation gains energy

When you 'clamp' the ball to the position on the edge, you're moving it up slightly.

Consider what happens when a bounce is triggered at the floor (this applies at any wall, but let's just think about the floor for now):

Suppose that when the bounce was triggered because self.y - RADIUS = -1 (which is < 0), we had self.vy = -10. At this point in time, self.y was RADIUS - 1. After you processed the bounce, you set self.vy = 10, and self.y = RADIUS. The ball was raised 1 unit above its previous position so it has more potential energy, but it still has the same kinetic energy as it had when it was lower. Here is your magical gain of energy.

To fix this, you just need to correctly calculate its velocity from energy conservation.

if self.y - RADIUS <= 0 or self.y + RADIUS >= HEIGHT:
    # Energy = potential.     + kinetic          ( per unit mass )
    energy = GRAVITY * self.y + 0.5 * self.vy**2
    # Update location
    self.y = max(RADIUS, min(HEIGHT - RADIUS, self.y))  # Clamp y within bounds
    # Direction of vy (1 = up, -1 = down)
    vy_direction = self.vy / abs(self.vy)
    # Recalculate velocity from energy and new location, flip sign
    self.vy = -vy_direction * (2 * (energy - GRAVITY * self.y))**0.5

Answer 4

I found that applying half of the gravity before the collision check and half of it after mostly fixed the problem on my end.

    def update(self, dt):
    # Update positions
    self.x += self.vx * dt
    self.y += self.vy * dt

    # Apply gravity step 1
    self.vy += GRAVITY/2 * dt

    # Bounce off walls
    if self.x - RADIUS < 0 or self.x + RADIUS > WIDTH:
        self.vx *= -1
        self.x = max(RADIUS, min(WIDTH - RADIUS, self.x))  # Clamp x within bounds
    if self.y - RADIUS < 0 or self.y + RADIUS > HEIGHT:
        self.vy *= -1
        self.y = max(RADIUS, min(HEIGHT - RADIUS, self.y))  # Clamp y within bounds

    # Apply gravity step 2
    self.vy += GRAVITY/2 * dt