From e14716ce79334e4e66a82d5ba6d107882cb8a642 Mon Sep 17 00:00:00 2001 From: Ryan Moulday Date: Tue, 11 Jun 2024 16:10:00 -0400 Subject: [PATCH] Adding taskflow w/o coloring --- include/deal.II/base/work_stream.h | 234 ++++++++++++++++++++++++----- tests/mpi/mesh_worker_01.cc | 6 +- tests/mpi/mesh_worker_02.cc | 6 +- tests/mpi/mesh_worker_04.cc | 6 +- 4 files changed, 207 insertions(+), 45 deletions(-) diff --git a/include/deal.II/base/work_stream.h b/include/deal.II/base/work_stream.h index a8b0b2bb1c..6379772321 100644 --- a/include/deal.II/base/work_stream.h +++ b/include/deal.II/base/work_stream.h @@ -164,6 +164,44 @@ namespace WorkStream */ namespace internal { + /** + * A structure that contains a pointer to a scratch data object + * along with a flag that indicates whether this object is currently + * in use. + */ + template + struct ScratchDataObject + { + std::unique_ptr scratch_data; + bool currently_in_use; + + /** + * Default constructor. + */ + ScratchDataObject() + : currently_in_use(false) + {} + + ScratchDataObject(std::unique_ptr &&p, const bool in_use) + : scratch_data(std::move(p)) + , currently_in_use(in_use) + {} + + ScratchDataObject(ScratchData *p, const bool in_use) + : scratch_data(p) + , currently_in_use(in_use) + {} + + // Provide a copy constructor that actually doesn't copy the + // internal state. This makes handling ScratchAndCopyDataObjects + // easier to handle with STL containers. + ScratchDataObject(const ScratchDataObject &) + : currently_in_use(false) + {} + + ScratchDataObject(ScratchDataObject &&o) noexcept = default; + }; + # ifdef DEAL_II_WITH_TBB /** * A namespace for the implementation of details of the WorkStream pattern @@ -196,44 +234,11 @@ namespace WorkStream */ struct ItemType { - /** - * A structure that contains a pointer to a scratch data object - * along with a flag that indicates whether this object is currently - * in use. - */ - struct ScratchDataObject - { - std::unique_ptr scratch_data; - bool currently_in_use; - - /** - * Default constructor. - */ - ScratchDataObject() - : currently_in_use(false) - {} - - ScratchDataObject(ScratchData *p, const bool in_use) - : scratch_data(p) - , currently_in_use(in_use) - {} - - // Provide a copy constructor that actually doesn't copy the - // internal state. This makes handling ScratchAndCopyDataObjects - // easier to handle with STL containers. - ScratchDataObject(const ScratchDataObject &) - : currently_in_use(false) - {} - - ScratchDataObject(ScratchDataObject &&o) noexcept = default; - }; - - /** * Typedef to a list of scratch data objects. The rationale for this * list is provided in the variables that use these lists. */ - using ScratchDataList = std::list; + using ScratchDataList = std::list>; /** * A list of iterators that need to be worked on. Only the first @@ -663,6 +668,152 @@ namespace WorkStream # endif // DEAL_II_WITH_TBB + +# ifdef DEAL_II_WITH_TASKFLOW + /** + * Mostly a copy of the 2nd implementation of the Workstream paper taking + * advantage of thread local lists for re-use. Uses taskflow for task + * scheduling rather than TBB. Currently does not support chunking. + */ + + + namespace taskflow_no_coloring + { + /** + * The main run function for the taskflow colorless implementation. The + * last two arguments in this function are for chunking support which + * currently does not exist but ideally will later. For now they are + * ignored but still here to permit existing programs to function + */ + template + void + run(const Iterator &begin, + const std_cxx20::type_identity_t &end, + Worker worker, + Copier copier, + const ScratchData &sample_scratch_data, + const CopyData &sample_copy_data, + const unsigned int /*queue_length*/ = 2 * + MultithreadInfo::n_threads(), + const unsigned int /*chunk_size*/ = 8) + + { + tf::Executor &executor = MultithreadInfo::get_taskflow_executor(); + tf::Taskflow taskflow; + + using ScratchDataList = std::list>; + + Threads::ThreadLocalStorage + thread_safe_scratch_data_list; + + tf::Task last_copier; + + // idx is used to connect each worker to its copier as communication + // between tasks is not supported. It does this by providing a unique + // index in the vector of pointers copy_datas at which the copy data + // object where the work done by work task #idx is stored + unsigned int idx = 0; + + std::vector> copy_datas; + + // Generate a static task graph. Here we generate a task for each cell + // that will be worked on. The tasks are not executed until all of them + // are created, this code runs sequentially. + for (Iterator it = begin; it != end; ++it, ++idx) + { + copy_datas.emplace_back(); + // Create a worker task. + auto worker_task = + taskflow + .emplace([it, + idx, + &thread_safe_scratch_data_list, + &sample_scratch_data, + &sample_copy_data, + ©_datas, + &worker]() { + ScratchData *scratch_data = nullptr; + + ScratchDataList &scratch_data_list = + thread_safe_scratch_data_list.get(); + // See if there is an unused object. if so, + // grab it and mark it as used. + for (auto &p : scratch_data_list) + { + if (p.currently_in_use == false) + { + scratch_data = p.scratch_data.get(); + p.currently_in_use = true; + break; + } + } + // If no element in the list was found, create + // one and mark it as used. + if (scratch_data == nullptr) + { + scratch_data_list.emplace_back( + std::make_unique(sample_scratch_data), + true); + scratch_data = + scratch_data_list.back().scratch_data.get(); + } + + // Create a unique copy data object where this + // worker's work will be stored. + auto © = copy_datas[idx]; + copy = std::make_unique(sample_copy_data); + worker(it, *scratch_data, *copy.get()); + + // Find our currently used scratch data and + // mark it as unused. + for (auto &p : scratch_data_list) + { + if (p.scratch_data.get() == scratch_data) + { + Assert(p.currently_in_use == true, + ExcInternalError()); + p.currently_in_use = false; + } + } + }) + .name("worker"); + + // Create a copier task. This task is a seperate object from the + // worker task. + tf::Task copier_task = taskflow + .emplace([idx, ©_datas, &copier]() { + copier(*copy_datas[idx].get()); + copy_datas[idx].reset(); + }) + .name("copy"); + + // Ensure the copy task runs after the worker task. + worker_task.precede(copier_task); + + // Ensure that only one copy task can run at a time. The code below + // makes each copy task wait until the previous one has finished + // before it can start + if (!last_copier.empty()) + last_copier.precede(copier_task); + + // Keep a handle to the last copier. Tasks in taskflow are + // basically handles to internally stored data, so this does not + // perform a copy: + last_copier = copier_task; + } + + // Now we run all the tasks in the task graph. They will be run in + // parallel and are eligible to run when their dependencies established + // above are met. + executor.run(taskflow).wait(); + } + } // namespace taskflow_no_coloring +# endif + /** * A reference implementation without using multithreading to be used if we * don't have multithreading support or if the user requests to run things @@ -1144,10 +1295,20 @@ namespace WorkStream if (MultithreadInfo::n_threads() > 1) { -# ifdef DEAL_II_WITH_TBB +# if defined(DEAL_II_WITH_TBB) || defined(DEAL_II_WITH_TASKFLOW) if (static_cast &>(copier)) { // If we have a copier, run the algorithm: +# if defined(DEAL_II_WITH_TASKFLOW) + internal::taskflow_no_coloring::run(begin, + end, + worker, + copier, + sample_scratch_data, + sample_copy_data, + queue_length, + chunk_size); +# elif defined(DEAL_II_WITH_TBB) internal::tbb_no_coloring::run(begin, end, worker, @@ -1156,6 +1317,7 @@ namespace WorkStream sample_copy_data, queue_length, chunk_size); +# endif } else { @@ -1191,7 +1353,7 @@ namespace WorkStream # endif } - // no TBB installed or we are requested to run sequentially: + // no TBB or Taskflow installed or we are requested to run sequentially: internal::sequential::run( begin, end, worker, copier, sample_scratch_data, sample_copy_data); } diff --git a/tests/mpi/mesh_worker_01.cc b/tests/mpi/mesh_worker_01.cc index 6016a26ecd..60b0650e8f 100644 --- a/tests/mpi/mesh_worker_01.cc +++ b/tests/mpi/mesh_worker_01.cc @@ -241,9 +241,9 @@ test() int main(int argc, char **argv) { - Utilities::MPI::MPI_InitFinalize mpi_initialization( - argc, argv, testing_max_num_threads()); - MPILogInitAll log; + // Disable multithreading so that text output order is consistent + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + MPILogInitAll log; test<2>(); } diff --git a/tests/mpi/mesh_worker_02.cc b/tests/mpi/mesh_worker_02.cc index 13a2247014..f384739bd5 100644 --- a/tests/mpi/mesh_worker_02.cc +++ b/tests/mpi/mesh_worker_02.cc @@ -200,9 +200,9 @@ test() int main(int argc, char **argv) { - Utilities::MPI::MPI_InitFinalize mpi_initialization( - argc, argv, testing_max_num_threads()); - MPILogInitAll log; + // Disable multithreading so that text output order is consistent + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + MPILogInitAll log; test<2>(); } diff --git a/tests/mpi/mesh_worker_04.cc b/tests/mpi/mesh_worker_04.cc index c05ac3aa77..d4b03cbb69 100644 --- a/tests/mpi/mesh_worker_04.cc +++ b/tests/mpi/mesh_worker_04.cc @@ -193,9 +193,9 @@ test() int main(int argc, char **argv) { - Utilities::MPI::MPI_InitFinalize mpi_initialization( - argc, argv, testing_max_num_threads()); - MPILogInitAll log; + // Disable multithreading so that text output order is consistent + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + MPILogInitAll log; test<2>(); } -- 2.39.5