diff options
Diffstat (limited to 'ell_grid_solve.c')
-rw-r--r-- | ell_grid_solve.c | 148 |
1 files changed, 88 insertions, 60 deletions
diff --git a/ell_grid_solve.c b/ell_grid_solve.c index 721dc11..61dc336 100644 --- a/ell_grid_solve.c +++ b/ell_grid_solve.c @@ -108,6 +108,8 @@ struct EGSInternal { DomainGeometry *dg; unsigned int local_component; + size_t fd_stencil; + int *sync_sendcounts; int *sync_senddispl; MPI_Datatype *sync_sendtypes; @@ -845,6 +847,74 @@ static int init_diff_coeffs_task(void *arg, unsigned int job_idx, unsigned int t return 0; } +static int init_mpi_sync(EGSContext *ctx) +{ + EGSInternal *priv = ctx->priv; + const DomainGeometry *dg = priv->dg; + Rect *overlaps_recv = NULL, *overlaps_send = NULL; + ptrdiff_t *lo; + int ret = 0; + + overlaps_recv = calloc(dg->nb_components, sizeof(*overlaps_recv)); + overlaps_send = calloc(dg->nb_components, sizeof(*overlaps_send)); + if (!overlaps_recv || !overlaps_send) { + ret = -ENOMEM; + goto fail; + } + + ret = mg2di_dg_edge_overlaps(overlaps_recv, overlaps_send, + dg, priv->local_component, FD_STENCIL_MAX); + if (ret < 0) + goto fail; + + priv->sync_sendtypes = calloc(dg->nb_components, sizeof(*priv->sync_sendtypes)); + priv->sync_senddispl = calloc(dg->nb_components, sizeof(*priv->sync_senddispl)); + priv->sync_sendcounts = calloc(dg->nb_components, sizeof(*priv->sync_sendcounts)); + priv->sync_recvtypes = calloc(dg->nb_components, sizeof(*priv->sync_recvtypes)); + priv->sync_recvdispl = calloc(dg->nb_components, sizeof(*priv->sync_recvdispl)); + priv->sync_recvcounts = calloc(dg->nb_components, sizeof(*priv->sync_recvcounts)); + if (!priv->sync_sendtypes || !priv->sync_senddispl || !priv->sync_sendcounts || + !priv->sync_recvtypes || !priv->sync_recvdispl || !priv->sync_recvcounts) + goto fail; + + lo = dg->components[priv->local_component].interior.start; + + /* construct the send/receive parameters */ + for (unsigned int i = 0; i < dg->nb_components; i++) { + if (i == priv->local_component) { + MPI_Type_dup(MPI_INT, &priv->sync_sendtypes[i]); + MPI_Type_dup(MPI_INT, &priv->sync_recvtypes[i]); + priv->sync_sendcounts[i] = 0; + priv->sync_recvcounts[i] = 0; + priv->sync_senddispl[i] = 0; + priv->sync_recvdispl[i] = 0; + + continue; + } + + /* receive */ + MPI_Type_vector(overlaps_recv[i].size[1], overlaps_recv[i].size[0], + ctx->u->stride[0], MPI_DOUBLE, &priv->sync_recvtypes[i]); + MPI_Type_commit(&priv->sync_recvtypes[i]); + priv->sync_recvcounts[i] = 1; + priv->sync_recvdispl[i] = ((overlaps_recv[i].start[1] - lo[1]) * ctx->u->stride[0] + + (overlaps_recv[i].start[0] - lo[0])) * sizeof(*ctx->u->data); + + /* send */ + MPI_Type_vector(overlaps_send[i].size[1], overlaps_send[i].size[0], + ctx->u->stride[0], MPI_DOUBLE, &priv->sync_sendtypes[i]); + MPI_Type_commit(&priv->sync_sendtypes[i]); + priv->sync_sendcounts[i] = 1; + priv->sync_senddispl[i] = ((overlaps_send[i].start[1] - lo[1]) * ctx->u->stride[0] + + (overlaps_send[i].start[0] - lo[0])) * sizeof(*ctx->u->data); + } + +fail: + free(overlaps_recv); + free(overlaps_send); + return ret; +} + int mg2di_egs_init(EGSContext *ctx, int flags) { EGSInternal *priv = ctx->priv; @@ -861,6 +931,22 @@ int mg2di_egs_init(EGSContext *ctx, int flags) goto finish; } + /* initialize the FD stencil-dependent state */ + if (!priv->fd_stencil) { + if (priv->dg->nb_components > 1) { + ret = init_mpi_sync(ctx); + if (ret < 0) + goto finish; + } + + priv->fd_stencil = ctx->fd_stencil; + } + if (priv->fd_stencil != ctx->fd_stencil) { + mg2di_log(&ctx->logger, MG2D_LOG_ERROR, "FD stencil changed\n"); + ret = -EINVAL; + goto finish; + } + if (r->relax_factor == 0.0) priv->r.relax_factor = relax_factors[ctx->fd_stencil - 1]; else @@ -1178,75 +1264,17 @@ EGSContext *mg2di_egs_alloc(const size_t domain_size[2]) EGSContext *mg2di_egs_alloc_mpi(MPI_Comm comm, const DomainGeometry *dg) { EGSContext *ctx = NULL; - Rect *overlaps_recv = NULL, *overlaps_send = NULL; - ptrdiff_t *lo; int local_component; - int ret; MPI_Comm_rank(comm, &local_component); - overlaps_recv = calloc(dg->nb_components, sizeof(*overlaps_recv)); - overlaps_send = calloc(dg->nb_components, sizeof(*overlaps_send)); - if (!overlaps_recv || !overlaps_send) - goto fail; - - ret = mg2di_dg_edge_overlaps(overlaps_recv, overlaps_send, - dg, local_component, FD_STENCIL_MAX); - if (ret < 0) - goto fail; - ctx = egs_alloc(dg, local_component); + if (!ctx) + return NULL; ctx->priv->comm = comm; - ctx->priv->sync_sendtypes = calloc(dg->nb_components, sizeof(*ctx->priv->sync_sendtypes)); - ctx->priv->sync_senddispl = calloc(dg->nb_components, sizeof(*ctx->priv->sync_senddispl)); - ctx->priv->sync_sendcounts = calloc(dg->nb_components, sizeof(*ctx->priv->sync_sendcounts)); - ctx->priv->sync_recvtypes = calloc(dg->nb_components, sizeof(*ctx->priv->sync_recvtypes)); - ctx->priv->sync_recvdispl = calloc(dg->nb_components, sizeof(*ctx->priv->sync_recvdispl)); - ctx->priv->sync_recvcounts = calloc(dg->nb_components, sizeof(*ctx->priv->sync_recvcounts)); - if (!ctx->priv->sync_sendtypes || !ctx->priv->sync_senddispl || !ctx->priv->sync_sendcounts || - !ctx->priv->sync_recvtypes || !ctx->priv->sync_recvdispl || !ctx->priv->sync_recvcounts) - goto fail; - - lo = dg->components[local_component].interior.start; - - /* construct the send/receive parameters */ - for (unsigned int i = 0; i < dg->nb_components; i++) { - if (i == local_component) { - MPI_Type_dup(MPI_INT, &ctx->priv->sync_sendtypes[i]); - MPI_Type_dup(MPI_INT, &ctx->priv->sync_recvtypes[i]); - ctx->priv->sync_sendcounts[i] = 0; - ctx->priv->sync_recvcounts[i] = 0; - ctx->priv->sync_senddispl[i] = 0; - ctx->priv->sync_recvdispl[i] = 0; - - continue; - } - - /* receive */ - MPI_Type_vector(overlaps_recv[i].size[1], overlaps_recv[i].size[0], - ctx->u->stride[0], MPI_DOUBLE, &ctx->priv->sync_recvtypes[i]); - MPI_Type_commit(&ctx->priv->sync_recvtypes[i]); - ctx->priv->sync_recvcounts[i] = 1; - ctx->priv->sync_recvdispl[i] = ((overlaps_recv[i].start[1] - lo[1]) * ctx->u->stride[0] + - (overlaps_recv[i].start[0] - lo[0])) * sizeof(*ctx->u->data); - - /* send */ - MPI_Type_vector(overlaps_send[i].size[1], overlaps_send[i].size[0], - ctx->u->stride[0], MPI_DOUBLE, &ctx->priv->sync_sendtypes[i]); - MPI_Type_commit(&ctx->priv->sync_sendtypes[i]); - ctx->priv->sync_sendcounts[i] = 1; - ctx->priv->sync_senddispl[i] = ((overlaps_send[i].start[1] - lo[1]) * ctx->u->stride[0] + - (overlaps_send[i].start[0] - lo[0])) * sizeof(*ctx->u->data); - } - return ctx; -fail: - free(overlaps_recv); - free(overlaps_send); - mg2di_egs_free(&ctx); - return NULL; } void mg2di_egs_free(EGSContext **pctx) |