diff --git a/src/partition.c b/src/partition.c
index 0f8eb3ebe334d71228510307dd9ccc4e56e234b3..7dbbb9552e603adee45097a379200f1493ce3349 100644
--- a/src/partition.c
+++ b/src/partition.c
@@ -927,7 +927,10 @@ static int check_complete(struct space *s, int verbose, int nregions) {
int failed = 0;
for (int i = 0; i < nregions; i++) present[i] = 0;
for (int i = 0; i < s->nr_cells; i++) {
- present[s->cells[i].nodeID]++;
+ if (s->cells[i].nodeID <= nregions)
+ present[s->cells[i].nodeID]++;
+ else
+ message("Bad nodeID: %d", s->cells[i].nodeID);
}
for (int i = 0; i < nregions; i++) {
if (!present[i]) {
@@ -938,3 +941,53 @@ static int check_complete(struct space *s, int verbose, int nregions) {
free(present);
return (!failed);
}
+
+
+/**
+ * @brief Partition a space of cells based on another space of cells.
+ *
+ * The two spaces are expected to be at different cell sizes, so what we'd
+ * like to do is assign the second space to geometrically closest nodes
+ * of the first, with the effect of minimizing particle movement when
+ * rebuilding the second space from the first.
+ *
+ * Since two spaces cannot exist simultaneously the old space is actually
+ * required in a decomposed state. These are the old cells sizes and counts
+ * per dimension, along with a list of the old nodeIDs. The old nodeIDs are
+ * indexed by the cellid (see cell_getid()), so should be stored that way.
+ *
+ * On exit the new space cells will have their nodeIDs assigned.
+ *
+ * @param oldh the cell dimensions of old space.
+ * @param oldcdim number of cells per dimension in old space.
+ * @param oldnodeIDs the nodeIDs of cells in the old space, indexed by old cellid.
+ * @param s the space to be partitioned.
+ *
+ * @return 1 if the new space contains nodeIDs from all nodes, 0 otherwise.
+ */
+int partition_space_to_space(double *oldh, double *oldcdim, int *oldnodeIDs,
+ struct space *s) {
+
+ /* Loop over all the new cells. */
+ int nr_nodes = 0;
+ for (int i = 0; i < s->cdim[0]; i++) {
+ for (int j = 0; j < s->cdim[1]; j++) {
+ for (int k = 0; k < s->cdim[2]; k++) {
+
+ /* Scale indices to old cell space. */
+ int ii = rint(i * s->ih[0] * oldh[0]);
+ int jj = rint(j * s->ih[1] * oldh[1]);
+ int kk = rint(k * s->ih[2] * oldh[2]);
+
+ int cid = cell_getid(s->cdim, i, j, k);
+ int oldcid = cell_getid(oldcdim, ii, jj, kk);
+ s->cells[cid].nodeID = oldnodeIDs[oldcid];
+
+ if (oldnodeIDs[oldcid] > nr_nodes) nr_nodes = oldnodeIDs[oldcid];
+ }
+ }
+ }
+
+ /* Check we have all nodeIDs present in the resample. */
+ return check_complete(s, 1, nr_nodes + 1);
+}
diff --git a/src/partition.h b/src/partition.h
index 3ab5f5a817bf5b77d45c7fb3313158b83d98e251..08906c765bbcf71b084829502e632e3159ebd0bf 100644
--- a/src/partition.h
+++ b/src/partition.h
@@ -56,4 +56,6 @@ void partition_repartition(enum repartition_type reparttype, int nodeID,
void partition_initial_partition(struct partition *initial_partition,
int nodeID, int nr_nodes, struct space *s);
+int partition_space_to_space(double *oldh, double *oldcdim, int *oldnodeID,
+ struct space *s);
#endif /* SWIFT_PARTITION_H */
diff --git a/src/space.c b/src/space.c
index 0785740335f6e782a875b7ab658a9a0564288aab..c25244b0c2c1d96b7ba2fd009586f54ec3f6f090 100644
--- a/src/space.c
+++ b/src/space.c
@@ -1,6 +1,7 @@
/*******************************************************************************
* This file is part of SWIFT.
* Copyright (c) 2012 Pedro Gonnet (pedro.gonnet@durham.ac.uk)
+* 2016 Peter W. Draper (p.w.draper@durham.ac.uk)
*
* This program is free software: you can redistribute it and/or modify
* it under the terms of the GNU Lesser General Public License as published
@@ -205,9 +206,33 @@ void space_regrid(struct space *s, double cell_max, int verbose) {
* global partition is recomputed and the particles redistributed.
* Be prepared to do that. */
#ifdef WITH_MPI
- int partition = 0;
- if (cdim[0] < s->cdim[0] || cdim[1] < s->cdim[1] || cdim[2] < s->cdim[2])
- partition = 1;
+ double oldh[3];
+ double oldcdim[3];
+ int *oldnodeIDs = NULL;
+ if (cdim[0] < s->cdim[0] || cdim[1] < s->cdim[1] || cdim[2] < s->cdim[2]) {
+
+ /* Capture state of current space. */
+ oldcdim[0] = s->cdim[0];
+ oldcdim[1] = s->cdim[1];
+ oldcdim[2] = s->cdim[2];
+ oldh[0] = s->h[0];
+ oldh[1] = s->h[1];
+ oldh[2] = s->h[2];
+
+ if ((oldnodeIDs = (int *)malloc(sizeof(int) * s->nr_cells)) == NULL)
+ error("Failed to allocate temporary nodeIDs.");
+
+ int cid = 0;
+ for (int i = 0; i < s->cdim[0]; i++) {
+ for (int j = 0; j < s->cdim[1]; j++) {
+ for (int k = 0; k < s->cdim[2]; k++) {
+ cid = cell_getid(oldcdim, i, j, k);
+ oldnodeIDs[cid] = s->cells[cid].nodeID;
+ }
+ }
+ }
+ }
+
#endif
/* Do we need to re-build the upper-level cells? */
@@ -268,29 +293,36 @@ void space_regrid(struct space *s, double cell_max, int verbose) {
fflush(stdout);
#ifdef WITH_MPI
- /* XXX create an engine_resplit() function to use the same method as the
- * initial partition. Fake for now. */
- if (partition) {
- if (s->e->nodeID == 0)
- message("cell dimensions have decreased. Recalculating the "
- "global partition.");
-
- /* Change the global partitioning. */
- int grid[3];
- factor(s->e->nr_nodes, &grid[0], &grid[1]);
- factor(s->e->nr_nodes / grid[1], &grid[0], &grid[2]);
- factor(grid[0] * grid[1], &grid[1], &grid[0]);
-
- /* Run through the cells and set their nodeID. */
- int ind[3];
- for (int k = 0; k < s->nr_cells; k++) {
- c = &s->cells[k];
- for (int j = 0; j < 3; j++) ind[j] = c->loc[j] / s->dim[j] * grid[j];
- c->nodeID = ind[0] + grid[0] * (ind[1] + grid[1] * ind[2]);
- }
+ if (oldnodeIDs != NULL) {
+ /* We have changed the top-level cell dimension, so need to redistribute
+ * cells around the nodes. We repartition using the old space node
+ * positions as a grid to resample. */
+ if (s->e->nodeID == 0)
+ message("basic cell dimensions have increased - recalculating the "
+ "global partition.");
+
+ if (!partition_space_to_space(oldh, oldcdim, oldnodeIDs, s) ) {
+
+ /* Failed, try another technique that requires no settings. */
+ message("Failed to get a new partition, trying less optimal method");
+ struct partition initial_partition;
+#ifdef HAVE_METIS
+ initial_partition.type = INITPART_METIS_NOWEIGHT;
+#else
+ initial_partition.type = INITPART_VECTORIZE;
+#endif
+ partition_initial_partition(&initial_partition, s->e->nodeID,
+ s->e->nr_nodes, s);
+ }
+
+ /* Re-distribute the particles to their new nodes. */
+ engine_redistribute(s->e);
+
+ /* Make the proxies. */
+ engine_makeproxies(s->e);
- /* Make the proxies. */
- engine_makeproxies(s->e);
+ /* Finished with these. */
+ free(oldnodeIDs);
}
#endif
} /* re-build upper-level cells? */