]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Use parallel_for instead of a pipeline in implementation3 of WorkStream. 298/head
authorBruno Turcksin <bruno.turcksin@gmail.com>
Fri, 12 Dec 2014 18:05:37 +0000 (12:05 -0600)
committerBruno Turcksin <bruno.turcksin@gmail.com>
Fri, 12 Dec 2014 18:08:19 +0000 (12:08 -0600)
include/deal.II/base/work_stream.h

index 0322e056a09bfd29c668256b924b2e129e45de66..4f86dec6cfba07713195183332cec4793009ff95 100644 (file)
@@ -728,48 +728,34 @@ namespace WorkStream
     namespace Implementation3
     {
       /**
-       * A class that creates a sequence of items from a range of iterators.
+       * A structure that contains a pointer to scratch and copy data objects along
+       * with a flag that indicates whether this object is currently in use.
        */
       template <typename Iterator,
                 typename ScratchData,
                 typename CopyData>
-      class IteratorRangeToItemStream : public tbb::filter
+      struct ScratchAndCopyDataObjects
       {
-      public:
+        std_cxx11::shared_ptr<ScratchData> scratch_data;
+        std_cxx11::shared_ptr<CopyData>    copy_data;
+        bool                               currently_in_use;
+
         /**
-         * A data type that we use to identify items to be worked on. This is
-         * the structure that is passed around between the different parts of
-         * the WorkStream implementation to identify what needs to be done by
-         * the various stages of the pipeline.
+         * Default constructor.
          */
-        struct ItemType
-        {
-          /**
-           * A structure that contains a pointer to scratch and copy data objects along
-           * with a flag that indicates whether this object is currently in use.
-           */
-          struct ScratchAndCopyDataObjects
-          {
-            std_cxx11::shared_ptr<ScratchData> scratch_data;
-            std_cxx11::shared_ptr<CopyData>    copy_data;
-            bool                               currently_in_use;
-
-            /**
-             * Default constructor.
-             */
-            ScratchAndCopyDataObjects ()
-              :
-              currently_in_use (false)
-            {}
+        ScratchAndCopyDataObjects ()
+          :
+          currently_in_use (false)
+        {}
 
-            ScratchAndCopyDataObjects (ScratchData *p,
-                                       CopyData *q,
-                                       const bool in_use)
-              :
-              scratch_data (p),
-              copy_data (q),
-              currently_in_use (in_use)
-            {}
+        ScratchAndCopyDataObjects (ScratchData *p,
+                                   CopyData *q,
+                                   const bool in_use)
+          :
+          scratch_data (p),
+          copy_data (q),
+          currently_in_use (in_use)
+        {}
 
 //TODO: when we push back an object to the list of scratch objects, in
 //      Worker::operator(), we first create an object and then copy
@@ -779,452 +765,69 @@ namespace WorkStream
 //      the original one goes out of scope right away again, but it's
 //      certainly awkward. one way to avoid this would be to use unique_ptr
 //      but we'd need to figure out a way to use it in non-C++11 mode
-            ScratchAndCopyDataObjects (const ScratchAndCopyDataObjects &o)
-              :
-              scratch_data (o.scratch_data),
-              copy_data (o.copy_data),
-              currently_in_use (o.currently_in_use)
-            {}
-          };
-
-
-          /**
-           * Typedef to a list of scratch data objects. The rationale for this
-           * list is provided in the variables that use these lists.
-           */
-          typedef std::list<ScratchAndCopyDataObjects> ScratchAndCopyDataList;
-
-          /**
-           * A list of iterators that need to be worked on. Only the first
-           * n_items are relevant.
-           */
-          std::vector<Iterator> work_items;
-
-          /**
-           * Number of items identified by the work_items array that the
-           * Worker and Copier pipeline stage need to work on. The maximum
-           * value of this variable will be chunk_size.
-           */
-          unsigned int          n_items;
-
-          /**
-           * Pointer to a thread local variable identifying the scratch
-           * and copy data objects this thread will use. The same considerations
-           * apply as documented in the Implementation2::IteratorRangeToItemStream
-           * class as well as in the paper by Turcksin, Kronbichler and Bangerth
-          * (see @ref workstream_paper).
-           */
-          Threads::ThreadLocalStorage<ScratchAndCopyDataList> *scratch_and_copy_data;
-
-          /**
-           * Pointer to a sample scratch data object, to be used to initialize
-           * the scratch data objects created for each individual thread.
-           */
-          const ScratchData *sample_scratch_data;
-
-          /**
-           * Pointer to a sample copy data object.
-           */
-          const CopyData *sample_copy_data;
-
-          /**
-           * Flag is true if the buffer is used and false if the buffer can be
-           * used.
-           */
-          bool currently_in_use;
-
-          /**
-           * Default constructor.
-           * Initialize everything that doesn't have a default constructor
-           * itself.
-           */
-          ItemType ()
-            :
-            n_items (0),
-            sample_scratch_data (0),
-            sample_copy_data (0),
-            currently_in_use (false)
-          {}
-        };
-
-
-        /**
-         * Constructor. Take range of iterators into an array of the kind of object we
-         * want to loop over, the size of a buffer that can
-         * hold items, and the sample additional data object that will be passed
-         * to each worker and copier function invokation.
-         */
-        IteratorRangeToItemStream (const typename std::vector<Iterator>::const_iterator &begin,
-                                   const typename std::vector<Iterator>::const_iterator &end,
-                                   const unsigned int    buffer_size,
-                                   const unsigned int    chunk_size,
-                                   const ScratchData    &sample_scratch_data,
-                                   const CopyData       &sample_copy_data)
+        ScratchAndCopyDataObjects (const ScratchAndCopyDataObjects &o)
           :
-          tbb::filter (/*is_serial=*/true),
-          remaining_iterator_range (begin, end),
-          item_buffer (buffer_size),
-          sample_scratch_data (sample_scratch_data),
-          sample_copy_data (sample_copy_data),
-          chunk_size (chunk_size)
-        {
-          Assert (begin != end, ExcMessage ("This class is not prepared to deal with empty ranges!"));
-          // initialize the elements of the item_buffer
-          for (unsigned int element=0; element<item_buffer.size(); ++element)
-            {
-              Assert (item_buffer[element].n_items == 0,
-                      ExcInternalError());
-
-              // resize the item_buffer. we have to initialize the new
-              // elements with something and because we can't rely on there
-              // being a default constructor for 'Iterator', we take the first
-              // element in the range [begin,end) pointed to.
-              item_buffer[element].work_items.resize (chunk_size,*begin);
-              item_buffer[element].scratch_and_copy_data = &thread_local_scratch_and_copy;
-              item_buffer[element].sample_scratch_data = &sample_scratch_data;
-              item_buffer[element].sample_copy_data = &sample_copy_data;
-              item_buffer[element].currently_in_use = false;
-            }
-        }
-
-
-        /**
-         * Create an item and return a
-         * pointer to it.
-         */
-        virtual void *operator () (void *)
-        {
-          // find first unused item. we know that there must be one
-          // because we have set the maximal number of tokens in flight
-          // and have set the ring buffer to have exactly this size. so
-          // if this function is called, we know that less than the
-          // maximal number of items in currently in flight
-          ItemType *current_item = 0;
-          for (unsigned int i=0; i<item_buffer.size(); ++i)
-            if (item_buffer[i].currently_in_use == false)
-              {
-                item_buffer[i].currently_in_use = true;
-                current_item = &item_buffer[i];
-                break;
-              }
-          Assert (current_item != 0, ExcMessage ("This can't be. There must be a free item!"));
-
-
-          // initialize the next item. it may
-          // consist of at most chunk_size
-          // elements
-          current_item->n_items = 0;
-          while ((remaining_iterator_range.first !=
-                  remaining_iterator_range.second)
-                 &&
-                 (current_item->n_items < chunk_size))
-            {
-              // initialize the iterators to work on with the elements
-              // of the vector that remaining_iterator_range
-              // points into
-              current_item->work_items[current_item->n_items]
-                = *remaining_iterator_range.first;
-
-              ++remaining_iterator_range.first;
-              ++current_item->n_items;
-            }
-
-          if (current_item->n_items == 0)
-            // there were no items
-            // left. terminate the pipeline
-            return 0;
-          else
-            return current_item;
-        }
-
-      private:
-        /**
-         * The interval of iterators still to
-         * be worked on. This range will shrink
-         * over time.
-         */
-        std::pair<typename std::vector<Iterator>::const_iterator,typename std::vector<Iterator>::const_iterator> remaining_iterator_range;
-
-        /**
-         * A buffer that will store items.
-         */
-        std::vector<ItemType>        item_buffer;
-
-        /**
-         * Pointer to a thread local variable identifying the scratch and
-         * copy data objects each thread will use. The same is true as
-         * discussed for the implementation in the
-         * Implementation2::IteratorRangeToItemStream class and the paper
-         * by Turcksin, Kronbichler and Bangerth (see @ref workstream_paper).
-         */
-        Threads::ThreadLocalStorage<typename ItemType::ScratchAndCopyDataList> thread_local_scratch_and_copy;
-
-        /**
-         * A reference to a sample scratch data that will be used to
-         * initialize the thread-local pointers to a scratch data object
-         * each of the worker tasks uses.
-         */
-        const ScratchData &sample_scratch_data;
-
-        /**
-         * A reference to a sample scratch data that will be used to
-         * initialize the thread-local pointers to a scratch data object
-         * each of the worker tasks uses.
-         */
-        const CopyData &sample_copy_data;
-
-        /**
-         * Number of elements of the
-         * iterator range that each
-         * thread should work on
-         * sequentially; a large number
-         * makes sure that each thread
-         * gets a significant amount of
-         * work before the next task
-         * switch happens, whereas a
-         * small number is better for
-         * load balancing.
-         */
-        const unsigned int           chunk_size;
+          scratch_data (o.scratch_data),
+          copy_data (o.copy_data),
+          currently_in_use (o.currently_in_use)
+        {}
+      };
 
-        /**
-         * Initialize the pointers and vector
-         * elements in the specified entry of
-         * the item_buffer.
-         */
-        void init_buffer_elements (const unsigned int element)
-        {
-          Assert (item_buffer[element].n_items == 0,
-                  ExcInternalError());
 
-          item_buffer[element].work_items
-          .resize (chunk_size, remaining_iterator_range.second);
-          item_buffer[element].scratch_and_copy_data
-            = &thread_local_scratch_and_copy;
-          item_buffer[element].sample_scratch_data
-            = &sample_scratch_data;
-          item_buffer[element].sample_copy_data
-            = &sample_copy_data;
-        }
-      };
 
 
 
       /**
-       * A class that manages calling the
-       * worker function on a number of
-       * parallel threads. Note that it is, in
-       * the TBB notation, a filter that can
-       * run in parallel.
+       * A class that manages calling the worker and copier functions. Unlike the
+       * other implementations, parallel_for is used instead of a pipeline.
        */
       template <typename Iterator,
                 typename ScratchData,
                 typename CopyData>
-      class WorkerAndCopier : public tbb::filter
+      class WorkerAndCopier
       {
       public:
         /**
-         * Constructor. Takes a
-         * reference to the object on
-         * which we will operate as
-         * well as a pointer to the
-         * function that will do the
-         * assembly.
+         * Constructor.
          */
         WorkerAndCopier (const std_cxx11::function<void (const Iterator &,
                                                          ScratchData &,
                                                          CopyData &)> &worker,
-                         const std_cxx11::function<void (const CopyData &)> &copier)
-          :
-          tbb::filter (/* is_serial= */ false),
-          worker (worker),
-          copier (copier)
-        {}
-
-
-        /**
-         * Work on an item.
-         */
-        void *operator () (void *item)
-        {
-          // first unpack the current item
-          typedef
-          typename IteratorRangeToItemStream<Iterator,ScratchData,CopyData>::ItemType
-          ItemType;
-
-          ItemType *current_item = static_cast<ItemType *> (item);
-
-          // 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 = current_item->scratch_and_copy_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(*current_item->sample_scratch_data);
-                copy_data    = new CopyData(*current_item->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 (unsigned int i=0; i<current_item->n_items; ++i)
-            {
-              try
-                {
-                  if (worker)
-                    worker (current_item->work_items[i],
-                            *scratch_data,
-                            *copy_data);
-                  if (copier)
-                    copier (*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 = current_item->scratch_and_copy_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;
-                }
-          }
-
-          // mark current item as usable again
-          current_item->currently_in_use = false;
-
-          // return an invalid item since we are at the end of the
-          // pipeline
-          return 0;
-        }
-
-
-      private:
-        /**
-         * Pointer to the function
-         * that does the assembling
-         * on the sequence of cells.
-         */
-        const std_cxx11::function<void (const Iterator &,
-                                        ScratchData &,
-                                        CopyData &)> worker;
-
-        /**
-         * Pointer to the function that does the copying from
-         * local contribution to global object.
-         */
-        const std_cxx11::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_cxx11::function<void (const Iterator &,
-                                                ScratchData &,
-                                                CopyData &)> &worker,
-                const ScratchData    &sample_scratch_data,
-                const CopyData       &sample_copy_data)
+                         const std_cxx11::function<void (const CopyData &)> &copier,
+                         const ScratchData    &sample_scratch_data,
+                         const CopyData       &sample_copy_data)
           :
           worker (worker),
+          copier (copier),
           sample_scratch_data (sample_scratch_data),
           sample_copy_data (sample_copy_data)
         {}
 
 
         /**
-         * The function that calls the worker function on a
+         * The function that calls the worker and the copier functions 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
+          // 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();
+            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
+            for (typename ScratchAndCopyDataList::iterator
                  p = scratch_and_copy_data_list.begin();
                  p != scratch_and_copy_data_list.end(); ++p)
               if (p->currently_in_use == false)
@@ -1242,25 +845,25 @@ namespace WorkStream
                 scratch_data = new ScratchData(sample_scratch_data);
                 copy_data    = new CopyData(sample_copy_data);
 
-                typename ItemType::ScratchAndCopyDataList::value_type
+                typename 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
+          // then call the worker and copier functions on each
+          // element of the chunk we were given.
           for (typename std::vector<Iterator>::const_iterator p=range.begin();
                p != range.end(); ++p)
             {
               try
                 {
-                  worker (*p,
-                          *scratch_data,
-                          *copy_data);
+                  if (worker)
+                    worker (*p,
+                            *scratch_data,
+                            *copy_data);
+                  if (copier)
+                    copier (*copy_data);
                 }
               catch (const std::exception &exc)
                 {
@@ -1276,10 +879,9 @@ namespace WorkStream
           // 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();
+            ScratchAndCopyDataList &scratch_and_copy_data_list = data.get();
 
-            for (typename ItemType::ScratchAndCopyDataList::iterator p =
+            for (typename ScratchAndCopyDataList::iterator p =
                    scratch_and_copy_data_list.begin(); p != scratch_and_copy_data_list.end();
                  ++p)
               if (p->scratch_data.get() == scratch_data)
@@ -1293,12 +895,14 @@ namespace WorkStream
 
       private:
         typedef
-        typename Implementation3::IteratorRangeToItemStream<Iterator,ScratchData,CopyData>::ItemType
-        ItemType;
+        typename Implementation3::ScratchAndCopyDataObjects<Iterator,ScratchData,CopyData>
+        ScratchAndCopyDataObjects;
 
-        typedef
-        typename ItemType::ScratchAndCopyDataList
-        ScratchAndCopyDataList;
+        /**
+         * Typedef to a list of scratch data objects. The rationale for this
+         * list is provided in the variables that use these lists.
+         */
+        typedef std::list<ScratchAndCopyDataObjects> ScratchAndCopyDataList;
 
         Threads::ThreadLocalStorage<ScratchAndCopyDataList> data;
 
@@ -1311,6 +915,12 @@ namespace WorkStream
                                         ScratchData &,
                                         CopyData &)> worker;
 
+        /**
+         * Pointer to the function that does the copying from
+         * local contribution to global object.
+         */
+        const std_cxx11::function<void (const CopyData &)> copier;
+
         /**
          * References to sample scratch and copy data for
          * when we need them.
@@ -1575,63 +1185,27 @@ namespace WorkStream
         for (unsigned int color=0; color<colored_iterators.size(); ++color)
           if (colored_iterators[color].size() > 0)
             {
-              if (static_cast<const std_cxx11::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);
-
+              typedef
+              internal::Implementation3::WorkerAndCopier<Iterator,ScratchData,CopyData>
+              WorkerAndCopier;
 
-                  internal::Implementation3::WorkerAndCopier<Iterator, ScratchData, CopyData>
-                  worker_and_copier_filter (worker, copier);
+              typedef
+              typename std::vector<Iterator>::const_iterator
+              RangeType;
 
-                  // 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);
+              WorkerAndCopier worker_and_copier (worker,
+                                                 copier,
+                                                 sample_scratch_data,
+                                                 sample_copy_data);
 
-                  assembly_line.clear ();
-                }
-              else
-                {
-                  // no copier function, we can implement things as a parallel for
-                  Assert (static_cast<const std_cxx11::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_cxx11::bind (&ParallelForWorker::operator(),
-                                                      std_cxx11::ref(parallel_for_worker),
-                                                      std_cxx11::_1),
-                                     tbb::auto_partitioner());
-                }
+              tbb::parallel_for (tbb::blocked_range<RangeType>
+                                 (colored_iterators[color].begin(),
+                                  colored_iterators[color].end(),
+                                  /*grain_size=*/chunk_size),
+                                 std_cxx11::bind (&WorkerAndCopier::operator(),
+                                                  std_cxx11::ref(worker_and_copier),
+                                                  std_cxx11::_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.