gravityMath.c++

#include <iostream>

using namespace std;

#define self (*this)

template <class T>
struct vec
{
    T x, y, z;

    vec operator + (vec r) {
        return { x + r.x, y + r.y, z + r.z };
    }

    vec operator - (vec r) {
        return { x - r.x, y - r.y, z - r.z };
    }

    vec operator * (vec r) {
        return { x * r.x, y * r.y, z * r.z };
    }

    vec operator / (vec r) {
        return { x / r.x, y / r.y, z / r.z };
    }

    static vec diag(T t) {
        return { t, t, t };
    }

    T sum() {
        return x + y + z;
    }
};

struct body {
    float mass;

    vec<float> position;
    vec<float> velocity;

    vec<float> accelerationTo(body other) {
        auto distance = self.position - other.position;
        auto lenSqr   = (distance * distance * distance).sum();
        auto accel    = other.mass / lenSqr;

        return distance * vec<float>::diag(accel);
    }
};

body sun = { 1, vec<float>::diag(0), vec<float>::diag(0) };

int main() {
    vec<float> a = vec<float>::diag(0);

    for (float i = 1; i < 1000000; i++) {
        body probe = 
            { 1
            , { i, 0, 0 }
            , vec<float>::diag(0)
            };

        a = a + probe.accelerationTo(sun);
    }

    cout << a.x << ", " << a.y << ", " << a.z << endl;
}