diff options
Diffstat (limited to 'Carpet/CarpetMask')
| -rw-r--r-- | Carpet/CarpetMask/interface.ccl | 2 | ||||
| -rw-r--r-- | Carpet/CarpetMask/src/mask_excluded.cc | 20 | ||||
| -rw-r--r-- | Carpet/CarpetMask/src/mask_surface.cc | 22 |
3 files changed, 40 insertions, 4 deletions
diff --git a/Carpet/CarpetMask/interface.ccl b/Carpet/CarpetMask/interface.ccl index 7297a1c46..9cd31f3bd 100644 --- a/Carpet/CarpetMask/interface.ccl +++ b/Carpet/CarpetMask/interface.ccl @@ -3,3 +3,5 @@ IMPLEMENTS: CarpetMask INHERITS: grid SphericalSurface + +USES INCLUDE HEADER: CarpetReduce_bits.h diff --git a/Carpet/CarpetMask/src/mask_excluded.cc b/Carpet/CarpetMask/src/mask_excluded.cc index 1b1c0beb5..aec5ffa6c 100644 --- a/Carpet/CarpetMask/src/mask_excluded.cc +++ b/Carpet/CarpetMask/src/mask_excluded.cc @@ -8,6 +8,8 @@ #include "mask_excluded.hh" +#include <CarpetReduce_bits.h> + namespace CarpetMask { @@ -26,15 +28,27 @@ namespace CarpetMask { DECLARE_CCTK_ARGUMENTS; DECLARE_CCTK_PARAMETERS; - CCTK_INT * const iweight = + CCTK_INT * restrict const iweight = static_cast <CCTK_INT *> (CCTK_VarDataPtr (cctkGH, 0, "CarpetReduce::iweight")); - if (not iweight) { + CCTK_REAL * restrict const excised_cells = + static_cast <CCTK_REAL *> + (CCTK_VarDataPtr (cctkGH, 0, "CarpetReduce::excised_cells")); + + if (not iweight or not excised_cells) { CCTK_WARN (CCTK_WARN_ABORT, "CarpetReduce is not active, or CarpetReduce::iweight does not have storage"); } + // Volume of a grid cell on this level, in terms of coarse grid + // cells + CCTK_REAL const cell_volume = + 1.0 / (cctk_levfac[0] * cctk_levfac[1] * cctk_levfac[2]); + + unsigned const bits = BMSK(cctk_dim); + CCTK_REAL const factor = 1.0 / bits; + for (int n = 0; n < 10; ++ n) { CCTK_REAL const r0 = excluded_radius[n]; @@ -60,6 +74,8 @@ namespace CarpetMask { dx2 + dy2 + dz2 >= r2 : dx2 + dy2 + dz2 <= r2) { + // Tally up the weight we are removing + * excised_cells += cell_volume * factor * BCNT(iweight[ind]); iweight[ind] = 0; } diff --git a/Carpet/CarpetMask/src/mask_surface.cc b/Carpet/CarpetMask/src/mask_surface.cc index 02d80f9f8..f2a592b3c 100644 --- a/Carpet/CarpetMask/src/mask_surface.cc +++ b/Carpet/CarpetMask/src/mask_surface.cc @@ -11,6 +11,8 @@ #include "mask_surface.hh" +#include <CarpetReduce_bits.h> + namespace CarpetMask { @@ -67,15 +69,27 @@ namespace CarpetMask { last_output.resize (num_excluded, -1); } - CCTK_INT * const iweight = + CCTK_INT * restrict const iweight = static_cast <CCTK_INT *> (CCTK_VarDataPtr (cctkGH, 0, "CarpetReduce::iweight")); - if (not iweight) { + CCTK_REAL * restrict const excised_cells = + static_cast <CCTK_REAL *> + (CCTK_VarDataPtr (cctkGH, 0, "CarpetReduce::excised_cells")); + + if (not iweight or not excised_cells) { CCTK_WARN (CCTK_WARN_ABORT, "CarpetReduce is not active, or CarpetReduce::iweight does not have storage"); } + // Volume of a grid cell on this level, in terms of coarse grid + // cells + CCTK_REAL const cell_volume = + 1.0 / (cctk_levfac[0] * cctk_levfac[1] * cctk_levfac[2]); + + unsigned const bits = BMSK(cctk_dim); + CCTK_REAL const factor = 1.0 / bits; + for (int n = 0; n < num_excluded; ++ n) { int const sn = excluded_surface[n]; @@ -127,6 +141,8 @@ namespace CarpetMask { sqrt (pow (dx, 2) + pow (dy, 2) + pow (dz, 2)); if (rho < 1.0e-12) { // Always excise the surface origin + // Tally up the weight we are removing + * excised_cells += cell_volume * factor * BCNT(iweight[ind]); iweight[ind] = 0; } else { CCTK_REAL theta = @@ -171,6 +187,8 @@ namespace CarpetMask { CCTK_REAL const dr = sf_radius[a + maxntheta * (b + maxnphi * sn)]; if (rho <= dr * shrink_factor) { + // Tally up the weight we are removing + * excised_cells += cell_volume * factor * BCNT(iweight[ind]); iweight[ind] = 0; } } |