2
$\begingroup$

I have set up Newtonian Gravity in my Game Engine, allowing me to simulate the gravitational attraction between celestial bodies.

I have the following variables defined:

btScalar whiteHoleMass = 5e+28;
btScalar whiteHoleRadius = 500'000.0;
btScalar sunRadius = 100'000;
btScalar sunOrbitRadius = 50000;
btScalar sunMass = 5e+24;
btScalar sunOrbitTime = 1800; // Orbit time is in seconds
btScalar planetRadius = 50'000;
btScalar planetMass = 5e+21;
btScalar planetOrbitTime = 900;

  • I add the white hole at position: {0, whiteHoleRadius, 0}

  • I add the sun at position: {sunOrbitRadius, whiteHoleRadius, 0}

  • I add the planet at position: {planetOrbitRadius + sunOrbitRadius, whiteHoleRadius, 0}

  • I add 10 orbiting satellites at position: {thisOrbitalHeight + planetOrbitRadius + sunOrbitRadius, whiteHoleRadius, 0}

sunOrbitRadius, planetOrbitRadius and thisOrbitalHeight are calculated with the following functions:

    btScalar sunOrbitRadius = 0;
    btScalar sunOrbitVelocity = 0;
    calculateOrbit(sunOrbitTime, whiteHoleMass, sunOrbitRadius, sunOrbitVelocity);
    btScalar planetOrbitRadius = 0;
    btScalar planetOrbitVelocity = 0;
    calculateOrbit(planetOrbitTime, sunMass, planetOrbitRadius, planetOrbitVelocity);
    btScalar orbitalSatelliteOrbitVelocity = 0;
    btScalar thisOrbitalHeight = 0;
    calculateOrbit(currentOrbitalLength, planetMass, thisOrbitalHeight, orbitalSatelliteOrbitVelocity);

and the calculateOrbit function is defined as:

void calculateOrbit(const btScalar &orbitTimeSecs, const btScalar &centralMassKg, btScalar &orbitRadius, btScalar &orbitVelocity)
{
    orbitRadius = std::cbrt((Common::G * centralMassKg * orbitTimeSecs * orbitTimeSecs) / (4 * Common::Pi * Common::Pi));
    orbitVelocity = std::sqrt(Common::G * centralMassKg / orbitRadius);
};

I then set the Linear Velocity of the Sun, Planet and Satellites like so:

sunPointer->rigidBody->setLinearVelocity({0, 0, sunOrbitVelocity});
planetPointer->rigidBody->setLinearVelocity({0, 0, planetOrbitVelocity + sunOrbitVelocity});
orbitalPointer->rigidBody->setLinearVelocity({0, 0, orbitalSatelliteOrbitVelocity + planetOrbitVelocity + sunOrbitVelocity});

The linear velocity is added in the Z dimension.

I would have thought this was sufficient to keep Satellite/Planet/Sun in constant orbit, however 45 seconds into running the simulation the satellites fall out of orbit and crash into the planet.

How can I set up orbits of multiple bodies so they are stable?

I calculate gravity between objects like so:

for (uInteger64 i = 1; i <= entityMap._size; i++)
{
    btVector3 accelerations(0, 0, 0);
    auto &entityA = *entityMap._data[i]->value.pointer;
    auto &bodyA = *entityA.rigidBody;
    auto bodyAMass = bodyA.getMass();
    auto &positionA = bodyA.getCenterOfMassPosition();
    if (bodyA.getActivationState() == ACTIVE_TAG)
    {
        for (uInteger64 j = 1; j <= entityMap._size; j++)
        {
            if (i == j)
            {
                continue; // Skip self
            }
            auto &entityB = *entityMap._data[j]->value.pointer;
            auto &bodyB = *entityB.rigidBody;
            auto bodyBMass = bodyB.getMass();
            auto &positionB = bodyB.getCenterOfMassPosition();
            btVector3 vectorAB = positionB - positionA;
            btScalar distanceAB = vectorAB.length();
            btVector3 directionAB = vectorAB.normalized();
            Floating128 forceMagnitude = (Floating128)Common::G * ((Floating128)bodyAMass * (Floating128)bodyBMass) / ((Floating128)distanceAB * (Floating128)distanceAB);
            btVector3 forceAB = directionAB * forceMagnitude;
            accelerations += forceAB / bodyAMass;
        }
        entityA.applyCustomAcceleration(accelerations);
    }
}

and applyCustomAcceleration is defined as:

void IEntity::applyCustomAcceleration(const btVector3 &acceleration)
{
    if (!rigidBody || rigidBody->isStaticOrKinematicObject())
    {
        return;
    }
    if (rigidBody)
    {
        btVector3 force = acceleration * rigidBody->getMass();
        rigidBody->applyCentralForce(force);
    }
};
Improve this question
$\endgroup$
3
  • 4
    $\begingroup$ Could you perhaps write some equations describing how you are updating your position and velocity at each time step? That is much easier for readers to parse than game engine code. $\endgroup$
    – whpowell96
    Commented Feb 28 at 17:47
  • 1
    $\begingroup$ Is the step size small enough that the satellite orbits are sufficiently sampled? If so, are the satellite orbits actually circular or narrow ellipses? If the periapsis is below the surface of the planet, you get a case of litho-braking. $\endgroup$
    – Lutz Lehmann
    Commented Feb 28 at 18:10
  • $\begingroup$ Crossposted from physics.stackexchange.com/q/804284/2451 $\endgroup$
    – Qmechanic
    Commented Feb 28 at 22:53

1 Answer 1

Reset to default
7
$\begingroup$

  1. Start with fewer bodies, can you simulate two masses with a spring or an earth/moon system?

  2. Track the total kinetic energy and potential energy of the system to verify that your integrator is (resonably) valid.

  3. Check your units. Are you using SI standard units everywhere, or did an astronomical unit slip into the calculation?

  4. Plot the objects positions in 2D/3D

Improve this answer
$\endgroup$

Not the answer you're looking for? Browse other questions tagged or ask your own question.