Commit bb9bde8e authored by Matthieu Schaller's avatar Matthieu Schaller
Browse files

Star particles can now be used over MPI.

parent db251a47
......@@ -4,7 +4,7 @@
if [ ! -e multiTypes.hdf5 ]
then
echo "Generating initial conditions for the multitype box example..."
python makeIC.py 50 60
python makeIC.py 17 24 12
fi
../swift -s -g -t 16 multiTypes.yml 2>&1 | tee output.log
../swift -s -g -S -t 1 multiTypes.yml 2>&1 | tee output.log
......@@ -196,6 +196,31 @@ int cell_link_gparts(struct cell *c, struct gpart *gparts) {
return c->gcount;
}
/**
* @brief Link the cells recursively to the given #spart array.
*
* @param c The #cell.
* @param sparts The #spart array.
*
* @return The number of particles linked.
*/
int cell_link_sparts(struct cell *c, struct spart *sparts) {
c->sparts = sparts;
/* Fill the progeny recursively, depth-first. */
if (c->split) {
int offset = 0;
for (int k = 0; k < 8; k++) {
if (c->progeny[k] != NULL)
offset += cell_link_sparts(c->progeny[k], &sparts[offset]);
}
}
/* Return the total number of linked particles. */
return c->scount;
}
/**
* @brief Pack the data of the given cell and all it's sub-cells.
*
......@@ -216,6 +241,7 @@ int cell_pack(struct cell *c, struct pcell *pc) {
pc->ti_old = c->ti_old;
pc->count = c->count;
pc->gcount = c->gcount;
pc->scount = c->scount;
c->tag = pc->tag = atomic_inc(&cell_next_tag) % cell_max_tag;
/* Fill in the progeny, depth-first recursion. */
......
......@@ -320,6 +320,7 @@ int cell_unpack_ti_ends(struct cell *c, integertime_t *ti_ends);
int cell_getsize(struct cell *c);
int cell_link_parts(struct cell *c, struct part *parts);
int cell_link_gparts(struct cell *c, struct gpart *gparts);
int cell_link_sparts(struct cell *c, struct spart *sparts);
void cell_convert_hydro(struct cell *c, void *data);
void cell_clean_links(struct cell *c, void *data);
int cell_are_neighbours(const struct cell *restrict ci,
......
This diff is collapsed.
......@@ -234,7 +234,9 @@ void engine_maketasks(struct engine *e);
void engine_split(struct engine *e, struct partition *initial_partition);
void engine_exchange_strays(struct engine *e, size_t offset_parts,
int *ind_part, size_t *Npart, size_t offset_gparts,
int *ind_gpart, size_t *Ngpart);
int *ind_gpart, size_t *Ngpart,
size_t offset_sparts, int *ind_spart,
size_t *Nspart);
void engine_rebuild(struct engine *e);
void engine_repartition(struct engine *e);
void engine_makeproxies(struct engine *e);
......
......@@ -204,9 +204,8 @@ void part_verify_links(struct part *parts, struct gpart *gparts,
MPI_Datatype part_mpi_type;
MPI_Datatype xpart_mpi_type;
MPI_Datatype gpart_mpi_type;
#endif
MPI_Datatype spart_mpi_type;
#ifdef WITH_MPI
/**
* @brief Registers MPI particle types.
*/
......@@ -233,5 +232,10 @@ void part_create_mpi_types() {
MPI_Type_commit(&gpart_mpi_type) != MPI_SUCCESS) {
error("Failed to create MPI type for gparts.");
}
if (MPI_Type_contiguous(sizeof(struct spart) / sizeof(unsigned char),
MPI_BYTE, &spart_mpi_type) != MPI_SUCCESS ||
MPI_Type_commit(&spart_mpi_type) != MPI_SUCCESS) {
error("Failed to create MPI type for sparts.");
}
}
#endif
......@@ -85,6 +85,7 @@ void part_verify_links(struct part *parts, struct gpart *gparts,
extern MPI_Datatype part_mpi_type;
extern MPI_Datatype xpart_mpi_type;
extern MPI_Datatype gpart_mpi_type;
extern MPI_Datatype spart_mpi_type;
void part_create_mpi_types();
#endif
......
......@@ -45,7 +45,6 @@
*
* @param p The #proxy.
*/
void proxy_cells_exch1(struct proxy *p) {
#ifdef WITH_MPI
......@@ -126,7 +125,6 @@ void proxy_cells_exch2(struct proxy *p) {
* @param p The #proxy.
* @param c The #cell.
*/
void proxy_addcell_in(struct proxy *p, struct cell *c) {
/* Check if the cell is already registered with the proxy. */
......@@ -155,7 +153,6 @@ void proxy_addcell_in(struct proxy *p, struct cell *c) {
* @param p The #proxy.
* @param c The #cell.
*/
void proxy_addcell_out(struct proxy *p, struct cell *c) {
/* Check if the cell is already registered with the proxy. */
......@@ -183,7 +180,6 @@ void proxy_addcell_out(struct proxy *p, struct cell *c) {
*
* @param p The #proxy.
*/
void proxy_parts_exch1(struct proxy *p) {
#ifdef WITH_MPI
......@@ -191,7 +187,8 @@ void proxy_parts_exch1(struct proxy *p) {
/* Send the number of particles. */
p->buff_out[0] = p->nr_parts_out;
p->buff_out[1] = p->nr_gparts_out;
if (MPI_Isend(p->buff_out, 2, MPI_INT, p->nodeID,
p->buff_out[2] = p->nr_sparts_out;
if (MPI_Isend(p->buff_out, 3, MPI_INT, p->nodeID,
p->mynodeID * proxy_tag_shift + proxy_tag_count, MPI_COMM_WORLD,
&p->req_parts_count_out) != MPI_SUCCESS)
error("Failed to isend nr of parts.");
......@@ -218,13 +215,22 @@ void proxy_parts_exch1(struct proxy *p) {
if (MPI_Isend(p->gparts_out, p->nr_gparts_out, gpart_mpi_type, p->nodeID,
p->mynodeID * proxy_tag_shift + proxy_tag_gparts,
MPI_COMM_WORLD, &p->req_gparts_out) != MPI_SUCCESS)
error("Failed to isend part data.");
error("Failed to isend gpart data.");
// message( "isent gpart data (%i) to node %i." , p->nr_parts_out ,
// p->nodeID ); fflush(stdout);
}
if (p->nr_sparts_out > 0) {
if (MPI_Isend(p->sparts_out, p->nr_sparts_out, spart_mpi_type, p->nodeID,
p->mynodeID * proxy_tag_shift + proxy_tag_sparts,
MPI_COMM_WORLD, &p->req_sparts_out) != MPI_SUCCESS)
error("Failed to isend spart data.");
// message( "isent gpart data (%i) to node %i." , p->nr_parts_out ,
// p->nodeID ); fflush(stdout);
}
/* Receive the number of particles. */
if (MPI_Irecv(p->buff_in, 2, MPI_INT, p->nodeID,
if (MPI_Irecv(p->buff_in, 3, MPI_INT, p->nodeID,
p->nodeID * proxy_tag_shift + proxy_tag_count, MPI_COMM_WORLD,
&p->req_parts_count_in) != MPI_SUCCESS)
error("Failed to irecv nr of parts.");
......@@ -241,8 +247,9 @@ void proxy_parts_exch2(struct proxy *p) {
/* Unpack the incomming parts counts. */
p->nr_parts_in = p->buff_in[0];
p->nr_gparts_in = p->buff_in[1];
p->nr_sparts_in = p->buff_in[2];
/* Is there enough space in the buffer? */
/* Is there enough space in the buffers? */
if (p->nr_parts_in > p->size_parts_in) {
do {
p->size_parts_in *= proxy_buffgrow;
......@@ -264,6 +271,15 @@ void proxy_parts_exch2(struct proxy *p) {
p->size_gparts_in)) == NULL)
error("Failed to re-allocate gparts_in buffers.");
}
if (p->nr_sparts_in > p->size_sparts_in) {
do {
p->size_sparts_in *= proxy_buffgrow;
} while (p->nr_sparts_in > p->size_sparts_in);
free(p->sparts_in);
if ((p->sparts_in = (struct spart *)malloc(sizeof(struct spart) *
p->size_sparts_in)) == NULL)
error("Failed to re-allocate sparts_in buffers.");
}
/* Receive the particle buffers. */
if (p->nr_parts_in > 0) {
......@@ -285,6 +301,14 @@ void proxy_parts_exch2(struct proxy *p) {
// message( "irecv gpart data (%i) from node %i." , p->nr_gparts_in ,
// p->nodeID ); fflush(stdout);
}
if (p->nr_sparts_in > 0) {
if (MPI_Irecv(p->sparts_in, p->nr_sparts_in, spart_mpi_type, p->nodeID,
p->nodeID * proxy_tag_shift + proxy_tag_sparts,
MPI_COMM_WORLD, &p->req_sparts_in) != MPI_SUCCESS)
error("Failed to irecv spart data.");
// message( "irecv gpart data (%i) from node %i." , p->nr_gparts_in ,
// p->nodeID ); fflush(stdout);
}
#else
error("SWIFT was not compiled with MPI support.");
......@@ -299,7 +323,6 @@ void proxy_parts_exch2(struct proxy *p) {
* @param xparts Pointer to an array of #xpart to send.
* @param N The number of parts.
*/
void proxy_parts_load(struct proxy *p, const struct part *parts,
const struct xpart *xparts, int N) {
......@@ -332,13 +355,12 @@ void proxy_parts_load(struct proxy *p, const struct part *parts,
}
/**
* @brief Load parts onto a proxy for exchange.
* @brief Load gparts onto a proxy for exchange.
*
* @param p The #proxy.
* @param gparts Pointer to an array of #gpart to send.
* @param N The number of parts.
* @param N The number of gparts.
*/
void proxy_gparts_load(struct proxy *p, const struct gpart *gparts, int N) {
/* Is there enough space in the buffer? */
......@@ -362,6 +384,36 @@ void proxy_gparts_load(struct proxy *p, const struct gpart *gparts, int N) {
p->nr_gparts_out += N;
}
/**
* @brief Load sparts onto a proxy for exchange.
*
* @param p The #proxy.
* @param sparts Pointer to an array of #spart to send.
* @param N The number of sparts.
*/
void proxy_sparts_load(struct proxy *p, const struct spart *sparts, int N) {
/* Is there enough space in the buffer? */
if (p->nr_sparts_out + N > p->size_sparts_out) {
do {
p->size_sparts_out *= proxy_buffgrow;
} while (p->nr_sparts_out + N > p->size_sparts_out);
struct spart *tp;
if ((tp = (struct spart *)malloc(sizeof(struct spart) *
p->size_sparts_out)) == NULL)
error("Failed to re-allocate sparts_out buffers.");
memcpy(tp, p->sparts_out, sizeof(struct spart) * p->nr_sparts_out);
free(p->sparts_out);
p->sparts_out = tp;
}
/* Copy the parts and xparts data to the buffer. */
memcpy(&p->sparts_out[p->nr_sparts_out], sparts, sizeof(struct spart) * N);
/* Increase the counters. */
p->nr_sparts_out += N;
}
/**
* @brief Initialize the given proxy.
*
......@@ -369,7 +421,6 @@ void proxy_gparts_load(struct proxy *p, const struct gpart *gparts, int N) {
* @param mynodeID The node this proxy is running on.
* @param nodeID The node with which this proxy will communicate.
*/
void proxy_init(struct proxy *p, int mynodeID, int nodeID) {
/* Set the nodeID. */
......@@ -427,4 +478,20 @@ void proxy_init(struct proxy *p, int mynodeID, int nodeID) {
error("Failed to allocate gparts_out buffers.");
}
p->nr_gparts_out = 0;
/* Allocate the spart send and receive buffers, if needed. */
if (p->sparts_in == NULL) {
p->size_sparts_in = proxy_buffinit;
if ((p->sparts_in = (struct spart *)malloc(sizeof(struct spart) *
p->size_sparts_in)) == NULL)
error("Failed to allocate sparts_in buffers.");
}
p->nr_sparts_in = 0;
if (p->sparts_out == NULL) {
p->size_sparts_out = proxy_buffinit;
if ((p->sparts_out = (struct spart *)malloc(sizeof(struct spart) *
p->size_sparts_out)) == NULL)
error("Failed to allocate sparts_out buffers.");
}
p->nr_sparts_out = 0;
}
......@@ -33,7 +33,8 @@
#define proxy_tag_parts 1
#define proxy_tag_xparts 2
#define proxy_tag_gparts 3
#define proxy_tag_cells 4
#define proxy_tag_sparts 4
#define proxy_tag_cells 5
/* Data structure for the proxy. */
struct proxy {
......@@ -55,13 +56,16 @@ struct proxy {
struct part *parts_in, *parts_out;
struct xpart *xparts_in, *xparts_out;
struct gpart *gparts_in, *gparts_out;
struct spart *sparts_in, *sparts_out;
int size_parts_in, size_parts_out;
int nr_parts_in, nr_parts_out;
int size_gparts_in, size_gparts_out;
int nr_gparts_in, nr_gparts_out;
int size_sparts_in, size_sparts_out;
int nr_sparts_in, nr_sparts_out;
/* Buffer to hold the incomming/outgoing particle counts. */
int buff_out[2], buff_in[2];
int buff_out[3], buff_in[3];
/* MPI request handles. */
#ifdef WITH_MPI
......@@ -69,6 +73,7 @@ struct proxy {
MPI_Request req_parts_out, req_parts_in;
MPI_Request req_xparts_out, req_xparts_in;
MPI_Request req_gparts_out, req_gparts_in;
MPI_Request req_sparts_out, req_sparts_in;
MPI_Request req_cells_count_out, req_cells_count_in;
MPI_Request req_cells_out, req_cells_in;
#endif
......@@ -79,6 +84,7 @@ void proxy_init(struct proxy *p, int mynodeID, int nodeID);
void proxy_parts_load(struct proxy *p, const struct part *parts,
const struct xpart *xparts, int N);
void proxy_gparts_load(struct proxy *p, const struct gpart *gparts, int N);
void proxy_sparts_load(struct proxy *p, const struct spart *sparts, int N);
void proxy_parts_exch1(struct proxy *p);
void proxy_parts_exch2(struct proxy *p);
void proxy_addcell_in(struct proxy *p, struct cell *c);
......
......@@ -521,18 +521,22 @@ void space_rebuild(struct space *s, int verbose) {
for (size_t k = 0; k < nr_parts;) {
if (cells_top[ind[k]].nodeID != local_nodeID) {
nr_parts -= 1;
/* Swap the particle */
const struct part tp = s->parts[k];
s->parts[k] = s->parts[nr_parts];
s->parts[nr_parts] = tp;
/* Swap the link with the gpart */
if (s->parts[k].gpart != NULL) {
s->parts[k].gpart->id_or_neg_offset = -k;
}
if (s->parts[nr_parts].gpart != NULL) {
s->parts[nr_parts].gpart->id_or_neg_offset = -nr_parts;
}
/* Swap the xpart */
const struct xpart txp = s->xparts[k];
s->xparts[k] = s->xparts[nr_parts];
s->xparts[nr_parts] = txp;
/* Swap the index */
const int t = ind[k];
ind[k] = ind[nr_parts];
ind[nr_parts] = t;
......@@ -556,20 +560,67 @@ void space_rebuild(struct space *s, int verbose) {
}
#endif
/* Move non-local sparts to the end of the list. */
for (size_t k = 0; k < nr_sparts;) {
if (cells_top[sind[k]].nodeID != local_nodeID) {
nr_sparts -= 1;
/* Swap the particle */
const struct spart tp = s->sparts[k];
s->sparts[k] = s->sparts[nr_sparts];
s->sparts[nr_sparts] = tp;
/* Swap the link with the gpart */
if (s->sparts[k].gpart != NULL) {
s->sparts[k].gpart->id_or_neg_offset = -k;
}
if (s->sparts[nr_parts].gpart != NULL) {
s->sparts[nr_parts].gpart->id_or_neg_offset = -nr_parts;
}
/* Swap the index */
const int t = sind[k];
sind[k] = sind[nr_sparts];
sind[nr_sparts] = t;
} else {
/* Increment when not exchanging otherwise we need to retest "k".*/
k++;
}
}
#ifdef SWIFT_DEBUG_CHECKS
/* Check that all sparts are in the correct place (untested). */
for (size_t k = 0; k < nr_sparts; k++) {
if (cells_top[sind[k]].nodeID != local_nodeID) {
error("Failed to move all non-local sparts to send list");
}
}
for (size_t k = nr_sparts; k < s->nr_sparts; k++) {
if (cells_top[sind[k]].nodeID == local_nodeID) {
error("Failed to remove local sparts from send list");
}
}
#endif
/* Move non-local gparts to the end of the list. */
for (size_t k = 0; k < nr_gparts;) {
if (cells_top[gind[k]].nodeID != local_nodeID) {
nr_gparts -= 1;
/* Swap the particle */
const struct gpart tp = s->gparts[k];
s->gparts[k] = s->gparts[nr_gparts];
s->gparts[nr_gparts] = tp;
if (s->gparts[k].id_or_neg_offset <= 0) {
/* Swap the link with part/spart */
if (s->gparts[k].type == swift_type_gas) {
s->parts[-s->gparts[k].id_or_neg_offset].gpart = &s->gparts[k];
} else if (s->gparts[k].type == swift_type_star) {
s->sparts[-s->gparts[k].id_or_neg_offset].gpart = &s->gparts[k];
}
if (s->gparts[nr_gparts].id_or_neg_offset <= 0) {
if (s->gparts[nr_gparts].type == swift_type_gas) {
s->parts[-s->gparts[nr_gparts].id_or_neg_offset].gpart =
&s->gparts[nr_gparts];
} else if (s->gparts[nr_gparts].type == swift_type_star) {
s->sparts[-s->gparts[nr_gparts].id_or_neg_offset].gpart =
&s->gparts[nr_gparts];
}
/* Swap the index */
const int t = gind[k];
gind[k] = gind[nr_gparts];
gind[nr_gparts] = t;
......@@ -597,14 +648,17 @@ void space_rebuild(struct space *s, int verbose) {
the parts arrays. */
size_t nr_parts_exchanged = s->nr_parts - nr_parts;
size_t nr_gparts_exchanged = s->nr_gparts - nr_gparts;
size_t nr_sparts_exchanged = s->nr_sparts - nr_sparts;
engine_exchange_strays(s->e, nr_parts, &ind[nr_parts], &nr_parts_exchanged,
nr_gparts, &gind[nr_gparts], &nr_gparts_exchanged);
nr_gparts, &gind[nr_gparts], &nr_gparts_exchanged,
nr_sparts, &sind[nr_sparts], &nr_sparts_exchanged);
/* Set the new particle counts. */
s->nr_parts = nr_parts + nr_parts_exchanged;
s->nr_gparts = nr_gparts + nr_gparts_exchanged;
s->nr_sparts = nr_sparts + nr_sparts_exchanged;
/* Re-allocate the index array if needed.. */
/* Re-allocate the index array for the parts if needed.. */
if (s->nr_parts + 1 > ind_size) {
int *ind_new;
if ((ind_new = (int *)malloc(sizeof(int) * (s->nr_parts + 1))) == NULL)
......@@ -614,10 +668,20 @@ void space_rebuild(struct space *s, int verbose) {
ind = ind_new;
}
/* Re-allocate the index array for the sparts if needed.. */
if (s->nr_sparts + 1 > ind_size) {
int *sind_new;
if ((sind_new = (int *)malloc(sizeof(int) * (s->nr_sparts + 1))) == NULL)
error("Failed to allocate temporary particle indices.");
memcpy(sind_new, sind, sizeof(int) * nr_sparts);
free(sind);
sind = sind_new;
}
const int cdim[3] = {s->cdim[0], s->cdim[1], s->cdim[2]};
const double ih[3] = {s->iwidth[0], s->iwidth[1], s->iwidth[2]};
/* Assign each particle to its cell. */
/* Assign each part to its cell. */
for (size_t k = nr_parts; k < s->nr_parts; k++) {
const struct part *const p = &s->parts[k];
ind[k] =
......@@ -630,6 +694,19 @@ void space_rebuild(struct space *s, int verbose) {
}
nr_parts = s->nr_parts;
/* Assign each spart to its cell. */
for (size_t k = nr_sparts; k < s->nr_sparts; k++) {
const struct spart *const sp = &s->sparts[k];
sind[k] =
cell_getid(cdim, sp->x[0] * ih[0], sp->x[1] * ih[1], sp->x[2] * ih[2]);
#ifdef SWIFT_DEBUG_CHECKS
if (cells_top[sind[k]].nodeID != local_nodeID)
error("Received part that does not belong to me (nodeID=%i).",
cells_top[sind[k]].nodeID);
#endif
}
nr_sparts = s->nr_sparts;
#endif /* WITH_MPI */
/* Sort the parts according to their cells. */
......@@ -685,7 +762,7 @@ void space_rebuild(struct space *s, int verbose) {
#ifdef WITH_MPI
/* Re-allocate the index array if needed.. */
/* Re-allocate the index array for the gparts if needed.. */
if (s->nr_gparts + 1 > gind_size) {
int *gind_new;
if ((gind_new = (int *)malloc(sizeof(int) * (s->nr_gparts + 1))) == NULL)
......@@ -695,7 +772,7 @@ void space_rebuild(struct space *s, int verbose) {
gind = gind_new;
}
/* Assign each particle to its cell. */
/* Assign each gpart to its cell. */
for (size_t k = nr_gparts; k < s->nr_gparts; k++) {
const struct gpart *const p = &s->gparts[k];
gind[k] =
......
......@@ -145,6 +145,10 @@ struct space {
struct gpart *gparts_foreign;
size_t nr_gparts_foreign, size_gparts_foreign;
/*! Buffers for g-parts that we will receive from foreign cells. */
struct spart *sparts_foreign;
size_t nr_sparts_foreign, size_sparts_foreign;
#endif
};
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment