ell_relax: fill in the corner ghostpoints

1 files changed, 34 insertions, 0 deletions

diff --git a/ell_relax.c b/ell_relax.c
index 679f340..9454e84 100644
--- a/ell_relax.c
+++ b/ell_relax.c
@@ -170,6 +170,40 @@ static void boundaries_apply(EllRelaxContext *ctx)
             break;
         }
     }
+
+    /* fill in the corner ghosts */
+    {
+        double *dst = ctx->u;
+        for (int i = 1; i <= FD_STENCIL_MAX; i++)
+            for (int j = 1; j <= FD_STENCIL_MAX; j++) {
+                const ptrdiff_t idx = j * strides[1] + i;
+                dst[-idx] = dst[idx];
+            }
+    }
+    {
+        double *dst = ctx->u + ctx->domain_size[0] - 1;
+        for (int i = 1; i <= FD_STENCIL_MAX; i++)
+            for (int j = 1; j <= FD_STENCIL_MAX; j++) {
+                const ptrdiff_t idx = j * strides[1] - i;
+                dst[-idx] = dst[idx];
+            }
+    }
+    {
+        double *dst = ctx->u + strides[1] * (ctx->domain_size[1] - 1);
+        for (int i = 1; i <= FD_STENCIL_MAX; i++)
+            for (int j = 1; j <= FD_STENCIL_MAX; j++) {
+                const ptrdiff_t idx = -j * strides[1] + i;
+                dst[-idx] = dst[idx];
+            }
+    }
+    {
+        double *dst = ctx->u + strides[1] * (ctx->domain_size[1] - 1) + ctx->domain_size[0] - 1;
+        for (int i = 1; i <= FD_STENCIL_MAX; i++)
+            for (int j = 1; j <= FD_STENCIL_MAX; j++) {
+                const ptrdiff_t idx = -j * strides[1] - i;
+                dst[-idx] = dst[idx];
+            }
+    }
 }
 
 int mg2di_ell_relax_step(EllRelaxContext *ctx)