Unexpected behavior of bodies affected by gravity in 2D

Question

I have a simple code that calculates the acceleration of a body every update of my simulation using the Newton's law of universal gravitation and the second law of motion.

def move_bodies():
    def get_dist(b1, b2):
        return math.hypot(b2['x'] - b1['x'], b2['y'] - b1['y'])

    for body in bodies:
        acceleration = 0
        angle = 0

        for other_body in bodies:
            dist = get_dist(body, other_body)

            if body != other_body:
                angle = math.atan2(other_body['y'] - body['y'], other_body['x'] - body['x'])
                acceleration = (GRAVITY_CONSTANT * body['m'] * other_body['m'] / dist ** 2) / body['m']

            body['x'] += math.cos(angle) * body['vx']
            body['vx'] += acceleration

            body['y'] += math.sin(angle) * body['vy']
            body['vy'] += acceleration

when these bodies interact, something unexpected happens when their x/y values are approaching one another -- the bodies move at extremely high speeds and seemingly disappear off the screen.

For example, if we have a body at (200, 100) and another at (200, 300), they do not gravitate towards one another as they would do in a situation where they have unique (x, y) values such as (200, 100) and (300, 400).

Looking at the code, as the y-values of body and other_body approach one another, values get exceedingly large(?)

I can't quite pin down the issue as I'm a beginner. Any suggestions/criticism?

Answer 1

I'm on my phone and can't provide a super detailed answer, but here goes.

This answer will not cover the physics calculation that you use, only the logic.

It looks like there is trash data in your loop.

When evaluating a body with itself, you do not update the angle and the acceleration, it remeains to what was used with a previous body. Whatever is in there is still used with the body influencing itself.

I would expect that bodies don't influence themselves in that loop, so perhaps it should look a bit more like this:

def move_bodies():
    def get_dist(b1, b2):
        return math.hypot(b2['x'] - b1['x'], b2['y'] - b1['y'])

    for body in bodies:
        for other_body in bodies:
            if body == other_body:
                continue

            dist = get_dist(body, other_body)
            angle = math.atan2(other_body['y'] - body['y'], other_body['x'] - body['x'])
            acceleration = (GRAVITY_CONSTANT * body['m'] * other_body['m'] / dist ** 2) / body['m']

            body['x'] += math.cos(angle) * body['vx']
            body['vx'] += acceleration

            body['y'] += math.sin(angle) * body['vy']
            body['vy'] += acceleration

As a final note: this will most likely give inacurate results. As you update the x and the y values of bodies during the loop, you change how the next bodies are influenced by it.

A more accurate approach would be to use the "current position" to do the calculations, and to store it somewhere else for further processing:

def move_bodies():
    def get_dist(b1, b2):
        return math.hypot(b2['x'] - b1['x'], b2['y'] - b1['y'])

    for body in bodies:
        for other_body in bodies:
            if body == other_body:
                continue

            dist = get_dist(body, other_body)
            angle = math.atan2(other_body['y'] - body['y'], other_body['x'] - body['x'])
            acceleration = (GRAVITY_CONSTANT * body['m'] * other_body['m'] / dist ** 2) / body['m']

            body['x_sum'] += math.cos(angle) * body['vx']
            body['vx_sum'] += acceleration

            body['y_sum'] += math.sin(angle) * body['vy']
            body['vy_sum'] += acceleration

    for body in bodies:
        body['x'] = body['x_sum'] / (len(bodies) - 1)
        body['vx'] = body['vx_sum'] / (len(bodies) - 1)

        body['y'] = body['y_sum'] / (len(bodies) - 1)
        body['vy'] = body['vy_sum'] / (len(bodies) - 1)

        body['x_sum'] = 0
        body['vx_sum'] = 0

        body['y_sum'] = 0
        body['vy_sum'] = 0