Merge pull request #556 from lanl/jdolence/iterative_tasks

Iterative Tasks

Author: jdolence

CHANGELOG.md

@@ -3,6 +3,7 @@
 ## Current develop
 
 ### Added (new features/APIs/variables/...)
+- [[PR 556]](https://github.com/lanl/parthenon/pull/556) Introduce iterative tasks and regionally dependent tasks
 - [[PR 578]](https://github.com/lanl/parthenon/pull/578) Add some profiling regions to tasks in particles example
 - [[PR 577]](https://github.com/lanl/parthenon/pull/577) Update invalid indices to allow for no-op loops
 - [[PR 564]](https://github.com/lanl/parthenon/pull/564) Add EstimateTimestep to particles example task list

CMakeLists.txt

@@ -67,9 +67,9 @@ include(cmake/Format.cmake)
 include(cmake/Lint.cmake)
 
 # regression test reference data
-set(REGRESSION_GOLD_STANDARD_VER 10 CACHE STRING "Version of gold standard to download and use")
+set(REGRESSION_GOLD_STANDARD_VER 11 CACHE STRING "Version of gold standard to download and use")
 set(REGRESSION_GOLD_STANDARD_HASH
-  "SHA512=59c0eec495bab38629edae9c99e1da51d68f7a2815087638d0b6b2fb7558e5d34171caae32b552fa52c10abe0fb3b720723c3b2a12645f8e144e7d134e7fe98b"
+  "SHA512=11fdcb9e42457f41f2f4d4cce4636a73cdab601eac9638fb63f9284c10dc90f5719b190e075a079c13804d8126e56a80590ffc71c6cc5eeb0d0afa72cbeaa9a1"
   CACHE STRING "Hash of default gold standard file to download")
 option(REGRESSION_GOLD_STANDARD_SYNC "Automatically sync gold standard files." ON)
 

docs/tasks.md

@@ -1,7 +1,7 @@
 # Tasks
 
 ## TaskList
-The `TaskList` class implements methods to build and execute a set of tasks with associated dependencies.  The main functionality of the class is implemented in two member functions:
+The `TaskList` class implements methods to build and execute a set of tasks with associated dependencies.  The class implements a few public facing member functions that provide useful functionality for downstream apps:
 
 ### AddTask
 `AddTask` is a templated variadic function that takes the task function to be executed, the task dependencies (see `TaskID` below), and the arguments to the task function as it's arguments.  All arguments are captured by value in a lambda for later execution.
@@ -10,14 +10,17 @@ When adding functions that are non-static class member functions, a slightly dif
 
 Examples of both `AddTask` calls can be found in the advection example [here](../example/advection/advection_driver.cpp).
 
+### AddIteration
+`AddIteration` provides a means of grouping a set of tasks together that will be executed repeatedly until stopping criteria are satisfied.  `AddIteration` returns an `IterativeTasks` object which provides overloaded `AddTask` functions as described above, but internally handles the bookkeeping necessary to maintain the association of all the tasks associated with the iterative process.  A special function `SetCompletionTask`, which behaves identically to `AddTask`, allows a task to be defined that evaluates the stopping criteria.  The maximum number of iterations can be controlled through the `SetMaxIterations` member function and the number of iterations between evaluating the stopping criteria can be set with the `SetCheckInterval` function.
+
 ### DoAvailable
-`DoAvailable` loops over the task list once, executing all tasks whose dependencies are satisfied.  The function returns either `TaskListStatus::complete` if all tasks have been executed (and the task list is therefore empty) or `TaskListStatus::running` if tasks remain to be completed.
+`DoAvailable` loops over the task list once, executing all tasks whose dependencies are satisfied.  Completed tasks are removed from the task list.
 
 ## TaskID
 The `TaskID` class implements methods that allow Parthenon to keep track of tasks, their dependencies, and what remains to be completed.  The main way application code will interact with this object is as a returned object from `TaskList::AddTask` and as an argument to subsequent calls to `TaskList::AddTask` as a dependency for other tasks.  When used as a dependency, `TaskID` objects can be combined with the bitwise or operator (`|`) to specify multiple dependencies.
 
  ## TaskRegion
- `TaskRegion` is defined via a simple using statement as a `std::vector<TaskList>`.  During task execution (described below), all task lists in a `TaskRegion` can be operated on concurrently.  For example, a `TaskRegion` can be used to construct independent task lists for each `MeshBlock`.
+ `TaskRegion` is a lightweight class that wraps `std::vector<TaskList>`, providing a little extra functionality.  During task execution (described below), all task lists in a `TaskRegion` can be operated on concurrently.  For example, a `TaskRegion` can be used to construct independent task lists for each `MeshBlock`.  Occasionally, it is useful to have a task not be considered complete until that task completes in all lists of a region.  For example, a global iterative solver cannot be considered complete until the stopping criteria are satisfied everywhere, which may require evaluating those criteria in tasks that live in different lists within a region.  An example of this use case is shown[here](../example/poisson/poisson_driver.cpp).
 
  ## TaskCollection
  A `TaskCollection` contains a `std::vector<TaskRegion>`, i.e. an ordered list of `TaskRegion`s.  Importantly, each `TaskRegion` will be executed to completion before subsequent `TaskRegion`s, introducing a notion of sequential execution and enabling flexibility in task granularity.  For example, the following code fragment uses the `TaskCollection` and `TaskRegion` abstractions to express work that can be done asynchronously across blocks, followed by a bulk synchronous task involving all blocks, and finally another round of asynchronous work.

example/CMakeLists.txt

@@ -19,3 +19,4 @@ add_subdirectory(kokkos_pi)
 add_subdirectory(particles)
 add_subdirectory(particle_leapfrog)
 add_subdirectory(particle_tracers)
+add_subdirectory(poisson)

example/poisson/CMakeLists.txt

@@ -0,0 +1,27 @@
+#=========================================================================================
+# (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.
+#=========================================================================================
+
+get_property(DRIVER_LIST GLOBAL PROPERTY DRIVERS_USED_IN_TESTS)
+if( "poisson-example" IN_LIST DRIVER_LIST OR NOT PARTHENON_DISABLE_EXAMPLES)
+  add_executable(
+      poisson-example
+          poisson_driver.cpp
+          poisson_driver.hpp
+          poisson_package.cpp
+          poisson_package.hpp
+          main.cpp
+          parthenon_app_inputs.cpp
+  )
+  target_link_libraries(poisson-example PRIVATE Parthenon::parthenon)
+  lint_target(poisson-example)
+endif()

example/poisson/main.cpp

@@ -0,0 +1,53 @@
+//========================================================================================
+// (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 "parthenon_manager.hpp"
+
+#include "poisson_driver.hpp"
+
+int main(int argc, char *argv[]) {
+  using parthenon::ParthenonManager;
+  using parthenon::ParthenonStatus;
+  ParthenonManager pman;
+
+  // Redefine parthenon defaults
+  pman.app_input->ProcessPackages = poisson_example::ProcessPackages;
+  pman.app_input->ProblemGenerator = poisson_example::ProblemGenerator;
+
+  // call ParthenonInit to initialize MPI and Kokkos, parse the input deck, and set up
+  auto manager_status = pman.ParthenonInit(argc, argv);
+  if (manager_status == ParthenonStatus::complete) {
+    pman.ParthenonFinalize();
+    return 0;
+  }
+  if (manager_status == ParthenonStatus::error) {
+    pman.ParthenonFinalize();
+    return 1;
+  }
+  // Now that ParthenonInit has been called and setup succeeded, the code can now
+  // make use of MPI and Kokkos
+
+  // Initialize the driver
+  poisson_example::PoissonDriver driver(pman.pinput.get(), pman.app_input.get(),
+                                        pman.pmesh.get());
+
+  // This line actually runs the simulation
+  auto driver_status = driver.Execute();
+
+  // call MPI_Finalize and Kokkos::finalize if necessary
+  pman.ParthenonFinalize();
+
+  // MPI and Kokkos can no longer be used
+
+  return (0);
+}

example/poisson/parthenon_app_inputs.cpp

@@ -0,0 +1,75 @@
+//========================================================================================
+// (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 <sstream>
+#include <string>
+
+#include <parthenon/package.hpp>
+
+#include "config.hpp"
+#include "defs.hpp"
+#include "poisson_package.hpp"
+#include "utils/error_checking.hpp"
+
+using namespace parthenon::package::prelude;
+using namespace parthenon;
+
+// *************************************************//
+// redefine some weakly linked parthenon functions *//
+// *************************************************//
+
+namespace poisson_example {
+
+void ProblemGenerator(MeshBlock *pmb, ParameterInput *pin) {
+  auto &data = pmb->meshblock_data.Get();
+
+  Real x0 = pin->GetOrAddReal("poisson", "x0", 0.0);
+  Real y0 = pin->GetOrAddReal("poisson", "y0", 0.0);
+  Real z0 = pin->GetOrAddReal("poisson", "z0", 0.0);
+  Real radius = pin->GetOrAddReal("poisson", "radius", 0.1);
+
+  auto cellbounds = pmb->cellbounds;
+  IndexRange ib = cellbounds.GetBoundsI(IndexDomain::entire);
+  IndexRange jb = cellbounds.GetBoundsJ(IndexDomain::entire);
+  IndexRange kb = cellbounds.GetBoundsK(IndexDomain::entire);
+
+  auto coords = pmb->coords;
+  PackIndexMap imap;
+  const std::vector<std::string> vars({"density", "potential"});
+  const auto &q = data->PackVariables(vars, imap);
+  const int irho = imap["density"].first;
+  const int iphi = imap["potential"].first;
+
+  pmb->par_for(
+      "Poisson::ProblemGenerator", kb.s, kb.e, jb.s, jb.e, ib.s, ib.e,
+      KOKKOS_LAMBDA(const int k, const int j, const int i) {
+        Real dist2 = std::pow(coords.x1v(i) - x0, 2) + std::pow(coords.x2v(j) - y0, 2) +
+                     std::pow(coords.x3v(k) - z0, 2);
+        if (dist2 < radius * radius) {
+          q(irho, k, j, i) = 1.0 / (4.0 / 3.0 * M_PI * std::pow(radius, 3));
+        } else {
+          q(irho, k, j, i) = 0.0;
+        }
+        q(iphi, k, j, i) = 0.0;
+      });
+}
+
+Packages_t ProcessPackages(std::unique_ptr<ParameterInput> &pin) {
+  Packages_t packages;
+  auto pkg = poisson_package::Initialize(pin.get());
+  packages.Add(pkg);
+
+  return packages;
+}
+
+} // namespace poisson_example

example/poisson/parthinput.poisson

@@ -0,0 +1,62 @@
+# ========================================================================================
+#  Athena++ astrophysical MHD code
+#  Copyright(C) 2014 James M. Stone <[email protected]> and other code contributors
+#  Licensed under the 3-clause BSD License, see LICENSE file for details
+# ========================================================================================
+#  (C) (or copyright) 21. 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.
+# ========================================================================================
+
+<parthenon/job>
+problem_id = poisson
+
+<parthenon/mesh>
+
+nx1 = 64
+x1min = -1.0
+x1max = 1.0
+ix1_bc = outflow
+ox1_bc = outflow
+
+nx2 = 64
+x2min = -1.0
+x2max = 1.0
+ix2_bc = outflow
+ox2_bc = outflow
+
+nx3 = 1
+x3min = -0.5
+x3max = 0.5
+ix3_bc = periodic
+ox3_bc = periodic
+
+<parthenon/meshblock>
+nx1 = 16
+nx2 = 16
+nx3 = 1
+
+<parthenon/time>
+#nlim = -1
+#tlim = 1.0
+#integrator = rk2
+#ncycle_out_mesh = -10000
+
+<parthenon/output0>
+file_type = hdf5
+dt = 0.05
+variables = potential, density
+
+<poisson>
+max_iterations = 100000
+check_interval = 100
+error_tolerance = 1.e-8
+fail_without_convergence = false
+warn_without_convergence = true

example/poisson/poisson_driver.cpp

@@ -0,0 +1,104 @@
+//========================================================================================
+// (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 <memory>
+#include <string>
+#include <vector>
+
+// Local Includes
+#include "bvals/cc/bvals_cc_in_one.hpp"
+#include "interface/metadata.hpp"
+#include "interface/update.hpp"
+#include "mesh/meshblock_pack.hpp"
+#include "mesh/refinement_cc_in_one.hpp"
+#include "parthenon/driver.hpp"
+#include "poisson_driver.hpp"
+#include "poisson_package.hpp"
+#include "refinement/refinement.hpp"
+
+using namespace parthenon::driver::prelude;
+
+namespace poisson_example {
+
+parthenon::DriverStatus PoissonDriver::Execute() {
+  pouts->MakeOutputs(pmesh, pinput);
+  ConstructAndExecuteTaskLists<>(this);
+  pouts->MakeOutputs(pmesh, pinput);
+  return DriverStatus::complete;
+}
+
+TaskCollection PoissonDriver::MakeTaskCollection(BlockList_t &blocks) {
+  using namespace parthenon::Update;
+  TaskCollection tc;
+  TaskID none(0);
+
+  for (int i = 0; i < blocks.size(); i++) {
+    auto &pmb = blocks[i];
+    auto &base = pmb->meshblock_data.Get();
+    pmb->meshblock_data.Add("delta", base);
+  }
+
+  int max_iters = pmesh->packages.Get("poisson_package")->Param<int>("max_iterations");
+  int check_interval =
+      pmesh->packages.Get("poisson_package")->Param<int>("check_interval");
+  bool fail_flag =
+      pmesh->packages.Get("poisson_package")->Param<bool>("fail_without_convergence");
+  bool warn_flag =
+      pmesh->packages.Get("poisson_package")->Param<bool>("warn_without_convergence");
+
+  const int num_partitions = pmesh->DefaultNumPartitions();
+  TaskRegion &solver_region = tc.AddRegion(num_partitions);
+
+  for (int i = 0; i < num_partitions; i++) {
+    // make/get a mesh_data container for the state
+    auto &md = pmesh->mesh_data.GetOrAdd("base", i);
+    auto &mdelta = pmesh->mesh_data.GetOrAdd("delta", i);
+
+    TaskList &tl = solver_region[i];
+
+    auto &solver = tl.AddIteration("poisson solver");
+    solver.SetMaxIterations(max_iters);
+    solver.SetCheckInterval(check_interval);
+    solver.SetFailWithMaxIterations(fail_flag);
+    solver.SetWarnWithMaxIterations(warn_flag);
+    auto start_recv = solver.AddTask(none, &MeshData<Real>::StartReceiving, md.get(),
+                                     BoundaryCommSubset::all);
+
+    auto update = solver.AddTask(none, poisson_package::UpdatePhi<MeshData<Real>>,
+                                 md.get(), mdelta.get());
+
+    auto send =
+        solver.AddTask(update, parthenon::cell_centered_bvars::SendBoundaryBuffers, md);
+
+    auto recv = solver.AddTask(
+        start_recv, parthenon::cell_centered_bvars::ReceiveBoundaryBuffers, md);
+
+    auto setb =
+        solver.AddTask(recv | update, parthenon::cell_centered_bvars::SetBoundaries, md);
+
+    auto clear = solver.AddTask(send | setb, &MeshData<Real>::ClearBoundary, md.get(),
+                                BoundaryCommSubset::all);
+
+    auto check = solver.SetCompletionTask(
+        clear, poisson_package::CheckConvergence<MeshData<Real>>, md.get(), mdelta.get());
+    // mark task so that dependent tasks (below) won't execute
+    // until all task lists have completed it
+    solver_region.AddRegionalDependencies(0, i, check);
+
+    auto print = tl.AddTask(check, poisson_package::PrintComplete);
+  }
+
+  return tc;
+}
+
+} // namespace poisson_example