Skip to content
Snippets Groups Projects
Commit ff9a358c authored by Matthieu Schaller's avatar Matthieu Schaller
Browse files

Corrected the N^2 interaction loop of the example.

parent 0d848b57
No related branches found
No related tags found
No related merge requests found
...@@ -64,7 +64,6 @@ accx_new = accx_new[rank] ...@@ -64,7 +64,6 @@ accx_new = accx_new[rank]
accy_new = accy_new[rank] accy_new = accy_new[rank]
accz_new = accz_new[rank] accz_new = accz_new[rank]
# # Read Gadget accelerations # # Read Gadget accelerations
# data=loadtxt("particle_dump_gadget.dat") # data=loadtxt("particle_dump_gadget.dat")
# id = data[:,0] # id = data[:,0]
......
...@@ -815,6 +815,7 @@ void legacy_comp_com(struct cell *c, int *countCoMs) { ...@@ -815,6 +815,7 @@ void legacy_comp_com(struct cell *c, int *countCoMs) {
c->legacy.com[2] = 0.0; c->legacy.com[2] = 0.0;
c->legacy.mass = 0.0; c->legacy.mass = 0.0;
} }
} }
/** /**
...@@ -968,14 +969,13 @@ void interact_exact(int N, struct part *parts, int monitor) { ...@@ -968,14 +969,13 @@ void interact_exact(int N, struct part *parts, int monitor) {
/* Some things to store locally. */ /* Some things to store locally. */
double pix[3] = {parts[i].x[0], parts[i].x[1], parts[i].x[2]}; double pix[3] = {parts[i].x[0], parts[i].x[1], parts[i].x[2]};
double mi = parts[i].mass; double mi = parts[i].mass;
double ai[3] = {0.0, 0.0, 0.0};
/* Loop over every other particle. */ /* Loop over every other particle. */
for (j = i + 1; j < N; ++j) { for (j = i + 1; j < N; ++j) {
/* Compute the pairwise distance. */ /* Compute the pairwise distance. */
for (r2 = 0.0, k = 0; k < 3; k++) { for (r2 = 0.0, k = 0; k < 3; k++) {
dx[k] = parts[i].x[k] - pix[k]; dx[k] = parts[j].x[k] - pix[k];
r2 += dx[k] * dx[k]; r2 += dx[k] * dx[k];
} }
...@@ -985,15 +985,11 @@ void interact_exact(int N, struct part *parts, int monitor) { ...@@ -985,15 +985,11 @@ void interact_exact(int N, struct part *parts, int monitor) {
for (k = 0; k < 3; k++) { for (k = 0; k < 3; k++) {
double wdx = w * dx[k]; double wdx = w * dx[k];
parts[j].a_exact[k] += wdx * mi; parts[j].a_exact[k] -= wdx * mi;
ai[k] -= wdx * parts[j].mass; parts[i].a_exact[k] += wdx * parts[j].mass;
} }
} }
/* Store the particle acceleration. */
for (k = 0; k < 3; k++) {
parts[i].a_exact[k] = ai[k];
}
} }
for (i = 0; i < N; ++i) for (i = 0; i < N; ++i)
...@@ -1098,8 +1094,6 @@ void test_bh(int N, int nr_threads, int runs, char *fileName) { ...@@ -1098,8 +1094,6 @@ void test_bh(int N, int nr_threads, int runs, char *fileName) {
root->loc[1] = 0.0; root->loc[1] = 0.0;
root->loc[2] = 0.0; root->loc[2] = 0.0;
root->h = 1.0; root->h = 1.0;
root->h = 1.0;
root->h = 1.0;
root->count = N; root->count = N;
root->parts = parts; root->parts = parts;
cell_split(root, &s); cell_split(root, &s);
...@@ -1109,7 +1103,7 @@ void test_bh(int N, int nr_threads, int runs, char *fileName) { ...@@ -1109,7 +1103,7 @@ void test_bh(int N, int nr_threads, int runs, char *fileName) {
/* Do a N^2 interactions calculation */ /* Do a N^2 interactions calculation */
tic_exact = getticks(); tic_exact = getticks();
// interact_exact( N , parts , ICHECK ); interact_exact( N , parts , ICHECK );
toc_exact = getticks(); toc_exact = getticks();
printf("Exact calculation (1 thread) took %lli (= %e) ticks\n", printf("Exact calculation (1 thread) took %lli (= %e) ticks\n",
...@@ -1210,7 +1204,7 @@ void test_bh(int N, int nr_threads, int runs, char *fileName) { ...@@ -1210,7 +1204,7 @@ void test_bh(int N, int nr_threads, int runs, char *fileName) {
/* Dump the particles to a file */ /* Dump the particles to a file */
file = fopen("particle_dump.dat", "w"); file = fopen("particle_dump.dat", "w");
fprintf(file, "# a_exact.x a_exact.y a_exact.z a_legacy.x " fprintf(file, "# x y z a_exact.x a_exact.y a_exact.z a_legacy.x "
"a_legacy.y a_legacy.z a_new.x a_new.y a_new.z\n"); "a_legacy.y a_legacy.z a_new.x a_new.y a_new.z\n");
for (k = 0; k < N; ++k) for (k = 0; k < N; ++k)
fprintf(file, "%d %e %e %e %e %e %e %e %e %e %e %e %e\n", parts[k].id, fprintf(file, "%d %e %e %e %e %e %e %e %e %e %e %e %e\n", parts[k].id,
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment