Merge pull request #39 from AFD-Illinois/feature/parthenon-bump

Bump Parthenon & Fix SMR

.gitignore

@@ -1,3 +1,7 @@
+# Memory dumps
+core*
+*.swp
+
 # Various script results/logs
 out-*.txt
 *.json
@@ -6,7 +10,6 @@ convergence.txt
 *.png
 *.mp4
 *.webm
-core.*
 frames_*/
 logs/
 *.log
@@ -77,6 +80,3 @@ make_args
 # Python files
 __pycache__/
 *.pyc
-
-# added by Hyerin
-*.swp

CMakeLists.txt

@@ -22,6 +22,9 @@ set(PARTHENON_LINT_DEFAULT OFF CACHE BOOL "KHARMA Override")
 set(PARTHENON_DISABLE_HDF5_COMPRESSION OFF CACHE BOOL "KHARMA Override")
 # Don't build sparse (selectively-allocated) variable support
 set(PARTHENON_DISABLE_SPARSE ON CACHE BOOL "KHARMA Override")
+# Set to move MPI buffers to host; slower but less crashy
+# Favor setting this per-machine in machines/
+set(PARTHENON_ENABLE_HOST_COMM_BUFFERS OFF CACHE BOOL "KHARMA Override")
 
 # Parthenon internal build options
 set(BUILD_TESTING OFF CACHE BOOL "KHARMA Override")

external/kokkos-kernels/KokkosBatched_Dot_Internal.hpp

@@ -5,6 +5,8 @@
 
 #include "KokkosBatched_Util.hpp"
 
+#define KOKKOS_IMPL_DO_NOT_USE_PRINTF(...) ::printf(__VA_ARGS__)
+
 namespace KokkosBatched {
 
 ///

external/patches/parthenon-use-gr-coordinates.patch

@@ -1,5 +1,5 @@
 diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt
-index 45566b0b..a9abdc1c 100644
+index aaeabd8a..e354ef3d 100644
 --- a/src/CMakeLists.txt
 +++ b/src/CMakeLists.txt
 @@ -90,7 +90,7 @@ set(COMPILED_WITH ${CMAKE_CXX_COMPILER})
@@ -11,13 +11,13 @@ index 45566b0b..a9abdc1c 100644
 
  configure_file(config.hpp.in generated/config.hpp @ONLY)
 
-@@ -309,6 +309,8 @@ lint_target(parthenon)
+@@ -331,6 +331,8 @@ lint_target(parthenon)
  target_include_directories(parthenon PUBLIC
    $<BUILD_INTERFACE:${CMAKE_CURRENT_SOURCE_DIR}>
    $<BUILD_INTERFACE:${CMAKE_CURRENT_BINARY_DIR}/generated>
 +  $<BUILD_INTERFACE:${CMAKE_CURRENT_SOURCE_DIR}/../../../kharma>
 +  $<BUILD_INTERFACE:${CMAKE_CURRENT_SOURCE_DIR}/../../variant/include>
-   $<INSTALL_INTERFACE:${CMAKE_INSTALL_INCLUDEDIR}/parthenon>
+   $<INSTALL_INTERFACE:${CMAKE_INSTALL_INCLUDEDIR}>
    )
 
 diff --git a/src/coordinates/coordinates.hpp b/src/coordinates/coordinates.hpp
@@ -32,4 +32,3 @@ index d1290dee..50bfc840 100644
 
  namespace parthenon {
 
-

kharma/b_cleanup/b_cleanup.cpp

@@ -62,9 +62,6 @@ void B_Cleanup::CleanupDivergence(std::shared_ptr<MeshData<Real>>& md) {}
 using namespace parthenon;
 using namespace parthenon::solvers;
 
-// TODO get the transport manager working later
-// Needs a call every X steps option, probably return a TaskList or TaskRegion
-
 std::shared_ptr<KHARMAPackage> B_Cleanup::Initialize(ParameterInput *pin, std::shared_ptr<Packages_t>& packages)
 {
     auto pkg = std::make_shared<KHARMAPackage>("B_Cleanup");

kharma/b_cleanup/bicgstab_solver.hpp

@@ -22,10 +22,11 @@
 #include "interface/meshblock_data.hpp"
 #include "interface/state_descriptor.hpp"
 #include "kokkos_abstraction.hpp"
-#include "solvers/solver_utils.hpp"
 #include "tasks/task_id.hpp"
 #include "tasks/task_list.hpp"
 
+#include "solver_utils.hpp"
+
 namespace parthenon {
 
 namespace solvers {

kharma/b_cleanup/solver_utils.hpp

@@ -0,0 +1,176 @@
+//========================================================================================
+// (C) (or copyright) 2021. Triad National Security, LLC. All rights reserved.
+//
+// This program was produced under U.S. Government contract 89233218CNA000001 for Los
+// Alamos National Laboratory (LANL), which is operated by Triad National Security, LLC
+// for the U.S. Department of Energy/National Nuclear Security Administration. All rights
+// in the program are reserved by Triad National Security, LLC, and the U.S. Department
+// of Energy/National Nuclear Security Administration. The Government is granted for
+// itself and others acting on its behalf a nonexclusive, paid-up, irrevocable worldwide
+// license in this material to reproduce, prepare derivative works, distribute copies to
+// the public, perform publicly and display publicly, and to permit others to do so.
+//========================================================================================
+#ifndef SOLVERS_SOLVER_UTILS_HPP_
+#define SOLVERS_SOLVER_UTILS_HPP_
+
+#include <string>
+#include <vector>
+
+#include "basic_types.hpp"
+#include "kokkos_abstraction.hpp"
+
+namespace parthenon {
+
+namespace solvers {
+
+struct SparseMatrixAccessor {
+  ParArray1D<int> ioff, joff, koff;
+  ParArray1D<int> ioff_inv, joff_inv, koff_inv;
+  ParArray1D<int> inv_entries;
+
+  const int nstencil;
+  int ndiag;
+  SparseMatrixAccessor() : nstencil(0), ndiag(0) {}
+  SparseMatrixAccessor(const SparseMatrixAccessor &sp)
+      : ioff(sp.ioff), joff(sp.joff), koff(sp.koff), nstencil(sp.nstencil),
+        ndiag(sp.ndiag), inv_entries(sp.inv_entries) {}
+  SparseMatrixAccessor(const std::string &label, const int n,
+                       std::vector<std::vector<int>> off)
+      : ioff(label + "_ioff", n), joff(label + "_joff", n), koff(label + "_koff", n),
+        ioff_inv(label + "_ioff_inv", n), joff_inv(label + "_joff_inv", n),
+        koff_inv(label + "_koff_inv", n), inv_entries(label + "_inv_ent", n),
+        nstencil(n) {
+    PARTHENON_REQUIRE_THROWS(off.size() == 3,
+                             "Offset array must have dimensions off[3][*]");
+    PARTHENON_REQUIRE_THROWS(off[0].size() >= n, "Offset array off[0][*] too small");
+    PARTHENON_REQUIRE_THROWS(off[1].size() >= n, "Offset array off[1][*] too small");
+    PARTHENON_REQUIRE_THROWS(off[2].size() >= n, "Offset array off[2][*] too small");
+    auto ioff_h = Kokkos::create_mirror_view(Kokkos::HostSpace(), ioff);
+    auto joff_h = Kokkos::create_mirror_view(Kokkos::HostSpace(), joff);
+    auto koff_h = Kokkos::create_mirror_view(Kokkos::HostSpace(), koff);
+
+    auto ioff_inv_h = Kokkos::create_mirror_view(Kokkos::HostSpace(), ioff_inv);
+    auto joff_inv_h = Kokkos::create_mirror_view(Kokkos::HostSpace(), joff_inv);
+    auto koff_inv_h = Kokkos::create_mirror_view(Kokkos::HostSpace(), koff_inv);
+
+    auto inv_ent_h = Kokkos::create_mirror_view(Kokkos::HostSpace(), inv_entries);
+
+    for (int i = 0; i < n; i++) {
+      ioff_h(i) = off[0][i];
+      joff_h(i) = off[1][i];
+      koff_h(i) = off[2][i];
+      // this is inverse.
+      ioff_inv_h(i) = -off[0][i];
+      joff_inv_h(i) = -off[1][i];
+      koff_inv_h(i) = -off[2][i];
+
+      if (off[0][i] == 0 && off[1][i] == 0 && off[2][i] == 0) {
+        ndiag = i;
+      }
+    }
+    for (int i = 0; i < n; i++) {
+      for (int j = 0; j < n; j++) {
+        if (ioff_h(i) == ioff_inv_h(j) && joff_h(i) == joff_inv_h(j) &&
+            koff_h(i) == koff_inv_h(j)) {
+          inv_entries(i) = j;
+          std::cout << "inv_entries:" << i << " " << j << std::endl;
+        }
+      } // j
+    }   // i
+
+    Kokkos::deep_copy(ioff, ioff_h);
+    Kokkos::deep_copy(joff, joff_h);
+    Kokkos::deep_copy(koff, koff_h);
+
+    Kokkos::deep_copy(inv_entries, inv_ent_h);
+  }
+
+  template <typename PackType>
+  KOKKOS_INLINE_FUNCTION Real MatVec(const PackType &spmat, const int imat_lo,
+                                     const int imat_hi, const PackType &v, const int iv,
+                                     const int b, const int k, const int j,
+                                     const int i) const {
+    Real matvec = 0.0;
+    for (int n = imat_lo; n <= imat_hi; n++) {
+      const int m = n - imat_lo;
+      matvec += spmat(b, n, k, j, i) * v(b, iv, k + koff(m), j + joff(m), i + ioff(m));
+    }
+    return matvec;
+  }
+
+  template <typename PackType>
+  KOKKOS_INLINE_FUNCTION Real Jacobi(const PackType &spmat, const int imat_lo,
+                                     const int imat_hi, const PackType &v, const int iv,
+                                     const int b, const int k, const int j, const int i,
+                                     const Real rhs) const {
+    const Real matvec = MatVec(spmat, imat_lo, imat_hi, v, iv, b, k, j, i);
+    return (rhs - matvec + spmat(b, imat_lo + ndiag, k, j, i) * v(b, iv, k, j, i)) /
+           spmat(b, imat_lo + ndiag, k, j, i);
+  }
+};
+
+template <typename T>
+struct Stencil {
+  ParArray1D<T> w;
+  ParArray1D<int> ioff, joff, koff;
+  const int nstencil;
+  int ndiag;
+  Stencil() : nstencil(0), ndiag(0) {}
+  Stencil(const Stencil<T> &st)
+      : w(st.w), ioff(st.ioff), joff(st.joff), koff(st.koff), nstencil(st.nstencil),
+        ndiag(st.ndiag) {}
+  Stencil(const std::string &label, const int n, std::vector<T> wgt,
+          std::vector<std::vector<int>> off)
+      : w(label + "_w", n), ioff(label + "_ioff", n), joff(label + "_joff", n),
+        koff(label + "_koff", n), nstencil(n) {
+    PARTHENON_REQUIRE_THROWS(off.size() == 3,
+                             "Offset array must have dimensions off[3][*]");
+    PARTHENON_REQUIRE_THROWS(wgt.size() >= n, "Weight array too small");
+    PARTHENON_REQUIRE_THROWS(off[0].size() >= n, "Offset array off[0][*] too small");
+    PARTHENON_REQUIRE_THROWS(off[1].size() >= n, "Offset array off[1][*] too small");
+    PARTHENON_REQUIRE_THROWS(off[2].size() >= n, "Offset array off[2][*] too small");
+    auto w_h = Kokkos::create_mirror_view(Kokkos::HostSpace(), w);
+    auto ioff_h = Kokkos::create_mirror_view(Kokkos::HostSpace(), ioff);
+    auto joff_h = Kokkos::create_mirror_view(Kokkos::HostSpace(), joff);
+    auto koff_h = Kokkos::create_mirror_view(Kokkos::HostSpace(), koff);
+
+    for (int i = 0; i < n; i++) {
+      w_h(i) = wgt[i];
+      ioff_h(i) = off[0][i];
+      joff_h(i) = off[1][i];
+      koff_h(i) = off[2][i];
+      if (off[0][i] == 0 && off[1][i] == 0 && off[2][i] == 0) {
+        ndiag = i;
+      }
+    }
+
+    Kokkos::deep_copy(w, w_h);
+    Kokkos::deep_copy(ioff, ioff_h);
+    Kokkos::deep_copy(joff, joff_h);
+    Kokkos::deep_copy(koff, koff_h);
+  }
+
+  template <typename PackType>
+  KOKKOS_INLINE_FUNCTION Real MatVec(const PackType &v, const int iv, const int b,
+                                     const int k, const int j, const int i) const {
+    Real matvec = 0.0;
+    for (int n = 0; n < nstencil; n++) {
+      matvec += w(n) * v(b, iv, k + koff(n), j + joff(n), i + ioff(n));
+    }
+    return matvec;
+  }
+
+  template <typename PackType>
+  KOKKOS_INLINE_FUNCTION Real Jacobi(const PackType &v, const int iv, const int b,
+                                     const int k, const int j, const int i,
+                                     const Real rhs) const {
+    const Real matvec = MatVec(v, iv, b, k, j, i);
+    return (rhs - matvec + w(ndiag) * v(b, iv, k, j, i)) / w(ndiag);
+  }
+};
+
+} // namespace solvers
+
+} // namespace parthenon
+
+#endif // SOLVERS_SOLVER_UTILS_HPP_

kharma/b_cleanup/solvers.cpp

@@ -0,0 +1,22 @@
+//========================================================================================
+// (C) (or copyright) 2021. Triad National Security, LLC. All rights reserved.
+//
+// This program was produced under U.S. Government contract 89233218CNA000001 for Los
+// Alamos National Laboratory (LANL), which is operated by Triad National Security, LLC
+// for the U.S. Department of Energy/National Nuclear Security Administration. All rights
+// in the program are reserved by Triad National Security, LLC, and the U.S. Department
+// of Energy/National Nuclear Security Administration. The Government is granted for
+// itself and others acting on its behalf a nonexclusive, paid-up, irrevocable worldwide
+// license in this material to reproduce, prepare derivative works, distribute copies to
+// the public, perform publicly and display publicly, and to permit others to do so.
+//========================================================================================
+
+#include "bicgstab_solver.hpp"
+
+namespace parthenon {
+namespace solvers {
+
+int BiCGStabCounter::global_num_bicgstab_solvers = 0;
+
+} // namespace solvers
+} // namespace parthenon