BOOLEAN do_reflux "Manage flux registers to ensure conservation"
{
} yes

array<int, dim> get_group_fluxes(const int gi) {
  assert(gi >= 0);
  const int tags = CCTK_GroupTagsTableI(gi);
  assert(tags >= 0);
  vector<char> fluxes_buf(1000);
  const int iret =
      Util_TableGetString(tags, fluxes_buf.size(), fluxes_buf.data(), "fluxes");
  if (iret == UTIL_ERROR_TABLE_NO_SUCH_KEY) {
    fluxes_buf[0] = '\0'; // default: empty (no fluxes)
  } else if (iret >= 0) {
    // do nothing
  } else {
    assert(0);
  }

  const string str(fluxes_buf.data());
  vector<string> strs;
  size_t end = 0;
  while (end < str.size()) {
    size_t begin = str.find_first_not_of(' ', end);
    if (begin == string::npos)
      break;
    end = str.find(' ', begin);
    strs.push_back(str.substr(begin, end - begin));
  }
  array<int, dim> fluxes;
  fluxes.fill(-1);
  if (strs.empty())
    return fluxes; // No fluxes specified
  assert(strs.size() == dim); // Check number of fluxes
  for (int d = 0; d < dim; ++d) {
    auto str1 = strs[d];
    if (str1.find(':') == string::npos) {
      const char *impl = CCTK_GroupImplementationI(gi);
      str1 = string(impl) + "::" + str1;
    }
    int gi1 = CCTK_GroupIndex(str1.c_str());
    assert(gi1 >= 0); // Check fluxes are valid groups
    fluxes[d] = gi1;
  }
  for (int d = 0; d < dim; ++d)
    for (int d1 = d + 1; d1 < dim; ++d1)
      assert(fluxes[d] != fluxes[d1]); // Check groups are all different
  return fluxes;
}

    }
    if (level > 0) {
      array<int, dim> fluxes = get_group_fluxes(groupdata.groupindex);
      const bool have_fluxes = fluxes[0] >= 0;
      if (have_fluxes) {
        assert((groupdata.indextype == array<int, dim>{1, 1, 1}));
        groupdata.freg = make_unique<FluxRegister>(
            gba, dm, ghext->amrcore->refRatio(level - 1), level,
            groupdata.numvars);
        groupdata.fluxes = fluxes;
      } else {
        groupdata.fluxes.fill(-1);
      }
    }
  }
  // Check flux register consistency
  for (int gi = 0; gi < numgroups; ++gi) {
    const auto &groupdata = leveldata.groupdata.at(gi);
    if (groupdata.freg) {
      for (int d = 0; d < dim; ++d) {
        assert(groupdata.fluxes[d] != groupdata.groupindex);
        const auto &flux_groupdata =
            leveldata.groupdata.at(groupdata.fluxes[d]);
        array<int, dim> flux_indextype{1, 1, 1};
        flux_indextype[d] = 0;
        assert(flux_groupdata.indextype == flux_indextype);
        assert(flux_groupdata.numvars == groupdata.numvars);
      }

      if (groupdata.freg)
        groupdata.freg = make_unique<FluxRegister>(
            gba, dm, ghext->amrcore->refRatio(level - 1), level,
            groupdata.numvars);

#include <AMReX_FluxRegister.H>

      // flux register between this and the next coarser level
      unique_ptr<FluxRegister> freg;
      // associated flux group indices
      array<int, dim> fluxes; // [dir]

void Reflux(int level);

    // Reflux
    assert(current_level == -1);
    for (int level = int(ghext->leveldata.size()) - 2; level >= 0; --level)
      Reflux(level);

}
void Reflux(int level) {
  DECLARE_CCTK_PARAMETERS;
  if (!do_reflux)
    return;
  CCTK_VINFO("Refluxing level %d", level);
#warning "TODO"
  bool didprint = false;
  static Timer timer("Reflux");
  Interval interval(timer);
  auto &leveldata = ghext->leveldata.at(level);
  const auto &fineleveldata = ghext->leveldata.at(level + 1);
  for (int gi = 0; gi < int(leveldata.groupdata.size()); ++gi) {
    const int tl = 0;
    auto &groupdata = leveldata.groupdata.at(gi);
    const auto &finegroupdata = fineleveldata.groupdata.at(gi);
    // If the group has associated fluxes
    if (finegroupdata.freg) {
      if (!didprint) {
        CCTK_VINFO("  found flux registers");
        didprint = true;
      }
      // Check coarse and fine data and fluxes are valid
      for (int vi = 0; vi < finegroupdata.numvars; ++vi) {
        assert(finegroupdata.valid.at(tl).at(vi).valid_int &&
               finegroupdata.valid.at(tl).at(vi).valid_bnd);
        assert(groupdata.valid.at(tl).at(vi).valid_int &&
               groupdata.valid.at(tl).at(vi).valid_bnd);
      }
      for (int d = 0; d < dim; ++d) {
        int flux_gi = finegroupdata.fluxes[d];
        const auto &flux_finegroupdata = fineleveldata.groupdata.at(flux_gi);
        const auto &flux_groupdata = leveldata.groupdata.at(flux_gi);
        for (int vi = 0; vi < finegroupdata.numvars; ++vi) {
          assert(flux_finegroupdata.valid.at(tl).at(vi).valid_int &&
                 flux_finegroupdata.valid.at(tl).at(vi).valid_bnd);
          assert(flux_groupdata.valid.at(tl).at(vi).valid_int &&
                 flux_groupdata.valid.at(tl).at(vi).valid_bnd);
        }
      }
#warning "TODO"
      const Geometry &geom0 = ghext->amrcore->Geom(level);
      const CCTK_REAL *restrict gdx = geom0.CellSize();
      CCTK_REAL dt = 0.25 * pow(0.5, 2) * gdx[0];
      CCTK_REAL dx = 1.0 * pow(0.5, level) * gdx[0];
      for (int d = 0; d < dim; ++d) {
        int flux_gi = finegroupdata.fluxes[d];
        const auto &flux_finegroupdata = fineleveldata.groupdata.at(flux_gi);
        const auto &flux_groupdata = leveldata.groupdata.at(flux_gi);
        finegroupdata.freg->CrseInit(*flux_groupdata.mfab.at(tl), d, 0, 0,
                                     flux_groupdata.numvars, -dt * dx * dx);
        finegroupdata.freg->FineAdd(*flux_finegroupdata.mfab.at(tl), d, 0, 0,
                                    flux_finegroupdata.numvars,
                                    dt * (0.5 * dx) * (0.5 * dx));
      }
      const Geometry &geom = ghext->amrcore->Geom(level);
      finegroupdata.freg->Reflux(*groupdata.mfab.at(tl), 1.0, 0, 0,
                                 groupdata.numvars, geom);
    }
  }

    // If there is more than one time level, then we don't restrict
    // the oldest.

    // If there is more than one time level, then we don't restrict the oldest.

      // Only restrict valid grid functions. Restriction only uses the interior.

      // Only restrict valid grid functions. Restriction only uses the
      // interior.