Векторы C ++, операторы перегрузки

Question

В попытке переучить немного C ++, я пытаюсь написать код, который объединяет проблему N-body. У меня уже был код, который работал просто отлично, но вместо того, чтобы перечислять 3d-позиции в 3 отдельных двойных переменных, я сейчас пытаюсь собрать их в один вектор, одинаковый для скоростей. Первый шаг - перевести всю систему так, чтобы ее центр масс стал равным 0, а также общий импульс. Для этого нужно рассчитать текущий центр масс, сложив, по существу, произведения масс и позиций. Поскольку позиции теперь являются векторами, это приводит к необходимости добавления векторов. Так что я уже перегружен оператором +, а также еще немного. Мой код продолжает падать, хотя.

#include <iostream>
#include <fstream>
#include <list>
#include <cmath>
#include <vector>
#include <string>

using namespace std;

//OVERLOAD SOME OPERATORS TO SUIT VECTORS
//scalar product
double operator*(vector<double> a, vector<double> b)
{
    double result = 0.0;
    for(int i=0; i<a.size(); i++)
    {
        result = result + a[i]*b[i];
    };
    return result;
};
vector<double> operator*(double a, vector<double> b)
{
    vector<double> result;
    for(int i=0; i<b.size(); i++)
    {
        result.push_back(a*b[i]);
    };
    return result;
};
vector<double> operator*(vector<double> a, double b)
{
    vector<double> result;
    for(int i=0; i<a.size(); i++)
    {
        result.push_back(a[i]*b);
    };
    return result;
};
vector<double> operator/(vector<double> a, double b)
{
    vector<double> result;
    for(int i=0; i<a.size(); i++)
    {
        result.push_back(a[i]/b);
    };
    return result;
};
vector<double> operator+(vector<double> a, vector<double> b)
{
    vector<double> result;
    for(int i=0; i<a.size(); i++)
    {
        result.push_back(a[i]+b[i]);
    };
    return result;
};
vector<double> operator-(vector<double> a, vector<double> b)
{
    vector<double> result;
    for(int i=0; i<a.size(); i++)
    {
        result.push_back(a[i]-b[i]);
    };
    return result;
};

//all the data of a planet
struct planet
{
    double mass;
    vector<double> r;
    vector<double> v;
    vector<double> f;
};

//collect all planets in a class, for data storage and calculations
class planets
{
    public:
        planets()
        {
        };
        //add a new planet to the list of planets
        void join(planet newpl)
        {
            this->pl.push_back(newpl);
        };
        //read in a startfile with positions and velocities of planets
        void read_start(string startfile)
        {
            ifstream f;
            f.open(startfile.c_str());
            planet readpl;
            vector<double> r(3);
            vector<double> v(3);
            readpl.mass = 0.0;
            while(f >> r[0] >> r[1] >> r[2] >> v[0] >> v[1] >> v[2] >> readpl.mass)
            {
                this->join(readpl);
            };
            f.close();
        };
        //translate the system such that the center of mass and the total momentum are both 0
        void cm()
        {
            cout << " cm test 0" << endl;
            list<planet>::iterator it_i;
            cout << " cm test 1" << endl;
            int i = 0;
            cout << " cm test 2" << endl;
            vector<double> cmr(3,0.0);
            cout << " cm test 3" << endl;
            vector<double> cmv(3,0.0);
            cout << " cm test 4" << endl;
            double totmass = 0.0;
            cout << " cm test 5" << endl;
            for(it_i=this->pl.begin(); it_i!=this->pl.end(); ++it_i)
            {
                cout << " cm test 6" << endl;
                totmass = totmass + it_i->mass;
                cout << " cm test 7" << endl;
                cmr = it_i->r*it_i->mass + cmr;
                cout << " cm test 8" << endl;
                cmv = cmv + it_i->v*it_i->mass;
                cout << " cm test 9" << endl;
            };
            for(it_i=this->pl.begin(); it_i!=this->pl.end(); ++it_i)
            {
                it_i->r = it_i->r - cmr/totmass;
                it_i->v = it_i->v - cmv/totmass;
            };
        };
    private:
        list<planet> pl;
};

int main()
{
    cout << "test 0" << endl;
    planets solar_system;
    cout << "test 1" << endl;
    solar_system.read_start("planets.start");
    cout << "test 2" << endl;
    solar_system.cm();
    cout << "test 3" << endl;
    return 0;
};

Я попытался уменьшить свой код до необходимых частей. Код, как он есть, воспроизводит мои ошибки. Файл "planets.start" содержит это:

0.0 0.0 0.0 0.0 0.0 0.0 1.989e30
147.09e9 0.0 0.0 0.0 30.29e3 0.0 5.972e24

Это просто данные об орбите Земли и Солнца.

Я включил некоторые тестовые данные, чтобы локализовать проблемный c фрагмент код, и кажется, что cmv = cmv + it_i->v*it_i->mass; это проблема. Теперь мне интересно, почему строка прямо перед cmr = it_i->r*it_i->mass + cmr; работает просто отлично, так как она по сути такая же, просто порядок суммы обратный. Я не реализовал ни конструкторов, ни копировальных конструкторов и чувствую, что это может быть связано с этим. Или мои реализации перегрузок операторов неверны. Если бы кто-то мог мне это объяснить, я был бы благодарен. Ура!