4

I'm trying to simulate Newtonian gravity with C++ and GLM. I'm using GLM vectors and this equation I found on Wikipedia

I translated this into code:

const static double g_const = 6.67430E-11;
void gravitate(object& obj1,object& obj2) {
    glm::vec2 distance = glm::abs(obj2.position-obj1.position);
    glm::vec2 unitvec = (obj2.position-obj1.position)/distance;
    glm::vec2 force = -g_const*obj1.mass*obj2.mass/(distance*distance)*unitvec;
    obj2.updateForce(force);
}

The object class is structured as this

class object {
public:
    double mass;
    glm::vec2 velocity;
    glm::vec2 position;

    object(double massin, glm::vec2 pos) : mass(massin), position(pos) {;}
    object(double massin, glm::vec2 pos, glm::vec2 vel) : mass(massin), position(pos), velocity(vel) {;}

    void updateForce(glm::vec2 force) {
        velocity+=force/mass;
    }

    void updatePos() {
        position+=velocity;
    }

    void printPos(); // Not used in this example
    void printVel(); // Not used in this example
}

The main function is run like this:

int main() {
    object test1(massOfEarth,
        glm::vec2(0, 0),
        glm::vec2(0,0));
    object test2(massOfmoon,
        glm::vec2(distanceFromEarth,
        distanceFromEarth),
        moonInitalVelocity);
    while (true) {
       gravitate(test1,test2);
       test2.updatePos();
    }
}

When I run this, however, I get this output: GIF Version

NOTE: Sorry for the low frame rate, it's a byproduct of converting a video to a GIF

Does anyone know what's going wrong here?

BTW the full code is available here: https://github.com/ProtoByter/GravityLaw

EDIT 1:

The y force placed on the object at the point that it suddenly speeds up is 3.5079064080620432e+26 N

whilst when it's moving at a constant speed the y force is -2.9308778146655032e+23 N

EDIT 2:

For clarification I'm using the vector form of newton's universal law of gravity, on wikipedia this is defined as a vector quantity whilst glm::distance returns a double so is not suitable

genpfault
  • 51,148
  • 11
  • 85
  • 139
CosmoΓammaByte
  • 78
  • 1
  • 5
  • 1
    _Does anyone know what's going wrong here?_ Maybe, a physics expert can see immediately what's going wrong here. For the other dummies (like me), could you explain what output do you expect contrary to what you get? In general, I would use the [debugger](https://stackoverflow.com/q/25385173/7478597) to compare what the program computes to compare this with the "numbers on paper" (i.e. what it should compute). – Scheff's Cat Nov 21 '20 at 08:55
  • I asked my physics teacher and he says what I've done should work – CosmoΓammaByte Nov 21 '20 at 09:06
  • _that it suddenly speeds up is 3.5079064080620432e+26 N_. So, pause you program as short as possible before that happens and then check the computation in [single-step debugging](https://stackoverflow.com/q/25385173/7478597). The computations with `double` are not completely like you might have learnt in your Math/Physics lessons for Real numbers. They are similar at best and only for a few values (compared to the infinite set of Real numbers). For very small or very large values, you might get strange effects. This is something I would search for second. – Scheff's Cat Nov 21 '20 at 09:17
  • _my physics teacher ... says what I've done should work_: First I would search for implementation errors. (I just recently scrolled the doxygen man. up and down because a reference to a group in my code wasn't working as expected. Finally, I found that I had written `Progess` instead of `Progress` and somehow this is not easy to spot...) – Scheff's Cat Nov 21 '20 at 09:19
  • 2
    You are dividing by zero. In line 2 of `gravitate` when the two objects are close then `distance` becomes zero or close to zero. As you know lim 1/zero -> a huge number. Also note that you are using Manhatan distance and perhaps for your calculations you should be using Euclidian distance (the function to compute it in std C++ is called hypot but I don't know an equivalent in the library you are using) – Aleksander Bobiński Nov 21 '20 at 09:25
  • @CosmoΓammaByte I'm guessing your physics teacher isn't a programmer. No doubt what you **intended** to write would have worked, what you actually wrote does not. You have to find out where the difference between what you intended and what you actually wrote is. This is normal in programming. – john Nov 21 '20 at 09:39
  • [tag:glm] != [tag:glm-math] – genpfault Oct 15 '21 at 16:18

1 Answers1

3

You using wrong function for Calculating the distance between the two objects, You have used glm::abs(gentype const &x) which is just for calculating the absolute value of x, i.e. if x < 0 then it returns -x and if x >= it simply returns as it is.

You can used glm::distance(point1, point2) for calculating the distance between the point1 and point2.

for more information you can check out this link: https://glm.g-truc.net/0.9.4/api/a00129.html and https://glm.g-truc.net/0.9.4/api/a00131.html

Ravi Prakash
  • 313
  • 1
  • 8