diff options
author | schnetter <schnetter@20f44201-0f4f-0410-9130-e5fc2714a787> | 2007-03-02 23:25:35 +0000 |
---|---|---|
committer | schnetter <schnetter@20f44201-0f4f-0410-9130-e5fc2714a787> | 2007-03-02 23:25:35 +0000 |
commit | 9cba03dae49af1e8e4127178a6af30089c5030db (patch) | |
tree | 77e2757cb165b984cd2c6e79e8dca637f14e1075 | |
parent | 655126e98267998083de5d436930af235ef79485 (diff) |
Don't use == to look for poison, use memcmp instead.
git-svn-id: http://svn.cactuscode.org/arrangements/CactusNumerical/RotatingSymmetry180/trunk@42 20f44201-0f4f-0410-9130-e5fc2714a787
-rw-r--r-- | param.ccl | 6 | ||||
-rw-r--r-- | src/rotatingsymmetry180.c | 66 |
2 files changed, 46 insertions, 26 deletions
@@ -27,7 +27,7 @@ BOOLEAN poison_boundaries "Fill the symmetry boundary with a poison value before { } "no" -CCTK_REAL poison_value "Poison value (must not occur naturally)" +CCTK_INT poison_value "Integer value (0..255) used to poison new timelevels (with memset)" STEERABLE=always { - *:* :: "" -} -424242.0 + 0:255 :: "Must fit into a byte. Use 0 for zero, 255 for nan, and e.g. 113 for a large value." +} 255 diff --git a/src/rotatingsymmetry180.c b/src/rotatingsymmetry180.c index 2779e7a..62b859d 100644 --- a/src/rotatingsymmetry180.c +++ b/src/rotatingsymmetry180.c @@ -32,7 +32,7 @@ int BndRot180VI (cGH const * restrict const cctkGH, void * restrict * restrict varptrs; int * restrict vartypes; - int (*paritiess)[3]; + int (* restrict paritiess)[3]; int global_bbox[6]; int global_lbnd[3], global_ubnd[3]; @@ -392,16 +392,18 @@ int BndRot180VI (cGH const * restrict const cctkGH, case CCTK_VARIABLE_INT: /* do nothing */ break; - case CCTK_VARIABLE_REAL: + case CCTK_VARIABLE_REAL: { + CCTK_REAL * restrict const varptr = varptrs[var]; for (k=imin[0]; k<imax[2]; ++k) { for (j=imin[1]; j<imax[1]; ++j) { for (i=imin[2]; i<imax[0]; ++i) { const int ind = CCTK_GFINDEX3D(cctkGH, i, j, k); - ((CCTK_REAL * restrict) varptrs[var]) [ind] = poison_value; + memset (&varptr[ind], poison_value, sizeof varptr[ind]); } } } break; + } case CCTK_VARIABLE_COMPLEX: /* do nothing */ break; @@ -479,20 +481,35 @@ int BndRot180VI (cGH const * restrict const cctkGH, case CCTK_VARIABLE_INT: /* do nothing */ break; - case CCTK_VARIABLE_REAL: + case CCTK_VARIABLE_REAL: { + int poison_found = 0; + CCTK_REAL const * restrict const varptr = varptrs[var]; + CCTK_REAL poison; + memset (&poison, poison_value, sizeof poison); for (k=imin[0]; k<imax[2]; ++k) { for (j=imin[1]; j<imax[1]; ++j) { for (i=imin[2]; i<imax[0]; ++i) { const int ind = CCTK_GFINDEX3D(cctkGH, i, j, k); - CCTK_REAL const val = - ((CCTK_REAL const * restrict) varptrs[var]) [ind]; - if (val == poison_value) { - CCTK_WARN (CCTK_WARN_ABORT, "Poison found in symmetry regions -- there is an error in this thorn"); + if (memcmp (&varptr[ind], &poison, sizeof poison) == 0) { + printf (" ijk=[%d,%d,%d]\n", i, j, k); + poison_found = 1; } } } } + if (poison_found) { + printf ("Poison found:\n"); + printf (" levfac=[%d,%d,%d]\n", cctk_levfac[0], cctk_levfac[1], cctk_levfac[2]); + printf (" origin_space=[%g,%g,%g]\n", cctk_origin_space[0], cctk_origin_space[1], cctk_origin_space[2]); + printf (" delta_space=[%g,%g,%g]\n", cctk_delta_space[0], cctk_delta_space[1], cctk_delta_space[2]); + printf (" lbnd=[%d,%d,%d]\n", cctk_lbnd[0], cctk_lbnd[1], cctk_lbnd[2]); + printf (" lsh=[%d,%d,%d]\n", cctk_lsh[0], cctk_lsh[1], cctk_lsh[2]); + printf (" gsh=[%d,%d,%d]\n", cctk_gsh[0], cctk_gsh[1], cctk_gsh[2]); + printf (" bbox=[%d,%d,%d,%d,%d,%d]\n", cctk_bbox[0], cctk_bbox[1], cctk_bbox[2], cctk_bbox[3], cctk_bbox[4], cctk_bbox[5]); + CCTK_WARN (CCTK_WARN_ABORT, "Poison found in symmetry regions -- there is an error in this thorn"); + } break; + } case CCTK_VARIABLE_COMPLEX: /* do nothing */ break; @@ -522,41 +539,44 @@ int BndRot180VI (cGH const * restrict const cctkGH, } assert (group.dim == 3); switch (group.vartype) { - case CCTK_VARIABLE_INT: + case CCTK_VARIABLE_INT: { + CCTK_INT * restrict const varptr = varptrs[var]; for (k=imin[0]; k<imax[2]; ++k) { for (j=imin[1]; j<imax[1]; ++j) { for (i=imin[2]; i<imax[0]; ++i) { const int ind = CCTK_GFINDEX3D(cctkGH, i, j, k); - ((CCTK_INT * restrict) varptrs[var]) [ind] *= -1; + varptr[ind] *= -1; } } } break; - case CCTK_VARIABLE_REAL: + } + case CCTK_VARIABLE_REAL: { + CCTK_REAL * restrict const varptr = varptrs[var]; for (k=imin[0]; k<imax[2]; ++k) { for (j=imin[1]; j<imax[1]; ++j) { for (i=imin[2]; i<imax[0]; ++i) { const int ind = CCTK_GFINDEX3D(cctkGH, i, j, k); - ((CCTK_REAL * restrict) varptrs[var]) [ind] *= -1; + varptr[ind] *= -1; } } } break; - case CCTK_VARIABLE_COMPLEX: - { - CCTK_COMPLEX const czero = CCTK_Cmplx (0.0, 0.0); - for (k=imin[0]; k<imax[2]; ++k) { - for (j=imin[1]; j<imax[1]; ++j) { - for (i=imin[2]; i<imax[0]; ++i) { - const int ind = CCTK_GFINDEX3D(cctkGH, i, j, k); - CCTK_COMPLEX * restrict const ptr = - & ((CCTK_COMPLEX * restrict) varptrs[var]) [ind]; - * ptr = CCTK_CmplxSub (czero, * ptr); - } + } + case CCTK_VARIABLE_COMPLEX: { + CCTK_COMPLEX const czero = CCTK_Cmplx (0.0, 0.0); + CCTK_COMPLEX * restrict const varptr = varptrs[var]; + for (k=imin[0]; k<imax[2]; ++k) { + for (j=imin[1]; j<imax[1]; ++j) { + for (i=imin[2]; i<imax[0]; ++i) { + const int ind = CCTK_GFINDEX3D(cctkGH, i, j, k); + CCTK_COMPLEX * restrict const ptr = & varptr[ind]; + * ptr = CCTK_CmplxSub (czero, * ptr); } } } break; + } default: assert (0); } /* switch grouptype */ |