From 60d1383db4d9f646ce5504edc9b2c437836696db Mon Sep 17 00:00:00 2001 From: Anton Khirnov Date: Sat, 29 Jun 2019 10:23:40 +0200 Subject: transfer: implement and use 1D interpolation Stop abusing "2D of y size 1" for this. --- mg2d.c | 38 ++++++++++++++++---------------------- 1 file changed, 16 insertions(+), 22 deletions(-) (limited to 'mg2d.c') diff --git a/mg2d.c b/mg2d.c index 113830b..853e9ab 100644 --- a/mg2d.c +++ b/mg2d.c @@ -640,18 +640,18 @@ static int diff_coeffs_fixup(MG2DContext *ctx, MG2DLevel *l) double *dst = a_dst->data + mg2d_bnd_is_upper(bnd_idx) * ((a_dst->shape[!ci] - 1) * a_dst->stride[!ci]); const ptrdiff_t stride_dst = a_dst->stride[ci]; - const size_t size_dst[] = { 1, l->solver->domain_size[!ci] }; + const size_t size_dst = l->solver->domain_size[!ci]; ret = mg2di_ndarray_slice(&dc_src, priv->dc_bnd_val[i][bnd_idx], - (Slice [2]){ SLICE_NULL, SLICE(FD_STENCIL_MAX, -FD_STENCIL_MAX, 1) }); + &SLICE(FD_STENCIL_MAX, -FD_STENCIL_MAX, 1)); if (ret < 0) goto fail; - ret = mg2di_ndarray_alloc(&dc_dst, 2, size_dst, 0); + ret = mg2di_ndarray_alloc(&dc_dst, 1, &size_dst, 0); if (ret < 0) goto fail; - gt_bnd = mg2di_gt_alloc(GRID_TRANSFER_LAGRANGE_1); + gt_bnd = mg2di_gt_alloc(1, GRID_TRANSFER_LAGRANGE_1); if (!gt_bnd) { ret = -ENOMEM; goto fail; @@ -661,18 +661,12 @@ static int diff_coeffs_fixup(MG2DContext *ctx, MG2DLevel *l) gt_bnd->cpuflags = priv->cpuflags; gt_bnd->src.start[0] = 0; - gt_bnd->src.start[1] = 0; - gt_bnd->src.size[0] = 1; - gt_bnd->src.size[1] = priv->dg->domain_size[1]; - gt_bnd->src.step[0] = ctx->step[0]; - gt_bnd->src.step[1] = ctx->step[1]; - - gt_bnd->dst.start[0] = 0; - gt_bnd->dst.start[1] = l->dg->components[priv->local_component].interior.start[!ci]; - gt_bnd->dst.size[0] = 1; - gt_bnd->dst.size[1] = l->dg->components[priv->local_component].interior.size[!ci]; - gt_bnd->dst.step[0] = l->solver->step[0]; - gt_bnd->dst.step[1] = l->solver->step[1]; + gt_bnd->src.size[0] = priv->dg->domain_size[1]; + gt_bnd->src.step[0] = ctx->step[1]; + + gt_bnd->dst.start[0] = l->dg->components[priv->local_component].interior.start[!ci]; + gt_bnd->dst.size[0] = l->dg->components[priv->local_component].interior.size[!ci]; + gt_bnd->dst.step[0] = l->solver->step[1]; ret = mg2di_gt_init(gt_bnd); if (ret < 0) @@ -682,7 +676,7 @@ static int diff_coeffs_fixup(MG2DContext *ctx, MG2DLevel *l) if (ret < 0) goto fail; - for (int idx = 0; idx < dc_dst->shape[1]; idx++) + for (int idx = 0; idx < dc_dst->shape[0]; idx++) dst[stride_dst * idx] += dc_dst->data[idx]; fail: @@ -730,7 +724,7 @@ static int mg_setup_restrict(MG2DContext *ctx, MG2DLevel *l, enum GridTransferOp GridTransferContext *gt; int ret; - gt = mg2di_gt_alloc(op); + gt = mg2di_gt_alloc(2, op); if (!gt) return -ENOMEM; l->transfer_restrict = gt; @@ -764,7 +758,7 @@ static int mg_setup_prolong(MG2DContext *ctx, MG2DLevel *l, enum GridTransferOpe GridTransferContext *gt; int ret; - gt = mg2di_gt_alloc(op); + gt = mg2di_gt_alloc(2, op); if (!gt) return -ENOMEM; l->transfer_prolong = gt; @@ -1642,8 +1636,8 @@ static MG2DContext *solver_alloc(DomainGeometry *dg, unsigned int local_componen for (int bnd_idx = 0; bnd_idx < ARRAY_ELEMS(dc->boundaries); bnd_idx++) { const int ci = mg2d_bnd_coord_idx(bnd_idx); - ret = mg2di_ndarray_alloc(&priv->dc_bnd_val[i][bnd_idx], 2, - (size_t [2]){ 1, dg->domain_size[!ci] + 2 * FD_STENCIL_MAX }, + ret = mg2di_ndarray_alloc(&priv->dc_bnd_val[i][bnd_idx], 1, + (size_t [1]){ dg->domain_size[!ci] + 2 * FD_STENCIL_MAX }, NDARRAY_ALLOC_ZERO); if (ret < 0) goto fail; @@ -1994,7 +1988,7 @@ int mg2d_init_guess(MG2DContext *ctx, const double *src, } if (!priv->transfer_init) { - priv->transfer_init = mg2di_gt_alloc(GRID_TRANSFER_LAGRANGE_3); + priv->transfer_init = mg2di_gt_alloc(2, GRID_TRANSFER_LAGRANGE_3); if (!priv->transfer_init) return -ENOMEM; -- cgit v1.2.3