From 73ccf43c3240c9b4637f355199e905ce084bcc2a Mon Sep 17 00:00:00 2001 From: David Wells Date: Wed, 27 May 2020 15:54:35 -0400 Subject: [PATCH] Convert ThreadGroup to TaskGroup in sparse vanka --- include/deal.II/lac/sparse_vanka.templates.h | 84 +++++++------------- 1 file changed, 27 insertions(+), 57 deletions(-) diff --git a/include/deal.II/lac/sparse_vanka.templates.h b/include/deal.II/lac/sparse_vanka.templates.h index 00fc4937d7..06c0f5d571 100644 --- a/include/deal.II/lac/sparse_vanka.templates.h +++ b/include/deal.II/lac/sparse_vanka.templates.h @@ -111,13 +111,15 @@ SparseVanka::compute_inverses() #else const size_type n_inverses = std::count(selected->begin(), selected->end(), true); - const std::size_t n_threads = MultithreadInfo::n_threads(); - const size_type n_inverses_per_thread = - std::max(n_inverses / n_threads, static_cast(1U)); + // somewhat arbitrarily set up an equal number of tasks as we have threads + const std::size_t n_tasks = MultithreadInfo::n_threads(); + const size_type n_inverses_per_task = + std::max(static_cast(n_inverses / n_tasks), + static_cast(1U)); // set up start and end index // for each of the - // threads. note that we have + // tasks. note that we have // to work somewhat to get this // appropriate, since the // indices for which inverses @@ -131,39 +133,34 @@ SparseVanka::compute_inverses() // consecutive, with other // consecutive regions where we // do not have to do something - std::vector> blocking(n_threads); + std::vector> blocking(n_tasks); unsigned int c = 0; - unsigned int thread = 0; + unsigned int task_n = 0; blocking[0].first = 0; - for (size_type i = 0; (i < matrix->m()) && (thread + 1 < n_threads); ++i) + for (size_type i = 0; (i < matrix->m()) && (task_n + 1 < n_tasks); ++i) { if ((*selected)[i] == true) ++c; - if (c == n_inverses_per_thread) + if (c == n_inverses_per_task) { - blocking[thread].second = i; - blocking[thread + 1].first = i; - ++thread; + blocking[task_n].second = i; + blocking[task_n + 1].first = i; + ++task_n; c = 0; } } - blocking[n_threads - 1].second = matrix->m(); - - using FunPtr = - void (SparseVanka::*)(const size_type, const size_type); - const FunPtr fun_ptr = &SparseVanka::compute_inverses; - - // Now spawn the threads - Threads::ThreadGroup<> threads; - for (unsigned int i = 0; i < n_threads; ++i) - threads += Threads::new_thread(fun_ptr, - *this, - blocking[i].first, - blocking[i].second); - threads.join_all(); + blocking[n_tasks - 1].second = matrix->m(); + + // Now spawn the tasks + Threads::TaskGroup<> tasks; + for (unsigned int i = 0; i < n_tasks; ++i) + tasks += Threads::new_task([&, i]() { + this->compute_inverses(blocking[i].first, blocking[i].second); + }); + tasks.join_all(); #endif } @@ -603,38 +600,11 @@ SparseBlockVanka::vmult(Vector & dst, { dst = 0; - // if no blocking is required, pass down to the underlying class - if (n_blocks == 1) - this->apply_preconditioner(dst, src); - else - // otherwise: blocking requested - { -#ifdef DEAL_II_WITH_THREADS - // spawn threads. since some compilers have trouble finding out which - // 'encapsulate' function to take of all those possible ones if we simply - // idrop in the address of an overloaded template member function, make it - // simpler for the compiler by giving it the correct type right away: - using mem_fun_p = - void (SparseVanka::*)(Vector &, - const Vector &, - const std::vector *const) const; - - const mem_fun_p comp = - &SparseVanka::template apply_preconditioner; - Threads::ThreadGroup<> threads; - for (unsigned int block = 0; block < n_blocks; ++block) - threads += - Threads::new_thread(comp, - *static_cast *>(this), - dst, - src, - &dof_masks[block]); - threads.join_all(); -#else - for (unsigned int block = 0; block < n_blocks; ++block) - this->apply_preconditioner(dst, src, &dof_masks[block]); -#endif - } + Threads::TaskGroup<> tasks; + for (unsigned int block = 0; block < n_blocks; ++block) + tasks += Threads::new_task( + [&, block] { this->apply_preconditioner(dst, src, &dof_masks[block]); }); + tasks.join_all(); } -- 2.39.5