diff --git a/examples/plot.py b/examples/plot.py
index 941a33ff03841be073f9f6f800f15faf2786bb98..008ed265b729d54cbebf839bf85132e998177a10 100644
--- a/examples/plot.py
+++ b/examples/plot.py
@@ -127,7 +127,7 @@ plot(id, errx_new , 'b.')
text(id[-1], 0.18, "B-H: $%5.3f\\pm%5.3f$\n QuickShed: $%5.3f\\pm%5.3f$"%(meanx_bh, stdx_bh, meanx_new, stdx_new), backgroundcolor="w", va="top", ha="right" )
-ylim(-0.2, 0.2)
+ylim(-0.02, 0.02)
xlim(0,id[-1])
grid()
@@ -138,7 +138,7 @@ plot(id, erry_new , 'b.')
text(id[-1], 0.18, "B-H: $%5.3f\\pm%5.3f$\n QuickShed: $%5.3f\\pm%5.3f$"%(meany_bh, stdy_bh, meany_new, stdy_new), backgroundcolor="w", va="top", ha="right" )
-ylim(-0.2, 0.2)
+ylim(-0.02, 0.02)
xlim(0,id[-1])
grid()
@@ -150,7 +150,7 @@ plot(id, errz_bh , 'rx', label="Legacy")
#text(id[-1], 0.18, "B-H: $%5.3f\\pm%5.3f$\n QuickShed: $%5.3f\\pm%5.3f$"%(meanz_bh, stdz_bh, meanz_new, stdz_new), backgroundcolor="w", va="top", ha="right" )
legend(loc="upper right")
-ylim(-0.2, 0.2)
+ylim(-0.02, 0.02)
xlim(0,id[-1])
grid()
diff --git a/examples/test_fmm.c b/examples/test_fmm.c
index 8c73cc8187fa7d55dd59ed14afb7d212d824dc5c..eb56a76104dfd89715d5d85bce72d3d103b9985c 100644
--- a/examples/test_fmm.c
+++ b/examples/test_fmm.c
@@ -37,7 +37,7 @@
/* Some local constants. */
#define cell_pool_grow 1000
-#define cell_maxparts 1
+#define cell_maxparts 128
#define task_limit 1
#define const_G 1 // 6.6738e-8
#define dist_min 0.5 /* Used for legacy walk only */
@@ -255,19 +255,15 @@ void comp_down(struct cell *c) {
double dx=0.f, dy=0.f, dz=0.f;
struct part *parts = c->parts;
-
-
if (c->split) {
- //message("Downpass: %f", c->field_x.F_000);
-
/* Loop over the siblings and propagate accelerations. */
for (cp = c->firstchild; cp != c->sibling; cp = cp->sibling) {
/* Distance over which to shift */
dx = cp->new.M_100 - c->new.M_100;
- dy = cp->new.M_100 - c->new.M_100;
- dz = cp->new.M_100 - c->new.M_100;
+ dy = cp->new.M_010 - c->new.M_010;
+ dz = cp->new.M_001 - c->new.M_001;
/* Propagate the tensors downwards */
shiftAndAddTensor(c->field_x, &cp->field_x, dx, dy, dz);
@@ -294,10 +290,6 @@ void comp_down(struct cell *c) {
ay = const_G * applyFieldAcceleration(c->field_y, dx, dy, dz);
az = const_G * applyFieldAcceleration(c->field_z, dx, dy, dz);
-
-
- //message("Apply to particle: %f %f", ax, c->field_x.F_000);
-
parts[k].a[0] += ax;
parts[k].a[1] += ay;
parts[k].a[2] += az;
@@ -569,12 +561,6 @@ static inline void iact_pair_pc(struct cell *ci, struct cell *cj) {
computeFieldTensors_x( -r_x, -r_y, -r_z, inv_r, cp->new, &cps->field_x );
computeFieldTensors_y( -r_x, -r_y, -r_z, inv_r, cp->new, &cps->field_y );
computeFieldTensors_z( -r_x, -r_y, -r_z, inv_r, cp->new, &cps->field_z );
-
- /* message("r= %f %f %f %f M_000= %f %f F_000= %f", r_x, r_y, r_z, sqrt(r2), */
- /* cps->new.M_000, */
- /* cps->new.M_000 * D_100( r_x, r_y, r_z, inv_r), */
- /* cp->field_x.F_000); */
-
}
}
}
@@ -1260,36 +1246,36 @@ void test_bh(int N, int nr_threads, int runs, char *fileName) {
/* If no input file was specified, generate random particle positions. */
if (fileName == NULL || fileName[0] == 0) {
- /* for (k = 0; k < N; k++) { */
- /* parts[k].id = k; */
- /* parts[k].x[0] = ((double)rand()) / RAND_MAX; */
- /* parts[k].x[1] = ((double)rand()) / RAND_MAX; */
- /* parts[k].x[2] = ((double)rand()) / RAND_MAX; */
- /* parts[k].mass = 1.f;//((double)rand()) / RAND_MAX; */
- /* parts[k].a_legacy[0] = 0.0; */
- /* parts[k].a_legacy[1] = 0.0; */
- /* parts[k].a_legacy[2] = 0.0; */
- /* } */
+ for (k = 0; k < N; k++) {
+ parts[k].id = k;
+ parts[k].x[0] = ((double)rand()) / RAND_MAX;
+ parts[k].x[1] = ((double)rand()) / RAND_MAX;
+ parts[k].x[2] = ((double)rand()) / RAND_MAX;
+ parts[k].mass = 1.f;//((double)rand()) / RAND_MAX;
+ parts[k].a_legacy[0] = 0.0;
+ parts[k].a_legacy[1] = 0.0;
+ parts[k].a_legacy[2] = 0.0;
+ }
- int ii, jj, kk;
- for (ii = 0; ii < 8; ++ii) {
- for (jj = 0; jj < 8; ++jj) {
- for (kk = 0; kk < 8; ++kk) {
+ /* int ii, jj, kk; */
+ /* for (ii = 0; ii < 8; ++ii) { */
+ /* for (jj = 0; jj < 8; ++jj) { */
+ /* for (kk = 0; kk < 8; ++kk) { */
- int index = 64*ii + 8*jj + kk;
-
- parts[index].id = index;
- parts[index].x[0] = ii*0.125 + 0.0625;
- parts[index].x[1] = jj*0.125 + 0.0625;
- parts[index].x[2] = kk*0.125 + 0.0625;
- parts[index].mass = 1.;
- parts[index].a_legacy[0] = 0.f;
- parts[index].a_legacy[1] = 0.f;
- parts[index].a_legacy[2] = 0.f;
-
- }
- }
- }
+ /* int index = 64*ii + 8*jj + kk; */
+
+ /* parts[index].id = index; */
+ /* parts[index].x[0] = ii*0.125 + 0.0625; */
+ /* parts[index].x[1] = jj*0.125 + 0.0625; */
+ /* parts[index].x[2] = kk*0.125 + 0.0625; */
+ /* parts[index].mass = 1.; */
+ /* parts[index].a_legacy[0] = 0.f; */
+ /* parts[index].a_legacy[1] = 0.f; */
+ /* parts[index].a_legacy[2] = 0.f; */
+
+ /* } */
+ /* } */
+ /* } */
/* Otherwise, read them from a file. */
} else {
@@ -1419,7 +1405,7 @@ void test_bh(int N, int nr_threads, int runs, char *fileName) {
message("[check] root CoMx= %f %f", root->legacy.com[0], root->new.M_100);
message("[check] root CoMy= %f %f", root->legacy.com[1], root->new.M_010);
message("[check] root CoMz= %f %f", root->legacy.com[2], root->new.M_001);
-
+#ifdef QUADRUPOLES
message("[check] | %f %f %f |", root->new.M_200, root->new.M_110, root->new.M_101);
message("[check] I = | %f %f %f |", root->new.M_110, root->new.M_020, root->new.M_011);
message("[check] | %f %f %f |", root->new.M_101, root->new.M_011, root->new.M_002);
@@ -1445,6 +1431,7 @@ void test_bh(int N, int nr_threads, int runs, char *fileName) {
message("[check] | %f %f %f |", M_200, M_110, M_101);
message("[check] I = | %f %f %f |", M_110, M_020, M_011);
message("[check] | %f %f %f |", M_101, M_011, M_002);
+#endif
#if ICHECK >= 0
for (i = 0; i < N; ++i)