eschnett/cactusamrex - Change M5R6KQLXLGYSVKHVAX5AJKD6NYE6IM5Z6WVTR3BTKPJDNNKF3ARAC

AMReX: add support for grid scalars first version that compiles, untested otherwise, has a couple of sub-optimal choices

Created by Roland Haas on October 14, 2019

M5R6KQLXLGYSVKHVAX5AJKD6NYE6IM5Z6WVTR3BTKPJDNNKF3ARAC

Dependencies

In channels

master

main

Change contents

Insertion in AMReX/src/driver.cxx at line 47 [6.15972]

[5.166]


// Local functions
void SetupLevel(int level, const BoxArray &ba, const DistributionMapping &dm);
void SetupGlobals();

Replacement in AMReX/src/driver.cxx at line 60 [6.15972]

B:BD[7.269] → [8.393:466]

    : AmrCore(rb, max_level_in, n_cell_in, coord, ref_ratios, is_per) {}

[7.269]

[8.466]

    : AmrCore(rb, max_level_in, n_cell_in, coord, ref_ratios, is_per) {
  SetupGlobals();
}

Replacement in AMReX/src/driver.cxx at line 68 [6.15972]

B:BD[7.613] → [8.533:606]

    : AmrCore(rb, max_level_in, n_cell_in, coord, ref_ratios, is_per) {}

[7.613]

[7.686]

    : AmrCore(rb, max_level_in, n_cell_in, coord, ref_ratios, is_per) {
  SetupGlobals();
}

Insertion in AMReX/src/driver.cxx at line 116 [6.15972]

[9.385]

[9.397]

void SetupGlobals() {
  DECLARE_CCTK_PARAMETERS;
  if (verbose)
    CCTK_VINFO("SetupGlobals");
  GHExt::GlobalData &globaldata = ghext->globaldata;
  const int numgroups = CCTK_NumGroups();
  globaldata.scalargroupdata.resize(numgroups);
  for (int gi = 0; gi < numgroups; ++gi) {
    cGroup group;
    int ierr = CCTK_GroupData(gi, &group);
    assert(!ierr);
    /* only grid functions live on levels (and the grid) */
    if(group.grouptype != CCTK_SCALAR and group.grouptype != CCTK_ARRAY)
      continue;
    assert(group.grouptype == CCTK_SCALAR);
    assert(group.vartype == CCTK_VARIABLE_REAL);
    assert(group.disttype == CCTK_DISTRIB_DEFAULT);
    assert(group.dim == dim);
    // TODO: vectors must be contiguous in memory and mine aren't
    assert(group.vectorlength == 1);
    GHExt::GlobalData::ScalarGroupData &scalargroupdata = globaldata.scalargroupdata.at(gi);
    scalargroupdata.groupindex = gi;
    scalargroupdata.firstvarindex = CCTK_FirstVarIndexI(gi);
    scalargroupdata.numvars = group.numvars;
    // Allocate data
    scalargroupdata.data.resize(group.numtimelevels);
    for (int tl = 0; tl < int(scalargroupdata.data.size()); ++tl) {
      scalargroupdata.data.at(tl).resize(scalargroupdata.numvars);
      for (int vi = 0; vi < scalargroupdata.numvars; ++vi) {
        // TODO: find out something that avoid new ?
        scalargroupdata.data.at(tl).at(vi) = new CCTK_REAL;
        // TODO: decide that valid_bnd == false always and rely on
        // initialization magic?
        scalargroupdata.valid.at(tl).at(vi).valid_int = false;
        scalargroupdata.valid.at(tl).at(vi).valid_bnd = true;
        // TODO: make poison_invalid and check_invalid virtual members of
        // CommonGroupData
        poison_invalid(scalargroupdata, vi, tl);
        check_valid(scalargroupdata, vi, tl);
      }
    }
  }
}

Insertion in AMReX/src/driver.cxx at line 257 [6.15972]

[10.1852]


    /* only grid functions live on levels (and the grid) */
    if(group.grouptype != CCTK_GF)
      continue;

Insertion in AMReX/src/driver.cxx at line 483 [6.15972]

[10.2992]

    cGroup group;
    int ierr = CCTK_GroupData(gi, &group);
    assert(!ierr);
    if (group.grouptype != CCTK_GF)
      continue;

Insertion in AMReX/src/driver.cxx at line 572 [6.15972]

[10.4368]

    cGroup group;
    int ierr = CCTK_GroupData(gi, &group);
    assert(!ierr);
    if (group.grouptype != CCTK_GF)
      continue;

Insertion in AMReX/src/driver.hxx at line 92 [6.15275]

[11.1600]

[12.1840]

  struct CommonGroupData {
    int groupindex;
    int firstvarindex;
    int numvars;
    vector<vector<valid_t> > valid; // [time level][var index]
    // TODO: add poison_invalid and check_valid functions
  };
  struct GlobalData {
    // all data that exists on all levels
    struct ScalarGroupData : public CommonGroupData {
      // TODO: find out how to do this without a pointer to a single double but
      // also working around const vector<CCTK_REAL> issues
      vector<vector<CCTK_REAL *> > data; // [time level][var index]
    };
    // TODO: right now this is sized for the total number of groups
    vector<ScalarGroupData> scalargroupdata; // [group index]
  };
  GlobalData globaldata;

Replacement in AMReX/src/driver.hxx at line 121 [6.15275]
B:BD[13.1] → [14.132:155]
B:BD[14.155] → [15.3683:3705]
∅:D[15.3705] → [14.155:199]
B:BD[14.155] → [14.155:199]
```
    struct GroupData {
      int groupindex;
      int firstvarindex;
      int numvars;
```
[14.34]
[16.179]
```
    struct GroupData : public CommonGroupData {
```
Deletion in AMReX/src/driver.hxx at line 126 [6.15275]
B:BD[17.403] → [18.2770:2839]
```
      vector<vector<valid_t> > valid;     // [time level][var index]
```

Replacement in AMReX/src/driver.hxx at line 132 [6.15275]

B:BD[14.252] → [14.252:285]

    vector<GroupData> groupdata;

[14.252]

[19.309]

    // TODO: right now this is sized for the total number of groups
    vector<GroupData> groupdata; // [group index]

Replacement in AMReX/src/driver.hxx at line 135 [6.15275]
B:BD[19.314] → [12.1862:1893]
```
  vector<LevelData> leveldata;
```
[19.314]
[6.15872]
```
  vector<LevelData> leveldata; // [reflevel]
```

Insertion in AMReX/src/io.cxx at line 38 [20.107]

[9.4498]

    cGroup group;
    int ierr = CCTK_GroupData(gi, &group);
    assert(!ierr);
    if (group.grouptype != CCTK_GF)
      continue;

Insertion in AMReX/src/io.cxx at line 194 [20.107]

[17.500]

[21.162]

    cGroup group;
    int ierr = CCTK_GroupData(gi, &group);
    assert(!ierr);
    if (group.grouptype != CCTK_GF)
      continue;

Insertion in AMReX/src/io.cxx at line 240 [20.107]

[22.331]

    cGroup group;
    int ierr = CCTK_GroupData(gi, &group);
    assert(!ierr);
    if (group.grouptype != CCTK_GF)
      continue;

Insertion in AMReX/src/io.cxx at line 281 [20.107]

[3.108]

[14.2278]

    }
  }
}
void OutputScalars(const cGH *restrict cctkGH) {
  DECLARE_CCTK_ARGUMENTS;
  DECLARE_CCTK_PARAMETERS;
  static Timer timer("OutputNorms");
  Interval interval(timer);
  const bool is_root = CCTK_MyProc(nullptr) == 0;
  if (!is_root)
    return;
  const int numgroups = CCTK_NumGroups();
  for (int gi = 0; gi < numgroups; ++gi) {
    cGroup group;
    int ierr = CCTK_GroupData(gi, &group);
    assert(!ierr);
    if (group.grouptype != CCTK_SCALAR)
      continue;
    string groupname = unique_ptr<char>(CCTK_GroupName(gi)).get();
    groupname = regex_replace(groupname, regex("::"), "-");
    for (auto &c : groupname)
      c = tolower(c);
    ostringstream buf;
    buf << out_dir << "/" << groupname;
    const string filename = buf.str();
    ofstream file(filename.c_str(), std::ofstream::out | std::ofstream::app);
    // TODO: write proper header (once)
#if 0
    // Output header
    file << "#"
         << "\t"
         << "iteration"
         << "\t"
         << "time";
    for (int vi = 0; vi < groupdata.numvars; ++vi) {
      file << "\t" << CCTK_VarName(group.firstvarindex + vi);
    file << "\n";
#endif
    // Output data
    const GHExt::GlobalData &restrict globaldata = ghext->globaldata;
    const GHExt::GlobalData::ScalarGroupData &restrict scalargroupdata =
        globaldata.scalargroupdata.at(gi);
    const int tl = 0;
    file << cctkGH->cctk_iteration << "\t" << cctkGH->cctk_time;
    for (int vi = 0; vi < scalargroupdata.numvars; ++vi) {
      file << "\t" << *scalargroupdata.data.at(tl).at(vi);

Insertion in AMReX/src/io.cxx at line 337 [20.107]
[14.2284]
[17.1350]
```
    file << "\n";
```
Insertion in AMReX/src/io.cxx at line 357 [20.107]
[16.2101]
[16.2101]
```
    OutputScalars(cctkGH);
```

Insertion in AMReX/src/reduction.cxx at line 72 [22.870]

[23.8587]

[22.1763]

  cGroup group;
  int ierr = CCTK_GroupData(gi, &group);
  assert(!ierr);
  assert(group.grouptype == CCTK_GF);

Insertion in AMReX/src/schedule.cxx at line 299 [24.531]

[2.797]

  }
}
// Set grid scalars to nan
void poison_invalid(const GHExt::GlobalData::ScalarGroupData &scalargroupdata, int vi,
                    int tl) {
  DECLARE_CCTK_PARAMETERS;
  if (!poison_undefined_values)
    return;
  const valid_t &valid = scalargroupdata.valid.at(tl).at(vi);
  if (valid.valid_int && valid.valid_bnd)
    return;
  // scalar have no boundary so we expect it to alway be valid
  assert(valid.valid_bnd);
  if (!valid.valid_int) {
    CCTK_REAL *restrict const ptr = scalargroupdata.data.at(tl).at(vi);
    *ptr = 0.0 / 0.0;

Insertion in AMReX/src/schedule.cxx at line 320 [24.531]

[2.801]

[25.2001]

}
// Ensure grid scalars are not nan
void check_valid(const GHExt::GlobalData::ScalarGroupData &scalargroupdata, int vi,
                 int tl) {
  DECLARE_CCTK_PARAMETERS;
  if (!poison_undefined_values)
    return;
  const valid_t &valid = scalargroupdata.valid.at(tl).at(vi);
  if (!valid.valid_int && !valid.valid_bnd)
    return;
  // scalar have no boundary so we expect it to alway be valid
  assert(valid.valid_bnd);
  atomic<bool> found_nan{false};
  if (valid.valid_int) {
    CCTK_REAL *restrict const ptr = scalargroupdata.data.at(tl).at(vi);
    if (CCTK_BUILTIN_EXPECT(isnan(*ptr), false)) {
      found_nan = true;
    }
  }
  if (CCTK_BUILTIN_EXPECT(found_nan, false)) {
    const char *where = valid.valid_int && valid.valid_bnd
                            ? "interior and boundary"
                            : valid.valid_int ? "interior" : "boundary";
    CCTK_VERROR("Grid Scalar \"%s\" has nans on time level %d; expected valid %s",
                CCTK_FullVarName(scalargroupdata.firstvarindex + vi),
                tl, where);
  }

Insertion in AMReX/src/schedule.cxx at line 443 [24.531]

[26.1954]

[24.2241]


  // Grid scalar pointers
  {
    auto &restrict globaldata = ghext->globaldata;
    const int num_groups = CCTK_NumGroups();
    for (int gi = 0; gi < num_groups; ++gi) {
      cGroup group;
      int ierr = CCTK_GroupData(gi, &group);
      assert(!ierr);
      if (group.grouptype != CCTK_SCALAR)
        continue;
      auto &restrict scalargroupdata = globaldata.scalargroupdata.at(gi);
      for (int tl = 0; tl < int(scalargroupdata.data.size()); ++tl) {
        const auto &restrict vars = scalargroupdata.data.at(tl);
        for (int vi = 0; vi < scalargroupdata.numvars; ++vi) {
          cctkGH->data[scalargroupdata.firstvarindex + vi][tl] = vars.at(vi);
        }
      }
    }
  }

Insertion in AMReX/src/schedule.cxx at line 467 [24.531]

[24.2290]

  // Grid scalar pointers
  {
    auto &restrict globaldata = ghext->globaldata;
    const int num_groups = CCTK_NumGroups();
    for (int gi = 0; gi < num_groups; ++gi) {
      cGroup group;
      int ierr = CCTK_GroupData(gi, &group);
      assert(!ierr);
      if (group.grouptype != CCTK_SCALAR)
        continue;
      auto &restrict scalargroupdata = globaldata.scalargroupdata.at(gi);
      for (int tl = 0; tl < int(scalargroupdata.data.size()); ++tl) {
        for (int vi = 0; vi < scalargroupdata.numvars; ++vi) {
          cctkGH->data[scalargroupdata.firstvarindex + vi][tl] = nullptr;
        }
      }
    }
  }

Replacement in AMReX/src/schedule.cxx at line 550 [24.531]

∅:D[23.9879] → [12.2429:2486]

B:BD[13.792] → [12.2429:2486]

  for (auto &restrict groupdata : leveldata.groupdata) {

[13.760]

[17.1400]

  const int num_groups = CCTK_NumGroups();
  for (int gi = 0; gi < num_groups; ++gi) {
    cGroup group;
    int ierr = CCTK_GroupData(gi, &group);
    assert(!ierr);
    if (group.grouptype != CCTK_GF)
      continue;
    auto &restrict groupdata = leveldata.groupdata.at(gi);

Replacement in AMReX/src/schedule.cxx at line 612 [24.531]

∅:D[27.15968] → [12.2598:2655]

B:BD[13.1200] → [12.2598:2655]

  for (auto &restrict groupdata : leveldata.groupdata) {

[27.15968]

[17.1652]

  const int num_groups = CCTK_NumGroups();
  for (int gi = 0; gi < num_groups; ++gi) {
    cGroup group;
    int ierr = CCTK_GroupData(gi, &group);
    assert(!ierr);
    if (group.grouptype != CCTK_GF)
      continue;
    auto &restrict groupdata = leveldata.groupdata.at(gi);

Replacement in AMReX/src/schedule.cxx at line 991 [24.531]

B:BD[15.5662] → [15.5662:5721]

    for (auto &restrict groupdata : leveldata.groupdata) {

[15.5662]

[15.5721]

    const int num_groups = CCTK_NumGroups();
    for (int gi = 0; gi < num_groups; ++gi) {
      cGroup group;
      int ierr = CCTK_GroupData(gi, &group);
      assert(!ierr);
      if (group.grouptype != CCTK_GF)
        continue;
      auto &restrict groupdata = leveldata.groupdata.at(gi);

Insertion in AMReX/src/schedule.cxx at line 1014 [24.531]

[15.6167]

[19.4940]


  // invalidate scalars
  {
    auto &restrict globaldata = ghext->globaldata;
    const int num_groups = CCTK_NumGroups();
    for (int gi = 0; gi < num_groups; ++gi) {
      cGroup group;
      int ierr = CCTK_GroupData(gi, &group);
      assert(!ierr);
      if (group.grouptype != CCTK_SCALAR)
        continue;
      auto &restrict scalargroupdata = globaldata.scalargroupdata.at(gi);
      const bool checkpoint = get_group_checkpoint_flag(scalargroupdata.groupindex);
      if (!checkpoint) {
        // Invalidate all time levels
        const int ntls = scalargroupdata.data.size();
        for (int tl = 0; tl < ntls; ++tl) {
          for (int vi = 0; vi < scalargroupdata.numvars; ++vi) {
            scalargroupdata.valid.at(tl).at(vi) = valid_t();
            poison_invalid(scalargroupdata, vi, tl);
          }
        }
      }
    }
  }

Replacement in AMReX/src/schedule.cxx at line 1051 [24.531]

B:BD[12.2710] → [12.2710:2769]

    for (auto &restrict groupdata : leveldata.groupdata) {

[12.2710]

[17.1770]

    const int num_groups = CCTK_NumGroups();
    for (int gi = 0; gi < num_groups; ++gi) {
      cGroup group;
      int ierr = CCTK_GroupData(gi, &group);
      assert(!ierr);
      if (group.grouptype != CCTK_GF)
        continue;
      auto &restrict groupdata = leveldata.groupdata.at(gi);

Insertion in AMReX/src/schedule.cxx at line 1079 [24.531]

[17.2076]

[19.5302]


  // cycle scalars
  {
    auto &restrict globaldata = ghext->globaldata;
    const int num_groups = CCTK_NumGroups();
    for (int gi = 0; gi < num_groups; ++gi) {
      cGroup group;
      int ierr = CCTK_GroupData(gi, &group);
      assert(!ierr);
      if (group.grouptype != CCTK_SCALAR)
        continue;
      auto &restrict scalargroupdata = globaldata.scalargroupdata.at(gi);
      const int ntls = scalargroupdata.data.size();
      // Rotate time levels and invalidate current time level
      if (ntls > 1) {
        auto tmp = move(scalargroupdata.data.at(ntls - 1));
        for (int tl = ntls - 1; tl > 0; --tl)
          scalargroupdata.data.at(tl) = scalargroupdata.data.at(tl - 1);
        scalargroupdata.data.at(0) = tmp;
        vector<valid_t> vtmp = move(scalargroupdata.valid.at(ntls - 1));
        for (int tl = ntls - 1; tl > 0; --tl)
          scalargroupdata.valid.at(tl) = move(scalargroupdata.valid.at(tl - 1));
        scalargroupdata.valid.at(0) = move(vtmp);
        for (int vi = 0; vi < scalargroupdata.numvars; ++vi)
          scalargroupdata.valid.at(0).at(vi) = valid_t();
        for (int vi = 0; vi < scalargroupdata.numvars; ++vi)
          poison_invalid(scalargroupdata, vi, 0);
      }
    }
  }

Insertion in AMReX/src/schedule.cxx at line 1222 [24.531]
[18.9041]
[28.510]
```
    const auto &restrict globaldata = ghext->globaldata;
```

Replacement in AMReX/src/schedule.cxx at line 1228 [24.531]

B:BD[18.9332] → [18.9332:9404]

        const auto &restrict groupdata = leveldata.groupdata.at(rd.gi);

[18.9332]

[18.9404]

        cGroup group;
        int ierr = CCTK_GroupData(rd.gi, &group);
        assert(!ierr);
        // TODO: something about this cast (https://stackoverflow.com/questions/6179314/casting-pointers-and-the-ternary-operator-have-i-reinvented-the-wheel)
        const GHExt::CommonGroupData * groupdata = group.grouptype == CCTK_GF ?
          static_cast<const GHExt::CommonGroupData *>(&leveldata.groupdata.at(rd.gi)) :
          static_cast<const GHExt::CommonGroupData *>(&globaldata.scalargroupdata.at(rd.gi));

Replacement in AMReX/src/schedule.cxx at line 1237 [24.531]

B:BD[18.9444] → [18.9444:9511]

        const valid_t &have = groupdata.valid.at(rd.tl).at(rd.vi);

[18.9444]

[18.9511]

        const valid_t &have = groupdata->valid.at(rd.tl).at(rd.vi);

Replacement in AMReX/src/schedule.cxx at line 1247 [24.531]

B:BD[18.9980] → [18.9980:10045]

              CCTK_FullVarName(groupdata.firstvarindex + rd.vi),

[18.9980]

[18.10045]

              CCTK_FullVarName(groupdata->firstvarindex + rd.vi),

Replacement in AMReX/src/schedule.cxx at line 1250 [24.531]

B:BD[18.10147] → [25.2024:2081]

        check_valid(leveldata, groupdata, rd.vi, rd.tl);

[18.10147]

        if (group.grouptype == CCTK_GF) {
          check_valid(leveldata, leveldata.groupdata.at(rd.gi), rd.vi, rd.tl);
        } else {
          check_valid(globaldata.scalargroupdata.at(rd.gi), rd.vi, rd.tl);
        }

Insertion in AMReX/src/schedule.cxx at line 1261 [24.531]
[15.6694]
[15.6694]
```
    auto &restrict globaldata = ghext->globaldata;
```

Replacement in AMReX/src/schedule.cxx at line 1283 [24.531]

B:BD[15.7476] → [15.7476:7603]

        auto &restrict groupdata = leveldata.groupdata.at(wr.gi);
        valid_t &have = groupdata.valid.at(wr.tl).at(wr.vi);

[15.7476]

[15.7603]

        cGroup group;
        int ierr = CCTK_GroupData(wr.gi, &group);
        assert(!ierr);
        GHExt::CommonGroupData * groupdata = group.grouptype == CCTK_GF ?
          static_cast<GHExt::CommonGroupData *>(&leveldata.groupdata.at(wr.gi)) :
          static_cast<GHExt::CommonGroupData *>(&globaldata.scalargroupdata.at(wr.gi));
        valid_t &have = groupdata->valid.at(wr.tl).at(wr.vi);

Insertion in AMReX/src/schedule.cxx at line 1341 [24.531]
[18.10435]
[28.563]
```
    auto &restrict globaldata = ghext->globaldata;
```

Replacement in AMReX/src/schedule.cxx at line 1347 [24.531]

B:BD[18.10724] → [18.10724:10790]

        auto &restrict groupdata = leveldata.groupdata.at(wr.gi);

[18.10724]

[18.10790]

        cGroup group;
        int ierr = CCTK_GroupData(wr.gi, &group);
        assert(!ierr);
        GHExt::CommonGroupData * groupdata = group.grouptype == CCTK_GF ?
          static_cast<GHExt::CommonGroupData *>(&leveldata.groupdata.at(wr.gi)) :
          static_cast<GHExt::CommonGroupData *>(&globaldata.scalargroupdata.at(wr.gi));

Replacement in AMReX/src/schedule.cxx at line 1355 [24.531]

B:BD[18.10834] → [18.10834:10895]

        valid_t &have = groupdata.valid.at(wr.tl).at(wr.vi);

[18.10834]

[18.10895]

        valid_t &have = groupdata->valid.at(wr.tl).at(wr.vi);

Replacement in AMReX/src/schedule.cxx at line 1359 [24.531]

B:BD[18.11024] → [25.2082:2139]

        check_valid(leveldata, groupdata, wr.vi, wr.tl);

[18.11024]

        if (group.grouptype == CCTK_GF) {
          check_valid(leveldata, leveldata.groupdata.at(wr.gi), wr.vi, wr.tl);
        } else {
          check_valid(globaldata.scalargroupdata.at(wr.gi), wr.vi, wr.tl);
        }

Insertion in AMReX/src/schedule.cxx at line 1400 [24.531]

[14.4154]

[12.3153]

      cGroup group;
      int ierr = CCTK_GroupData(gi, &group);
      assert(!ierr);
      if (group.grouptype != CCTK_GF)
        continue;

Insertion in AMReX/src/schedule.cxx at line 1510 [24.531]

[4.3525]

    cGroup group;
    int ierr = CCTK_GroupData(gi, &group);
    assert(!ierr);
    if (group.grouptype != CCTK_GF)
      continue;

Insertion in AMReX/src/schedule.cxx at line 1582 [24.531]

[29.631]

[10.8625]

    cGroup group;
    int ierr = CCTK_GroupData(gi, &group);
    assert(!ierr);
    if (group.grouptype != CCTK_GF)
      continue;

Insertion in AMReX/src/schedule.hxx at line 77 [24.5992]

[27.18749]

[24.6506]

void poison_invalid(const GHExt::GlobalData::ScalarGroupData &groupdata, int vi,
                    int tl);
void check_valid(const GHExt::GlobalData::ScalarGroupData &groupdata, int vi,
                 int tl);