N-body problem: Difference between revisions

Content added Content deleted
Line 428: Line 428:
===Further work===
===Further work===
And as in every experiment report/journal paper/whatever, time for some philosophizing, a more rigorous implementation will require the bodies to be realistic, this means they will not be point particles and hence the dynamics will be much more complex. Tidal effects, collisions, orbits will require much more computation. When I say realistic, I don't mean rigid bodies because even rigid bodies are an idealization of reality. And all of this, even with classical mechanics, employing the Newtonian version of Gravitation will be a worthy task for any competent Physicist, let alone a Programmer. Introduce the General Theory of Relativity and everything changes. I can go on and on, but this is a coding, not a rambling site. I will end by saying that the Inverse Square Law, a name often used for Newtonian Gravitation, is just one class of forces studied in Dynamics. The force relation can take any form, the N-body problem can therefore be solved in (countably/uncountably ?) infinite ways.
And as in every experiment report/journal paper/whatever, time for some philosophizing, a more rigorous implementation will require the bodies to be realistic, this means they will not be point particles and hence the dynamics will be much more complex. Tidal effects, collisions, orbits will require much more computation. When I say realistic, I don't mean rigid bodies because even rigid bodies are an idealization of reality. And all of this, even with classical mechanics, employing the Newtonian version of Gravitation will be a worthy task for any competent Physicist, let alone a Programmer. Introduce the General Theory of Relativity and everything changes. I can go on and on, but this is a coding, not a rambling site. I will end by saying that the Inverse Square Law, a name often used for Newtonian Gravitation, is just one class of forces studied in Dynamics. The force relation can take any form, the N-body problem can therefore be solved in (countably/uncountably ?) infinite ways.

=={{header|C++}}==
{{trans|Java}}
<lang cpp>#include <algorithm>
#include <iomanip>
#include <iostream>
#include <fstream>
#include <string>
#include <vector>

class Vector {
private:
double px, py, pz;

public:
Vector() : px(0.0), py(0.0), pz(0.0) {
// empty
}

Vector(double x, double y, double z) : px(x), py(y), pz(z) {
// empty
}

double mod() const {
return sqrt(px*px + py * py + pz * pz);
}

Vector operator+(const Vector& rhs) const {
return Vector(px + rhs.px, py + rhs.py, pz + rhs.pz);
}

Vector operator-(const Vector& rhs) const {
return Vector(px - rhs.px, py - rhs.py, pz - rhs.pz);
}

Vector operator*(double s) const {
return Vector(px*s, py*s, pz*s);
}

bool operator==(const Vector& rhs) const {
return px == rhs.px
&& py == rhs.py
&& pz == rhs.pz;
}

friend std::istream& operator>>(std::istream&, Vector&);
friend std::ostream& operator<<(std::ostream&, Vector&);
};

std::istream& operator>>(std::istream& in, Vector& v) {
return in >> v.px >> v.py >> v.pz;
}

std::ostream& operator<<(std::ostream& out, Vector& v) {
auto precision = out.precision();
auto width = out.width();
out << std::fixed << std::setw(width) << std::setprecision(precision) << v.px << " ";
out << std::fixed << std::setw(width) << std::setprecision(precision) << v.py << " ";
out << std::fixed << std::setw(width) << std::setprecision(precision) << v.pz;
return out;
}

const Vector ORIGIN{ 0.0, 0.0, 0.0 };

class NBody {
private:
double gc;
int bodies;
int timeSteps;
std::vector<double> masses;
std::vector<Vector> positions;
std::vector<Vector> velocities;
std::vector<Vector> accelerations;

void resolveCollisions() {
for (int i = 0; i < bodies; ++i) {
for (int j = i + 1; j < bodies; ++j) {
if (positions[i] == positions[j]) {
std::swap(velocities[i], velocities[j]);
}
}
}
}

void computeAccelerations() {
for (int i = 0; i < bodies; ++i) {
accelerations[i] = ORIGIN;
for (int j = 0; j < bodies; ++j) {
if (i != j) {
double temp = gc * masses[j] / pow((positions[i] - positions[j]).mod(), 3);
accelerations[i] = accelerations[i] + (positions[j] - positions[i]) * temp;
}
}
}
}

void computeVelocities() {
for (int i = 0; i < bodies; ++i) {
velocities[i] = velocities[i] + accelerations[i];
}
}

void computePositions() {
for (int i = 0; i < bodies; ++i) {
positions[i] = positions[i] + velocities[i] + accelerations[i] * 0.5;
}
}

public:
NBody(std::string& fileName) {
using namespace std;

ifstream ifs(fileName);
if (!ifs.is_open()) {
throw runtime_error("Could not open file.");
}

ifs >> gc >> bodies >> timeSteps;

masses.resize(bodies);
positions.resize(bodies);
fill(positions.begin(), positions.end(), ORIGIN);
velocities.resize(bodies);
fill(velocities.begin(), velocities.end(), ORIGIN);
accelerations.resize(bodies);
fill(accelerations.begin(), accelerations.end(), ORIGIN);

for (int i = 0; i < bodies; ++i) {
ifs >> masses[i] >> positions[i] >> velocities[i];
}

cout << "Contents of " << fileName << '\n';
cout << gc << ' ' << bodies << ' ' << timeSteps << '\n';
for (int i = 0; i < bodies; ++i) {
cout << masses[i] << '\n';
cout << positions[i] << '\n';
cout << velocities[i] << '\n';
}
cout << "\nBody : x y z | vx vy vz\n";
}

int getTimeSteps() {
return timeSteps;
}

void simulate() {
computeAccelerations();
computePositions();
computeVelocities();
resolveCollisions();
}

friend std::ostream& operator<<(std::ostream&, NBody&);
};

std::ostream& operator<<(std::ostream& out, NBody& nb) {
for (int i = 0; i < nb.bodies; ++i) {
out << "Body " << i + 1 << " : ";
out << std::setprecision(6) << std::setw(9) << nb.positions[i];
out << " | ";
out << std::setprecision(6) << std::setw(9) << nb.velocities[i];
out << '\n';
}
return out;
}

int main() {
std::string fileName = "nbody.txt";
NBody nb(fileName);

for (int i = 0; i < nb.getTimeSteps(); ++i) {
std::cout << "\nCycle " << i + 1 << '\n';
nb.simulate();
std::cout << nb;
}

return 0;
}</lang>
{{out}}
<pre>Contents of nbody.txt
0.01 3 20
1
0.000000 0.000000 0.000000
0.010000 0.000000 0.000000
0.100000
1.000000 1.000000 0.000000
0.000000 0.000000 0.020000
0.001000
0.000000 1.000000 1.000000
0.010000 -0.010000 -0.010000

Body : x y z | vx vy vz

Cycle 1
Body 1 : 0.010177 0.000179 0.000002 | 0.010354 0.000357 0.000004
Body 2 : 0.998230 0.998232 0.020002 | -0.003539 -0.003536 0.020004
Body 3 : 0.010177 0.988232 0.988055 | 0.010354 -0.013536 -0.013889

Cycle 2
Body 1 : 0.020709 0.000718 0.000011 | 0.010710 0.000721 0.000014
Body 2 : 0.992907 0.992896 0.039971 | -0.007109 -0.007138 0.019935
Body 3 : 0.020717 0.972888 0.972173 | 0.010727 -0.017153 -0.017876

Cycle 3
Body 1 : 0.031600 0.001625 0.000034 | 0.011072 0.001094 0.000033
Body 2 : 0.983985 0.983910 0.059834 | -0.010735 -0.010835 0.019790
Body 3 : 0.031643 0.953868 0.952235 | 0.011125 -0.020886 -0.021999

Cycle 4
Body 1 : 0.042858 0.002913 0.000081 | 0.011443 0.001481 0.000060
Body 2 : 0.971393 0.971163 0.079509 | -0.014448 -0.014659 0.019561
Body 3 : 0.042981 0.931039 0.928087 | 0.011552 -0.024772 -0.026299

Cycle 5
Body 1 : 0.054492 0.004595 0.000160 | 0.011826 0.001884 0.000097
Body 2 : 0.955030 0.954509 0.098909 | -0.018278 -0.018649 0.019238
Body 3 : 0.054766 0.904225 0.899522 | 0.012018 -0.028857 -0.030829

Cycle 6
Body 1 : 0.066517 0.006691 0.000281 | 0.012224 0.002308 0.000145
Body 2 : 0.934759 0.933760 0.117931 | -0.022265 -0.022849 0.018806
Body 3 : 0.067040 0.873197 0.866280 | 0.012530 -0.033199 -0.035655

Cycle 7
Body 1 : 0.078950 0.009225 0.000456 | 0.012642 0.002759 0.000206
Body 2 : 0.910400 0.908677 0.136456 | -0.026454 -0.027316 0.018244
Body 3 : 0.079856 0.837662 0.828023 | 0.013101 -0.037871 -0.040861

Cycle 8
Body 1 : 0.091815 0.012227 0.000702 | 0.013086 0.003245 0.000284
Body 2 : 0.881722 0.878958 0.154340 | -0.030902 -0.032122 0.017523
Body 3 : 0.093281 0.797239 0.784313 | 0.013749 -0.042975 -0.046559

Cycle 9
Body 1 : 0.105140 0.015737 0.001035 | 0.013564 0.003775 0.000383
Body 2 : 0.848429 0.844216 0.171401 | -0.035684 -0.037362 0.016600
Body 3 : 0.107405 0.751427 0.734579 | 0.014498 -0.048649 -0.052908

Cycle 10
Body 1 : 0.118964 0.019805 0.001481 | 0.014085 0.004362 0.000509
Body 2 : 0.810137 0.803953 0.187408 | -0.040900 -0.043166 0.015414
Body 3 : 0.122346 0.699554 0.678056 | 0.015384 -0.055097 -0.060138

Cycle 11
Body 1 : 0.133337 0.024498 0.002071 | 0.014662 0.005025 0.000672
Body 2 : 0.766343 0.757509 0.202050 | -0.046687 -0.049720 0.013868
Body 3 : 0.138268 0.640690 0.613685 | 0.016460 -0.062633 -0.068603

Cycle 12
Body 1 : 0.148327 0.029907 0.002851 | 0.015317 0.005792 0.000888
Body 2 : 0.716377 0.703998 0.214889 | -0.053246 -0.057302 0.011810
Body 3 : 0.155406 0.573482 0.539941 | 0.017816 -0.071782 -0.078886

Cycle 13
Body 1 : 0.164025 0.036157 0.003887 | 0.016079 0.006709 0.001184
Body 2 : 0.659310 0.642172 0.225282 | -0.060887 -0.066351 0.008976
Body 3 : 0.174112 0.495836 0.454475 | 0.019596 -0.083511 -0.092045

Cycle 14
Body 1 : 0.180564 0.043437 0.005286 | 0.017000 0.007852 0.001613
Body 2 : 0.593807 0.570186 0.232208 | -0.070119 -0.077621 0.004875
Body 3 : 0.194929 0.404136 0.353320 | 0.022038 -0.099890 -0.110265

Cycle 15
Body 1 : 0.198150 0.052049 0.007234 | 0.018171 0.009372 0.002283
Body 2 : 0.517817 0.485100 0.233878 | -0.081861 -0.092550 -0.001535
Body 3 : 0.218605 0.290860 0.228583 | 0.025314 -0.126661 -0.139210

Cycle 16
Body 1 : 0.217126 0.062542 0.010117 | 0.019781 0.011614 0.003484
Body 2 : 0.427899 0.381659 0.226654 | -0.097974 -0.114332 -0.012913
Body 3 : 0.244268 0.131956 0.057562 | 0.026013 -0.191148 -0.202831

Cycle 17
Body 1 : 0.238346 0.076539 0.015221 | 0.022658 0.016380 0.006723
Body 2 : 0.317489 0.248502 0.200967 | -0.122846 -0.151982 -0.038461
Body 3 : 0.075592 -0.559591 -0.487315 | -0.363366 -1.191945 -0.886924

Cycle 18
Body 1 : 0.263123 0.097523 0.026918 | 0.026898 0.025587 0.016672
Body 2 : 0.173428 0.050424 0.112716 | -0.165275 -0.244174 -0.138041
Body 3 : -0.286241 -1.745597 -1.369528 | -0.360299 -1.180066 -0.877501

Cycle 19
Body 1 : 0.270854 0.113045 0.061923 | -0.011436 0.005457 0.053339
Body 2 : 0.199821 -0.093105 -0.208666 | 0.218061 -0.042882 -0.504723
Body 3 : -0.646318 -2.924909 -2.246453 | -0.359856 -1.178559 -0.876350

Cycle 20
Body 1 : 0.258572 0.116046 0.112038 | -0.013129 0.000544 0.046890
Body 2 : 0.426346 -0.111425 -0.681150 | 0.234987 0.006241 -0.440245
Body 3 : -1.006089 -4.103186 -3.122591 | -0.359686 -1.177995 -0.875924</pre>


=={{header|C#|C sharp}}==
=={{header|C#|C sharp}}==
Retrieved from "https://rosettacode.org/wiki/N-body_problem"