diff --git a/src/fof.c b/src/fof.c
index dcae2e1569180064e1c965d8640310041889257b..a3282615869e422cb213954b0cd95083241aa70e 100644
--- a/src/fof.c
+++ b/src/fof.c
@@ -288,7 +288,7 @@ static void rec_fof_search_pair_foreign(struct cell *ci, struct cell *cj, struct
   /* Perform FOF search between pairs of cells that are within the linking
    * length and not the same cell. */
   else if (ci != cj)
-    fof_search_pair_cells_foreign_new(s, ci, cj, send_count, fof_send);
+    fof_search_pair_cells_foreign(s, ci, cj, send_count, fof_send);
   else if (ci == cj) error("Pair FOF called on same cell!!!");
 
 }
@@ -511,163 +511,6 @@ void fof_search_pair_cells_foreign(struct space *s, struct cell *ci, struct cell
   const int periodic = s->periodic;
   double diff[3];
   
-  if(ci_local) {
-
-    for (int k = 0; k < 3; k++) {
-      diff[k] = cj->loc[k] - ci->loc[k];
-      if (periodic && diff[k] < -dim[k] / 2)
-        shift[k] = dim[k];
-      else if (periodic && diff[k] > dim[k] / 2)
-        shift[k] = -dim[k];
-      else
-        shift[k] = 0.0;
-      diff[k] += shift[k];
-    }
-
-    /* Loop over particles and find which particles belong in the same group. */
-    for (size_t i = 0; i < count_i; i++) {
-
-      struct gpart *pi = &gparts_i[i];
-      const double pix = pi->x[0] - shift[0];
-      const double piy = pi->x[1] - shift[1];
-      const double piz = pi->x[2] - shift[2];
-
-      /* Find the root of pi. */
-      int root_i = fof_find(offset_i[i] - node_offset, group_index);
-
-      for (size_t j = 0; j < count_j; j++) {
-
-        struct gpart *pj = &gparts_j[j];
-        const double pjx = pj->x[0];
-        const double pjy = pj->x[1];
-        const double pjz = pj->x[2];
-
-        /* Compute pairwise distance, remembering to account for boundary
-         * conditions. */
-        float dx[3], r2 = 0.0f;
-        dx[0] = pix - pjx;
-        dx[1] = piy - pjy;
-        dx[2] = piz - pjz;
-
-        for (int k = 0; k < 3; k++) r2 += dx[k] * dx[k];
-
-        /* Hit or miss? */
-        if (r2 < l_x2) {
-          //message("Particle %lld in range of foreign particle %lld. Group ID: %d root_i: %d", pi->id_or_neg_offset, pj->id_or_neg_offset, group_index[root_i - node_offset], root_i);
-
-          /* Store the particle IDs and group ID for communication. */
-          fof_send[*send_count].local_pid = pi->id_or_neg_offset;        
-          fof_send[*send_count].foreign_pid = pj->id_or_neg_offset;        
-          fof_send[*send_count].root_i = root_i;
-          (*send_count)++;
-          local_send_count++;
-
-          //fprintf(fof_file, "  %7zu %7zu %7d\n", pi->id_or_neg_offset, pj->id_or_neg_offset, root_i);
-        }
-      }
-    }
-
-  }
-
-  if(cj_local) {
-
-    for (int k = 0; k < 3; k++) {
-      diff[k] = ci->loc[k] - cj->loc[k];
-      if (periodic && diff[k] < -dim[k] / 2)
-        shift[k] = dim[k];
-      else if (periodic && diff[k] > dim[k] / 2)
-        shift[k] = -dim[k];
-      else
-        shift[k] = 0.0;
-      diff[k] += shift[k];
-    }
-
-    /* Loop over particles and find which particles belong in the same group. */
-    for (size_t j = 0; j < count_j; j++) {
-
-      struct gpart *pj = &gparts_j[j];
-      const double pjx = pj->x[0] - shift[0];
-      const double pjy = pj->x[1] - shift[1];
-      const double pjz = pj->x[2] - shift[2];
-
-      /* Find the root of pj. */
-      int root_j = fof_find(offset_j[j] - node_offset, group_index);
-
-      for (size_t i = 0; i < count_i; i++) {
-
-        struct gpart *pi = &gparts_i[i];
-        const double pix = pi->x[0];
-        const double piy = pi->x[1];
-        const double piz = pi->x[2];
-
-        /* Compute pairwise distance, remembering to account for boundary
-         * conditions. */
-        float dx[3], r2 = 0.0f;
-        dx[0] = pjx - pix;
-        dx[1] = pjy - piy;
-        dx[2] = pjz - piz;
-
-        for (int k = 0; k < 3; k++) r2 += dx[k] * dx[k];
-
-        /* Hit or miss? */
-        if (r2 < l_x2) {
-
-          //message("Particle %lld in range of foreign particle %lld. Group ID: %d root_j: %d", pj->id_or_neg_offset, pi->id_or_neg_offset, group_index[root_j - node_offset], root_j);
-          /* Store the particle IDs and group ID for communication. */
-          fof_send[*send_count].local_pid = pj->id_or_neg_offset;        
-          fof_send[*send_count].foreign_pid = pi->id_or_neg_offset; 
-          fof_send[*send_count].root_i = root_j;
-          (*send_count)++;
-          local_send_count++;
-          //fprintf(fof_file, "  %7zu %7zu %7d\n", pj->id_or_neg_offset, pi->id_or_neg_offset, root_j);
-        }
-      }
-    }
-
-  }
-
-  //if(local_send_count > 0) {
-
-  //  if(ci_local) 
-  //    message("Rank: %d sending %zu links to rank %d for testing. ci: (%lf,%lf,%lf) -> cj: (%lf,%lf,%lf).", engine_rank, local_send_count, cj->nodeID, ci->loc[0], ci->loc[1], ci->loc[2], cj->loc[0], cj->loc[1], cj->loc[2]);
-  //  
-  //  if(cj_local) 
-  //    message("Rank: %d sending %zu links to rank %d for testing. cj: (%lf,%lf,%lf) -> ci: (%lf,%lf,%lf).", engine_rank, local_send_count, ci->nodeID, cj->loc[0], cj->loc[1], cj->loc[2], ci->loc[0], ci->loc[1], ci->loc[2]);
-  //}
-
-}
-
-/* Perform a FOF search between a local and foreign cell using the Union-Find algorithm. 
- * Communicate any links found between particles to the appropriate node.*/
-void fof_search_pair_cells_foreign_new(struct space *s, struct cell *ci, struct cell *cj, size_t *send_count, struct fof_mpi *fof_send ) {
-
-  const size_t count_i = ci->gcount;
-  const size_t count_j = cj->gcount;
-  struct gpart *gparts_i = ci->gparts;
-  struct gpart *gparts_j = cj->gparts;
-  const double dim[3] = {s->dim[0], s->dim[1], s->dim[2]};
-  const double l_x2 = s->l_x2;
-  int *group_index = s->fof_data.group_index;
-  
-  /* Make a list of particle offsets into the global gparts array. */
-  int *const offset_i = group_index + (ptrdiff_t)(gparts_i - s->gparts);
-  int *const offset_j = group_index + (ptrdiff_t)(gparts_j - s->gparts);
-
-  /* Account for boundary conditions.*/
-  double shift[3] = {0.0, 0.0, 0.0};
-
-  /* Check whether cells are local to the node. */
-  const int ci_local = (ci->nodeID == engine_rank);
-  const int cj_local = (cj->nodeID == engine_rank);
-
-  if((ci_local && cj_local) || (!ci_local && !cj_local)) error("FOF search of foreign cells called on two local cells or two foreign cells.");
-
-  size_t local_send_count = 0;
-
-  /* Get the relative distance between the pairs, wrapping. */
-  const int periodic = s->periodic;
-  double diff[3];
-  
   if(ci_local) {
 
     for (int k = 0; k < 3; k++) {
@@ -984,6 +827,8 @@ void fof_search_tree_mapper(void *map_data, int num_elements,
  */
 void fof_search_foreign_cells(struct space *s) {
 
+#ifdef WITH_MPI
+
   struct engine *e = s->e;
   struct cell *cells = s->cells_top;
   int *group_index = s->fof_data.group_index;
@@ -1067,8 +912,6 @@ void fof_search_foreign_cells(struct space *s) {
         /* Skip empty cells. */
         if(foreign_cell->gcount == 0) continue;
         
-        //rec_fof_search_pair_foreign(local_cell, foreign_cell, s, dim, search_r2, &send_count, fof_send);
-
         /* Check if local cell has already been added to the local list of cells. */
         if(!found) {
           const double r2 = cell_min_dist(local_cell, foreign_cell, dim);
@@ -1108,17 +951,21 @@ void fof_search_foreign_cells(struct space *s) {
 
   for(int i=0; i<sched->nr_tasks; i++) {
   
-    if(tasks[i].type == task_type_send) {
-      scheduler_activate(sched,&tasks[i]);
+    struct task *t = &tasks[i];
+
+    if(t->type == task_type_send && t->subtype == task_subtype_gpart) {
+      scheduler_activate(sched,t);
       nr_sends++;
     }
 
-    if(tasks[i].type == task_type_recv) {
-      scheduler_activate(sched,&tasks[i]);
+    if(t->type == task_type_recv && t->subtype == task_subtype_gpart) {
+      scheduler_activate(sched,t);
       nr_recvs++;
     } 
   }
- 
+
+  message("Rank %d. No. of sends activated: %d, no. of receives activated: %d", engine_rank, nr_sends, nr_recvs);
+
   engine_launch(e);
 
   for(int i=0; i<e->nr_proxies; i++) {
@@ -1127,14 +974,23 @@ void fof_search_foreign_cells(struct space *s) {
 
       struct cell *restrict foreign_cell = e->proxies[i].cells_in[j];
       struct gpart *gparts = foreign_cell->gparts;
-  
+     
+      if(foreign_cell->nodeID == engine_rank) error("Rank %d. Foreign cell is actually local.", engine_rank); 
+
+      if(gparts[0].root >= node_offset && gparts[0].root < node_offset + s->nr_gparts) {
+          //error("Rank %d received foreign particle with local root: %d", engine_rank, gparts[k].root);
+          message("Rank %d received foreign particle %lld from rank %d with a local root: %d. i=%d,j=%d.", engine_rank, gparts[0].id_or_neg_offset, foreign_cell->nodeID, gparts[0].root, i, j);
+      }
       for(int k=0; k<foreign_cell->gcount; k++) {
-        if(gparts[k].root >= node_offset && gparts[k].root < node_offset + s->nr_gparts) error("Rank %d received foreign particle with local root: %d", engine_rank, gparts[k].root);
+        if(gparts[k].root >= node_offset && gparts[k].root < node_offset + s->nr_gparts) {
+          //error("Rank %d received foreign particle with local root: %d", engine_rank, gparts[k].root);
+          //message("Rank %d received foreign particle %lld from rank %d with a local root: %d. i=%d,j=%d,k=%d.", engine_rank, gparts[k].id_or_neg_offset, foreign_cell->nodeID, gparts[k].root, i, j, k);
+        }
       }
     }
   }
 
-  size_t link_count_new = 0;
+  size_t link_count = 0;
 
   for(int i=0; i<interface_cell_count; i++) {
 
@@ -1152,27 +1008,29 @@ void fof_search_foreign_cells(struct space *s) {
         /* Skip empty cells. */
         if(foreign_cell->gcount == 0) continue;
 
-        rec_fof_search_pair_foreign(local_cell, foreign_cell, s, dim, search_r2, &link_count_new, fof_send);
+        rec_fof_search_pair_foreign(local_cell, foreign_cell, s, dim, search_r2, &link_count, fof_send);
 
       }
     }
   }
 
-  message("Rank %d found %zu links between local and foreign groups.", engine_rank, link_count_new);
+  message("Rank %d found %zu links between local and foreign groups.", engine_rank, link_count);
 
   struct fof_mpi *fof_list;
   int list_count = 0;
   
   if (posix_memalign((void**)&fof_list, SWIFT_STRUCT_ALIGNMENT,
-                       link_count_new * sizeof(struct fof_mpi)) != 0)
+                       link_count * sizeof(struct fof_mpi)) != 0)
     error("Error while allocating memory for FOF MPI communication");
 
-  for(int i=0; i<link_count_new; i++) {
+  for(int i=0; i<link_count; i++) {
     
     int found = 0;
     int local_root = fof_send[i].local_pid;
     int foreign_root = fof_send[i].foreign_pid;
 
+    if(local_root == foreign_root) error("Particles have same root. local_root: %d, foreign_root: %d", local_root, foreign_root);
+
     for(int j=0; j<list_count; j++) {
       if(fof_list[j].local_pid == local_root && fof_list[j].foreign_pid == foreign_root) {
         found = 1;
@@ -1193,6 +1051,13 @@ void fof_search_foreign_cells(struct space *s) {
   int global_list_count = 0;
 
   MPI_Allreduce(&list_count, &global_list_count, 1, MPI_INT, MPI_SUM, MPI_COMM_WORLD);
+    
+  int *list_counts = NULL;
+  if (posix_memalign((void**)&list_counts, SWIFT_STRUCT_ALIGNMENT,
+                       e->nr_nodes * sizeof(int)) != 0)
+    error("Error while allocating memory for FOF MPI communication");
+
+  MPI_Allgather(&list_count, 1, MPI_INT, list_counts, 1, MPI_INT, MPI_COMM_WORLD);
  
   message("Global list count: %d", global_list_count);
 
@@ -1204,13 +1069,23 @@ void fof_search_foreign_cells(struct space *s) {
                        global_list_count * sizeof(struct fof_mpi)) != 0)
     error("Error while allocating memory for FOF MPI communication");
  
-  MPI_Allgather(fof_list, list_count, fof_mpi_type, global_fof_list, list_count, fof_mpi_type, MPI_COMM_WORLD);
+  int *displ = NULL;
+  if (posix_memalign((void**)&displ, SWIFT_STRUCT_ALIGNMENT,
+                       e->nr_nodes * sizeof(int)) != 0)
+    error("Error while allocating memory for FOF MPI communication");
+
+  displ[0] = 0;
+  for(int i=1; i<e->nr_nodes; i++) displ[i] = displ[i-1] + list_counts[i-1];
+
+  MPI_Allgatherv(fof_list, list_count, fof_mpi_type, global_fof_list, list_counts, displ, fof_mpi_type, MPI_COMM_WORLD);
 
   for(int i=0; i<global_list_count; i++) {
     fprintf(fof_file, "  %7lld <-> %7lld\n", global_fof_list[i].local_pid, global_fof_list[i].foreign_pid);
     //message("Rank %d. Global FOF list: %lld <-> %lld", engine_rank, global_fof_list[i].local_pid, global_fof_list[i].foreign_pid);
   }
   
+  fclose(fof_file);
+  
   for(int i=0; i<global_list_count; i++) {
     global_fof_list_offset[i].local_pid = i;
     global_fof_list_offset[i].foreign_pid = i;
@@ -1219,50 +1094,49 @@ void fof_search_foreign_cells(struct space *s) {
   for(int i=0; i<global_list_count; i++) {
     int group_i = global_fof_list[i].local_pid;
     int group_j = global_fof_list[i].foreign_pid;
+    
+    //if((group_j >= node_offset && group_j < node_offset + s->nr_gparts)) {
+    //  int root = fof_find(group_j - node_offset, group_index);
 
-    int min_root = min(group_i, group_j);
+    //  group_j = root;
 
-    //if(group_j < group_i) {
-      
-      //min_root = group_j;
+    //}
+    //
+    //if((group_i >= node_offset && group_i < node_offset + s->nr_gparts)) {
+    //  int root = fof_find(group_i - node_offset, group_index);
 
-      for(int j=i + 1; j<global_list_count; j++) {
-        int root_i = global_fof_list[j].local_pid;
-        int root_j = global_fof_list[j].foreign_pid;
+    //  group_i = root;
 
-        if(root_i == group_j) {
-          if(root_j < min_root) min_root = root_j;
-        }
-        else if(root_j == group_j) {
-          if(root_i < min_root) min_root = root_i;
-        }
-      }
+    //}
 
-      //if( group_i >= node_offset && group_i < node_offset + s->nr_gparts) {
-      //  message("Rank %d. Group %d links to %d.", engine_rank,group_i, min_root);
-      //}
+    if(group_i == group_j) error("Particles have same root. local_root: %d, foreign_root: %d", group_i, group_j);
 
-    //}
-    //else {
-      //min_root = group_i;
-      
-      for(int j=i + 1; j<global_list_count; j++) {
-        int root_i = global_fof_list[j].local_pid;
-        int root_j = global_fof_list[j].foreign_pid;
+    int min_root = min(group_i, group_j);
 
-        if(root_i == group_i) {
-          if(root_j < min_root) min_root = root_j;
-        }
-        else if(root_j == group_i) {
-          if(root_i < min_root) min_root = root_i;
-        }
+    for(int j=i + 1; j<global_list_count; j++) {
+      int root_i = global_fof_list[j].local_pid;
+      int root_j = global_fof_list[j].foreign_pid;
+
+      if(root_i == group_j) {
+        if(root_j < min_root) min_root = root_j;
+      }
+      else if(root_j == group_j) {
+        if(root_i < min_root) min_root = root_i;
+      }
+    }
+
+    for(int j=i + 1; j<global_list_count; j++) {
+      int root_i = global_fof_list[j].local_pid;
+      int root_j = global_fof_list[j].foreign_pid;
+
+      if(root_i == group_i) {
+        if(root_j < min_root) min_root = root_j;
+      }
+      else if(root_j == group_i) {
+        if(root_i < min_root) min_root = root_i;
       }
+    }
 
-      //if( group_j >= node_offset && group_j < node_offset + s->nr_gparts) {
-      //  message("Rank %d. Group %d links to %d.", engine_rank, group_j, min_root);
-      //}
-    //}
-   
     if(((group_j >= node_offset && group_j < node_offset + s->nr_gparts) ||
        (group_i >= node_offset && group_i < node_offset + s->nr_gparts) ) &&
        (group_i != min_root && group_j != min_root) ) {
@@ -1272,21 +1146,60 @@ void fof_search_foreign_cells(struct space *s) {
     if((group_j >= node_offset && group_j < node_offset + s->nr_gparts) &&
        (group_j != min_root)) {
       int root = fof_find(group_j - node_offset, group_index);
+      if(root != group_j) message("Rank %d. The root of group_j has changed from: %d to %d", engine_rank, group_j, root);
+      //message("Rank %d. Updating local root: %d to: %d", engine_rank, root, min_root);
+
+      if(root < min_root) message("Rank %d. Group_j: %d Trying to set a root: %d to a larger value: %d", engine_rank, group_j, root, min_root);
+     
+      /* Check root is local to node. */ 
+      if(root >= node_offset && root < node_offset + s->nr_gparts) {
+        if(root > min_root) {
+          group_index[root - node_offset] = min_root;
+        }
+        else if(min_root >= node_offset && min_root < node_offset + s->nr_gparts){
+
+          int root_new = fof_find(min_root - node_offset, group_index);
+          //if(min_root != root_new) error("Wrong min_root used. min_root in fof_list: %d and actual root found in group_index: %d", min_root, root_new);
+
+          group_index[root_new - node_offset] = root;
 
-      group_index[root - node_offset] = min_root;
+        }
+      }
     }
     
     if((group_i >= node_offset && group_i < node_offset + s->nr_gparts) &&
        (group_i != min_root)) {
       int root = fof_find(group_i - node_offset, group_index);
+      if(root != group_i) message("Rank %d. The root of group_i has changed from: %d to %d", engine_rank, group_i, root);
+
+      //message("Rank %d. Updating local root: %d to: %d", engine_rank, root, min_root);
+      
+      if(root < min_root) message("Rank %d. Group_i: %d Trying to set a root: %d to a larger value: %d", engine_rank, group_i, root, min_root);
+      
+      /* Check root is local to node. */ 
+      if(root >= node_offset && root < node_offset + s->nr_gparts) {
+        if(root > min_root) {
+          group_index[root - node_offset] = min_root;
+        }
+        else if(min_root >= node_offset && min_root < node_offset + s->nr_gparts){
+
+          int root_new = fof_find(min_root - node_offset, group_index);
+          //if(min_root != root_new) error("Wrong min_root used. min_root in fof_list: %d and actual root found in group_index: %d", min_root, root_new);
 
-      group_index[root - node_offset] = min_root;
+          group_index[root_new - node_offset] = root;
+
+        }
+      }
     }
+
+    //if(engine_rank == 0) message("Rank 0. Updating local roots. %d of %d....", i, global_list_count);
     //global_fof_list[i].local_pid = min_root;
     //global_fof_list[i].foreign_pid = min_root;
   }
 
-  fclose(fof_file);
+  message("Rank %d finished linking local roots to foreign roots.", engine_rank);
+
+#endif /* WITH_MPI */
 
 }
 
@@ -1321,9 +1234,11 @@ void fof_search_tree(struct space *s) {
 
     bzero(global_nr_gparts, s->e->nr_nodes * sizeof(int));
 
-    MPI_Allgather(&s->nr_gparts, 1, MPI_INT, global_nr_gparts, s->e->nr_nodes - 1, MPI_INT, MPI_COMM_WORLD);
+    MPI_Allgather(&s->nr_gparts, 1, MPI_INT, global_nr_gparts, 1, MPI_INT, MPI_COMM_WORLD);
 
-    message("No. of particles: [%d, %d]", global_nr_gparts[0], global_nr_gparts[1]);
+    printf("No. of particles: [%d", global_nr_gparts[0]);
+    for(int i = 1; i<s->e->nr_nodes; i++) printf(", %d", global_nr_gparts[i]);
+    printf("]\n");
 
     for(int i = 0; i<engine_rank; i++) node_offset += global_nr_gparts[i];
   }
@@ -1418,7 +1333,13 @@ void fof_search_tree(struct space *s) {
     bzero(global_group_size, s->e->total_nr_gparts * sizeof(int));
     bzero(global_group_id, s->e->total_nr_gparts * sizeof(long long));
 
-    int displ[2] = {0, nr_gparts};
+    int *displ = NULL;
+    if (posix_memalign((void**)&displ, SWIFT_STRUCT_ALIGNMENT,
+          s->e->nr_nodes * sizeof(int)) != 0)
+      error("Error while allocating memory for FOF MPI communication");
+
+    displ[0] = 0;
+    for(int i=1; i<s->e->nr_nodes; i++) displ[i] = displ[i-1] + global_nr_gparts[i-1];
 
     MPI_Gatherv(group_index, nr_gparts, MPI_INT, global_group_index, global_nr_gparts, displ, MPI_INT, 0, MPI_COMM_WORLD);
     MPI_Gatherv(group_id, nr_gparts, MPI_LONG_LONG, global_group_id, global_nr_gparts, displ, MPI_LONG_LONG, 0, MPI_COMM_WORLD);
diff --git a/src/fof.h b/src/fof.h
index 2d748bfd6cb17d1cf23a723b0d4c5c694caa33ba..50d1edaf0ca171d999df4a1608a4fd719e82e2aa 100644
--- a/src/fof.h
+++ b/src/fof.h
@@ -47,7 +47,6 @@ void fof_search_serial(struct space *s);
 void fof_search_cell(struct space *s, struct cell *c);
 void fof_search_pair_cells(struct space *s, struct cell *ci, struct cell *cj);
 void fof_search_pair_cells_foreign(struct space *s, struct cell *ci, struct cell *cj, size_t *send_count, struct fof_mpi *fof_send);
-void fof_search_pair_cells_foreign_new(struct space *s, struct cell *ci, struct cell *cj, size_t *send_count, struct fof_mpi *fof_send);
 void fof_search_tree_serial(struct space *s);
 void fof_search_tree(struct space *s);
 void fof_dump_group_data(char *out_file, const size_t nr_gparts, int *group_index, int *num_in_groups, long long *group_id);