transfer: implement and use 1D interpolation

Stop abusing "2D of y size 1" for this.

1 files changed, 16 insertions, 22 deletions

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;