TestLoopControl: Correct errors in testing ghost zones

1 files changed, 52 insertions, 34 deletions

diff --git a/CarpetExtra/TestLoopControl/src/TestLoopPointwise.c b/CarpetExtra/TestLoopControl/src/TestLoopPointwise.c
index 8ac6a26bd..fe9baeef5 100644
--- a/CarpetExtra/TestLoopControl/src/TestLoopPointwise.c
+++ b/CarpetExtra/TestLoopControl/src/TestLoopPointwise.c
@@ -24,7 +24,7 @@ void TestLoopControlPointwise_All(CCTK_ARGUMENTS)
     for (int j=0; j<cctk_lsh[1]; ++j) {
       for (int i=0; i<cctk_lsh[0]; ++i) {
         int const ind3d = CCTK_GFINDEX3D(cctkGH, i,j,k);
-        pointtypes[ind3d] = 1;
+        pointtypes[ind3d] = 100;
       }
     }
   }
@@ -32,7 +32,7 @@ void TestLoopControlPointwise_All(CCTK_ARGUMENTS)
 #pragma omp parallel
   CCTK_LOOP3_ALL(loop3_all, cctkGH, i,j,k) {
     int const ind3d = CCTK_GFINDEX3D(cctkGH, i,j,k);
-    pointtypes[ind3d] -= 1;
+    pointtypes[ind3d] += 1;
   } CCTK_ENDLOOP3_ALL(loop3_all);
   
   int num_errors = 0;
@@ -40,7 +40,7 @@ void TestLoopControlPointwise_All(CCTK_ARGUMENTS)
     for (int j=0; j<cctk_lsh[1]; ++j) {
       for (int i=0; i<cctk_lsh[0]; ++i) {
         int const ind3d = CCTK_GFINDEX3D(cctkGH, i,j,k);
-        num_errors += pointtypes[ind3d] != 0;
+        num_errors += pointtypes[ind3d] != 101;
       }
     }
   }
@@ -65,9 +65,8 @@ void TestLoopControlPointwise_Int(CCTK_ARGUMENTS)
   
   int imin[3], imax[3];
   for (int d=0; d<3; ++d) {
-    imin[d] = cctk_bbox[2*d  ] ? bndsize[2*d  ] : cctk_nghostzones[d];
-    imax[d] = (cctk_lsh[d] -
-               (cctk_bbox[2*d+1] ? bndsize[2*d+1] : cctk_nghostzones[d]));
+    imin[d] = bndsize[2*d];
+    imax[d] = cctk_lsh[d] - bndsize[2*d+1];
   }
   
   for (int k=0; k<cctk_lsh[2]; ++k) {
@@ -83,7 +82,7 @@ void TestLoopControlPointwise_Int(CCTK_ARGUMENTS)
     for (int j=imin[1]; j<imax[1]; ++j) {
       for (int i=imin[0]; i<imax[0]; ++i) {
         int const ind3d = CCTK_GFINDEX3D(cctkGH, i,j,k);
-        pointtypes[ind3d] = 1;
+        pointtypes[ind3d] = 100;
       }
     }
   }
@@ -91,7 +90,7 @@ void TestLoopControlPointwise_Int(CCTK_ARGUMENTS)
 #pragma omp parallel
   CCTK_LOOP3_INT(loop3_int, cctkGH, i,j,k) {
     int const ind3d = CCTK_GFINDEX3D(cctkGH, i,j,k);
-    pointtypes[ind3d] -= 1;
+    pointtypes[ind3d] += 1;
   } CCTK_ENDLOOP3_INT(loop3_int);
   
   int num_errors = 0;
@@ -99,7 +98,7 @@ void TestLoopControlPointwise_Int(CCTK_ARGUMENTS)
     for (int j=0; j<cctk_lsh[1]; ++j) {
       for (int i=0; i<cctk_lsh[0]; ++i) {
         int const ind3d = CCTK_GFINDEX3D(cctkGH, i,j,k);
-        num_errors += pointtypes[ind3d] != 0;
+        num_errors += pointtypes[ind3d] != 0 && pointtypes[ind3d] != 101;
       }
     }
   }
@@ -148,7 +147,7 @@ void TestLoopControlPointwise_Bnd(CCTK_ARGUMENTS)
           for (int j=imin[1]; j<imax[1]; ++j) {
             for (int i=imin[0]; i<imax[0]; ++i) {
               int const ind3d = CCTK_GFINDEX3D(cctkGH, i,j,k);
-              pointtypes[ind3d] = 1;
+              pointtypes[ind3d] = 100;
             }
           }
         }
@@ -159,7 +158,7 @@ void TestLoopControlPointwise_Bnd(CCTK_ARGUMENTS)
 #pragma omp parallel
   CCTK_LOOP3_BND(loop3_bnd, cctkGH, i,j,k, ni,nj,nk) {
     int const ind3d = CCTK_GFINDEX3D(cctkGH, i,j,k);
-    pointtypes[ind3d] -= 1;
+    pointtypes[ind3d] += 1;
   } CCTK_ENDLOOP3_BND(loop3_bnd);
   
   int num_errors = 0;
@@ -167,10 +166,7 @@ void TestLoopControlPointwise_Bnd(CCTK_ARGUMENTS)
     for (int j=0; j<cctk_lsh[1]; ++j) {
       for (int i=0; i<cctk_lsh[0]; ++i) {
         int const ind3d = CCTK_GFINDEX3D(cctkGH, i,j,k);
-        num_errors += pointtypes[ind3d] != 0;
-        if (pointtypes[ind3d] != 0) {
-          printf("[%d %d %d] %d\n", i,j,k, (int)pointtypes[ind3d]);
-        }
+        num_errors += pointtypes[ind3d] != 0 && pointtypes[ind3d] != 101;
       }
     }
   }
@@ -202,26 +198,48 @@ void TestLoopControlPointwise_IntBnd(CCTK_ARGUMENTS)
     }
   }
   
-  for (int dir=0; dir<3; ++dir) {
-    for (int face=0; face<2; ++face) {
-      if (is_physbnd[2*dir+face] &&
-          !is_ghostbnd[2*dir+face] && !is_symbnd[2*dir+face])
-      {
-        int imin[3], imax[3];
+  int idir[3];
+  for (idir[2]=-1; idir[2]<=+1; ++idir[2]) {
+    for (idir[1]=-1; idir[1]<=+1; ++idir[1]) {
+      for (idir[0]=-1; idir[0]<=+1; ++idir[0]) {
+        int is_any_physbnd  = 0;
+        int is_any_ghostbnd = 0;
+        int is_any_symbnd   = 0;
         for (int d=0; d<3; ++d) {
-          imin[d] = 0;
-          imax[d] = cctk_lsh[d];
-        }
-        if (face==0) {
-          imax[dir] = bndsize[2*dir];
-        } else {
-          imin[dir] = cctk_lsh[dir] - bndsize[2*dir+1];
+          if (idir[d] < 0) {
+            is_any_physbnd  |= is_physbnd [2*d];
+            is_any_ghostbnd |= is_ghostbnd[2*d];
+            is_any_symbnd   |= is_symbnd  [2*d];
+          } else if (idir[d] == 0) {
+            /* do nothing */
+          } else if (idir[d] > 0) {
+            is_any_physbnd  |= is_physbnd [2*d+1];
+            is_any_ghostbnd |= is_ghostbnd[2*d+1];
+            is_any_symbnd   |= is_symbnd  [2*d+1];
+          }
         }
-        for (int k=imin[2]; k<imax[2]; ++k) {
-          for (int j=imin[1]; j<imax[1]; ++j) {
-            for (int i=imin[0]; i<imax[0]; ++i) {
-              int const ind3d = CCTK_GFINDEX3D(cctkGH, i,j,k);
-              pointtypes[ind3d] = 1;
+        /* Iterate over this region if it is a physical boundary, and
+           if it is neither a ghost nor a symmetry boundary */
+        if (is_any_physbnd && !is_any_ghostbnd && !is_any_symbnd) {
+          int imin[3], imax[3];
+          for (int d=0; d<3; ++d) {
+            if (idir[d] < 0) {
+              imin[d] = 0;
+              imax[d] = bndsize[2*d];
+            } else if (idir[d] == 0) {
+              imin[d] = bndsize[2*d];
+              imax[d] = cctk_lsh[d] - bndsize[2*d+1];
+            } else if (idir[d] > 0) {
+              imin[d] = cctk_lsh[d] - bndsize[2*d+1];
+              imax[d] = cctk_lsh[d];
+            }
+          }
+          for (int k=imin[2]; k<imax[2]; ++k) {
+            for (int j=imin[1]; j<imax[1]; ++j) {
+              for (int i=imin[0]; i<imax[0]; ++i) {
+                int const ind3d = CCTK_GFINDEX3D(cctkGH, i,j,k);
+                pointtypes[ind3d] += 100;
+              }
             }
           }
         }
@@ -240,7 +258,7 @@ void TestLoopControlPointwise_IntBnd(CCTK_ARGUMENTS)
     for (int j=0; j<cctk_lsh[1]; ++j) {
       for (int i=0; i<cctk_lsh[0]; ++i) {
         int const ind3d = CCTK_GFINDEX3D(cctkGH, i,j,k);
-        num_errors += pointtypes[ind3d] != 0 && pointtypes[ind3d] != 2;
+        num_errors += pointtypes[ind3d] != 0 && pointtypes[ind3d] != 101;
       }
     }
   }