diff options
Diffstat (limited to 'ell_relax.c')
| -rw-r--r-- | ell_relax.c | 61 |
1 files changed, 27 insertions, 34 deletions
diff --git a/ell_relax.c b/ell_relax.c index 982a95a..4ee2747 100644 --- a/ell_relax.c +++ b/ell_relax.c @@ -28,6 +28,8 @@ #include "cpu.h" #include "ell_relax.h" #include "log.h" +#include "mg2d_boundary.h" +#include "mg2d_boundary_internal.h" #include "mg2d_constants.h" static const double relax_factors[FD_STENCIL_MAX] = { @@ -42,8 +44,6 @@ struct EllRelaxInternal { double *residual_base; double *diff_coeffs_base[MG2D_DIFF_COEFF_NB]; - double *boundary_val_base[4]; - double *residual_max; size_t residual_max_size; @@ -298,25 +298,25 @@ static void boundaries_apply(EllRelaxContext *ctx) const ptrdiff_t dst_strides[] = { stride_boundary, (boundary_def[i].is_upper ? 1 : -1) * stride_offset }; - switch (ctx->boundaries[i].type) { + switch (ctx->boundaries[i]->type) { case MG2D_BC_TYPE_FIXVAL: - boundaries_apply_fixval(dst, dst_strides, ctx->boundaries[i].val, - ctx->boundaries[i].val_stride, size_boundary); + boundaries_apply_fixval(dst, dst_strides, ctx->boundaries[i]->val, + ctx->boundaries[i]->val_stride, size_boundary); break; case MG2D_BC_TYPE_FIXDIFF: - boundaries_apply_fixdiff(dst, dst_strides, ctx->boundaries[i].val, - ctx->boundaries[i].val_stride, size_boundary); + boundaries_apply_fixdiff(dst, dst_strides, ctx->boundaries[i]->val, + ctx->boundaries[i]->val_stride, size_boundary); break; } } /* fill in the corner ghosts */ - if (ctx->boundaries[MG2D_BOUNDARY_0L].type == MG2D_BC_TYPE_FIXDIFF || - ctx->boundaries[MG2D_BOUNDARY_1L].type == MG2D_BC_TYPE_FIXDIFF) { + if (ctx->boundaries[MG2D_BOUNDARY_0L]->type == MG2D_BC_TYPE_FIXDIFF || + ctx->boundaries[MG2D_BOUNDARY_1L]->type == MG2D_BC_TYPE_FIXDIFF) { double *dst = ctx->u; int fact_x = -1, fact_z = -1; - if (ctx->boundaries[MG2D_BOUNDARY_0L].type == MG2D_BC_TYPE_FIXDIFF) + if (ctx->boundaries[MG2D_BOUNDARY_0L]->type == MG2D_BC_TYPE_FIXDIFF) fact_z *= -1; else fact_x *= -1; @@ -328,12 +328,12 @@ static void boundaries_apply(EllRelaxContext *ctx) dst[idx_dst] = dst[idx_src]; } } - if (ctx->boundaries[MG2D_BOUNDARY_0L].type == MG2D_BC_TYPE_FIXDIFF || - ctx->boundaries[MG2D_BOUNDARY_1U].type == MG2D_BC_TYPE_FIXDIFF) { + if (ctx->boundaries[MG2D_BOUNDARY_0L]->type == MG2D_BC_TYPE_FIXDIFF || + ctx->boundaries[MG2D_BOUNDARY_1U]->type == MG2D_BC_TYPE_FIXDIFF) { double *dst = ctx->u + ctx->domain_size[0] - 1; int fact_x = 1, fact_z = -1; - if (ctx->boundaries[MG2D_BOUNDARY_0L].type == MG2D_BC_TYPE_FIXDIFF) + if (ctx->boundaries[MG2D_BOUNDARY_0L]->type == MG2D_BC_TYPE_FIXDIFF) fact_z *= -1; else fact_x *= -1; @@ -345,12 +345,12 @@ static void boundaries_apply(EllRelaxContext *ctx) dst[idx_dst] = dst[idx_src]; } } - if (ctx->boundaries[MG2D_BOUNDARY_1L].type == MG2D_BC_TYPE_FIXDIFF || - ctx->boundaries[MG2D_BOUNDARY_0U].type == MG2D_BC_TYPE_FIXDIFF) { + if (ctx->boundaries[MG2D_BOUNDARY_1L]->type == MG2D_BC_TYPE_FIXDIFF || + ctx->boundaries[MG2D_BOUNDARY_0U]->type == MG2D_BC_TYPE_FIXDIFF) { double *dst = ctx->u + strides[1] * (ctx->domain_size[1] - 1); int fact_x = -1, fact_z = 1; - if (ctx->boundaries[MG2D_BOUNDARY_0U].type == MG2D_BC_TYPE_FIXDIFF) + if (ctx->boundaries[MG2D_BOUNDARY_0U]->type == MG2D_BC_TYPE_FIXDIFF) fact_z *= -1; else fact_x *= -1; @@ -362,12 +362,12 @@ static void boundaries_apply(EllRelaxContext *ctx) dst[idx_dst] = dst[idx_src]; } } - if (ctx->boundaries[MG2D_BOUNDARY_0U].type == MG2D_BC_TYPE_FIXDIFF || - ctx->boundaries[MG2D_BOUNDARY_1U].type == MG2D_BC_TYPE_FIXDIFF) { + if (ctx->boundaries[MG2D_BOUNDARY_0U]->type == MG2D_BC_TYPE_FIXDIFF || + ctx->boundaries[MG2D_BOUNDARY_1U]->type == MG2D_BC_TYPE_FIXDIFF) { double *dst = ctx->u + strides[1] * (ctx->domain_size[1] - 1) + ctx->domain_size[0] - 1; int fact_x = 1, fact_z = 1; - if (ctx->boundaries[MG2D_BOUNDARY_0U].type == MG2D_BC_TYPE_FIXDIFF) + if (ctx->boundaries[MG2D_BOUNDARY_0U]->type == MG2D_BC_TYPE_FIXDIFF) fact_z *= -1; else fact_x *= -1; @@ -475,13 +475,13 @@ int mg2di_ell_relax_init(EllRelaxContext *ctx) priv->residual_calc_offset = 0; priv->residual_calc_size[0] = ctx->domain_size[0]; priv->residual_calc_size[1] = ctx->domain_size[1]; - if (ctx->boundaries[MG2D_BOUNDARY_0L].type == MG2D_BC_TYPE_FIXVAL) + if (ctx->boundaries[MG2D_BOUNDARY_0L]->type == MG2D_BC_TYPE_FIXVAL) priv->residual_calc_offset += ctx->residual_stride; - if (ctx->boundaries[MG2D_BOUNDARY_1L].type == MG2D_BC_TYPE_FIXVAL) + if (ctx->boundaries[MG2D_BOUNDARY_1L]->type == MG2D_BC_TYPE_FIXVAL) priv->residual_calc_offset++; for (int i = 0; i < ARRAY_ELEMS(ctx->boundaries); i++) { - EllRelaxBoundary *bnd = &ctx->boundaries[i]; + MG2DBoundary *bnd = ctx->boundaries[i]; if (bnd->type == MG2D_BC_TYPE_FIXDIFF) { for (int k = 0; k < ctx->domain_size[boundary_def[i].stride_idx]; k++) if (bnd->val[k] != 0.0) { @@ -553,16 +553,9 @@ static int ell_relax_arrays_alloc(EllRelaxContext *ctx, const size_t domain_size priv->stride = stride; for (int i = 0; i < ARRAY_ELEMS(ctx->boundaries); i++) { - const size_t boundary_len0 = domain_size[boundary_def[i].stride_idx]; - const size_t boundary_len_max = boundary_len0 + 2 * (FD_STENCIL_MAX - 1); - - ret = posix_memalign((void**)&priv->boundary_val_base[i], 32, - sizeof(*priv->boundary_val_base[i]) * boundary_len_max * FD_STENCIL_MAX); - if (ret != 0) - return -ret; - - ctx->boundaries[i].val = priv->boundary_val_base[i] + FD_STENCIL_MAX - 1; - ctx->boundaries[i].val_stride = boundary_len_max; + ctx->boundaries[i] = mg2di_bc_alloc(domain_size[boundary_def[i].stride_idx]); + if (!ctx->boundaries[i]) + return -ENOMEM; } return 0; @@ -613,8 +606,8 @@ void mg2di_ell_relax_free(EllRelaxContext **pctx) free(ctx->priv->residual_max); for (int i = 0; i < ARRAY_ELEMS(ctx->priv->diff_coeffs_base); i++) free(ctx->priv->diff_coeffs_base[i]); - for (int i = 0; i < ARRAY_ELEMS(ctx->priv->boundary_val_base); i++) - free(ctx->priv->boundary_val_base[i]); + for (int i = 0; i < ARRAY_ELEMS(ctx->boundaries); i++) + mg2di_bc_free(&ctx->boundaries[i]); tp_free(&ctx->priv->tp_internal); free(ctx->priv); |