]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Implement a parallel for approach for WorkStream::run when the copier function is...
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 26 Dec 2013 16:03:28 +0000 (16:03 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 26 Dec 2013 16:03:28 +0000 (16:03 +0000)
git-svn-id: https://svn.dealii.org/trunk@32124 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/include/deal.II/base/work_stream.h

index ccbfba6bc884429617cb8d5ed6018e4018084360..e398c01f87cf74c1fc2d2d8e9491a2a389a863e7 100644 (file)
@@ -85,11 +85,20 @@ inconvenience this causes.
 
 <ol>
 
+  <li> Improved: When you call WorkStream::run with an empty function object
+  for the copier, operations on individual cells are essentially all independent.
+  In other words, you have a massively parallel collection of jobs. In this
+  case, a parallel for loop over all elements is better suited than the
+  pipeline approach currently used. This has now been implemented.
+  <br>
+  (Wolfgang Bangerth, 2013/12/26)
+  </li>
+
   <li> New: The new function VectorTools::interpolate_based_on_material_id()
   can be used to interpolate several functions onto a mesh, based on the
   material id of each cell individually.
   <br>
-  (Valentin Zingan, 2013/12/24)
+  (Valentin Zingan, 2013/12/26)
   </li>
 
   <li> New: A new reinit() method has been introduced to
index 4c3eb912cfc82f35ff437b278f3ae1ac2b0cc82d..fd30a26371278e0f52ff06e2dd398676516a354a 100644 (file)
@@ -26,6 +26,7 @@
 #include <deal.II/base/std_cxx1x/function.h>
 #include <deal.II/base/std_cxx1x/bind.h>
 #include <deal.II/base/thread_local_storage.h>
+#include <deal.II/base/parallel.h>
 
 #ifdef DEAL_II_WITH_THREADS
 #  include <deal.II/base/thread_management.h>
@@ -1166,8 +1167,163 @@ namespace WorkStream
         const std_cxx1x::function<void (const CopyData &)> copier;
       };
     }
+
+
+    /**
+     * A namespace for functions used in the implementation of
+     * Implementation3 whenever we don't actually have a copier
+     * function. In that case, we can simply use parallel_for, but we
+     * have to be careful with the use of thread-local objects as we
+     * had to be with Implementation3 as well.
+     */
+    namespace ParallelFor
+    {
+      template <typename Iterator,
+               typename ScratchData,
+               typename CopyData>
+      class Worker
+      {
+      public:
+       /**
+        * Constructor.
+        */
+       Worker (const std_cxx1x::function<void (const Iterator &,
+                                               ScratchData &,
+                                               CopyData &)> &worker,
+               const ScratchData    &sample_scratch_data,
+               const CopyData       &sample_copy_data)
+       :
+       worker (worker),
+       sample_scratch_data (sample_scratch_data),
+       sample_copy_data (sample_copy_data)
+         {}
+
+
+       /**
+        * The function that calls the worker function on a
+        * range of items denoted by the two arguments.
+        */
+       void operator() (const tbb::blocked_range<typename std::vector<Iterator>::const_iterator> &range)
+         {
+           // we need to find an unused scratch and corresponding copy
+           // data object in the list that
+           // corresponds to the current thread and then mark it as used. if
+           // we can't find one, create one
+           //
+           // as discussed in the discussion of the documentation of the
+           // IteratorRangeToItemStream::scratch_data variable, there is no
+           // need to synchronize access to this variable using a mutex
+           // as long as we have no yield-point in between. this means that
+           // we can't take an iterator into the list now and expect it to
+           // still be valid after calling the worker, but we at least do
+           // not have to lock the following section
+           ScratchData *scratch_data = 0;
+           CopyData    *copy_data    = 0;
+           {
+             typename ItemType::ScratchAndCopyDataList &
+               scratch_and_copy_data_list = data.get();
+
+             // see if there is an unused object. if so, grab it and mark
+             // it as used
+             for (typename ItemType::ScratchAndCopyDataList::iterator
+                    p = scratch_and_copy_data_list.begin();
+                  p != scratch_and_copy_data_list.end(); ++p)
+               if (p->currently_in_use == false)
+                 {
+                   scratch_data = p->scratch_data.get();
+                   copy_data    = p->copy_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 == 0)
+               {
+                 Assert (copy_data==0, ExcInternalError());
+                 scratch_data = new ScratchData(sample_scratch_data);
+                 copy_data    = new CopyData(sample_copy_data);
+
+                 typename ItemType::ScratchAndCopyDataList::value_type
+                   new_scratch_object (scratch_data, copy_data, true);
+                 scratch_and_copy_data_list.push_back (new_scratch_object);
+               }
+           }
+
+           // then call the worker and copier function on each
+           // element of the chunk we were given. since these
+           // functions are called on separate threads, nothing good
+           // can happen if they throw an exception and we are best
+           // off catching it and showing an error message
+           for (typename std::vector<Iterator>::const_iterator p=range.begin();
+                p != range.end(); ++p)
+             {
+               try
+                 {
+                   worker (*p,
+                           *scratch_data,
+                           *copy_data);
+                 }
+               catch (const std::exception &exc)
+                 {
+                   Threads::internal::handle_std_exception (exc);
+                 }
+               catch (...)
+                 {
+                   Threads::internal::handle_unknown_exception ();
+                 }
+             }
+
+           // finally mark the scratch object as unused again. as above, there
+           // is no need to lock anything here since the object we work on
+           // is thread-local
+           {
+             typename ItemType::ScratchAndCopyDataList &
+               scratch_and_copy_data_list = data.get();
+
+             for (typename ItemType::ScratchAndCopyDataList::iterator p =
+                    scratch_and_copy_data_list.begin(); p != scratch_and_copy_data_list.end();
+                  ++p)
+               if (p->scratch_data.get() == scratch_data)
+                 {
+                   Assert(p->currently_in_use == true, ExcInternalError());
+                   p->currently_in_use = false;
+                 }
+           }
+
+         }
+
+      private:
+       typedef
+       typename Implementation3::IteratorRangeToItemStream<Iterator,ScratchData,CopyData>::ItemType
+       ItemType;
+
+       typedef
+       typename ItemType::ScratchAndCopyDataList
+       ScratchAndCopyDataList;
+
+       Threads::ThreadLocalStorage<ScratchAndCopyDataList> data;
+
+        /**
+         * Pointer to the function
+         * that does the assembling
+         * on the sequence of cells.
+         */
+        const std_cxx1x::function<void (const Iterator &,
+                                        ScratchData &,
+                                        CopyData &)> worker;
+
+       /**
+        * References to sample scratch and copy data for
+        * when we need them.
+        */
+       const ScratchData    &sample_scratch_data;
+       const CopyData       &sample_copy_data;
+      };
+    }
+
   }
 
+
 #endif // DEAL_II_WITH_THREADS
 
 
@@ -1419,30 +1575,65 @@ namespace WorkStream
         // loop over the various colors of what we're given
         for (unsigned int color=0; color<colored_iterators.size(); ++color)
           if (colored_iterators[color].size() > 0)
-            {
-              // create the three stages of the pipeline
-              internal::Implementation3::IteratorRangeToItemStream<Iterator,ScratchData,CopyData>
-              iterator_range_to_item_stream (colored_iterators[color].begin(),
-                                             colored_iterators[color].end(),
-                                             queue_length,
-                                             chunk_size,
-                                             sample_scratch_data,
-                                             sample_copy_data);
-
-
-              internal::Implementation3::WorkerAndCopier<Iterator, ScratchData, CopyData>
-              worker_and_copier_filter (worker, copier);
-
-              // now create a pipeline from these stages
-              tbb::pipeline assembly_line;
-              assembly_line.add_filter (iterator_range_to_item_stream);
-              assembly_line.add_filter (worker_and_copier_filter);
-
-              // and run it
-              assembly_line.run (queue_length);
-
-              assembly_line.clear ();
-            }
+           {
+             if (static_cast<const std_cxx1x::function<void (const CopyData &)>& >(copier))
+               {
+                 // there is a copier function, so we have to go with
+                 // the full three-stage design of the pipeline
+                 internal::Implementation3::IteratorRangeToItemStream<Iterator,ScratchData,CopyData>
+                   iterator_range_to_item_stream (colored_iterators[color].begin(),
+                                                  colored_iterators[color].end(),
+                                                  queue_length,
+                                                  chunk_size,
+                                                  sample_scratch_data,
+                                                  sample_copy_data);
+
+
+                 internal::Implementation3::WorkerAndCopier<Iterator, ScratchData, CopyData>
+                   worker_and_copier_filter (worker, copier);
+
+                 // now create a pipeline from these stages
+                 tbb::pipeline assembly_line;
+                 assembly_line.add_filter (iterator_range_to_item_stream);
+                 assembly_line.add_filter (worker_and_copier_filter);
+
+                 // and run it
+                 assembly_line.run (queue_length);
+
+                 assembly_line.clear ();
+               }
+             else
+               {
+                 // no copier function, we can implement things as a parallel for
+                 Assert (static_cast<const std_cxx1x::function<void (const Iterator &,
+                                                                     ScratchData &,
+                                                                     CopyData &)>& >(worker),
+                         ExcMessage ("It makes no sense to call this function with "
+                                     "empty functions for both the worker and the "
+                                     "copier!"));
+
+                 typedef
+                   internal::ParallelFor::Worker<Iterator,ScratchData,CopyData>
+                   ParallelForWorker;
+
+                 typedef
+                   typename std::vector<Iterator>::const_iterator
+                   RangeType;
+
+                 ParallelForWorker parallel_for_worker (worker,
+                                                        sample_scratch_data,
+                                                        sample_copy_data);
+
+                 tbb::parallel_for (tbb::blocked_range<RangeType>
+                                    (colored_iterators[color].begin(),
+                                     colored_iterators[color].end(),
+                                     /*grain_size=*/chunk_size),
+                                    std_cxx1x::bind (&ParallelForWorker::operator(),
+                                                     std_cxx1x::ref(parallel_for_worker),
+                                                     std_cxx1x::_1),
+                                    tbb::auto_partitioner());
+               }
+           }
       }
 #endif
   }

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.