]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add a mechanism for thread local storage to the library. Add corresponding tests.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 18 Jul 2011 20:37:24 +0000 (20:37 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 18 Jul 2011 20:37:24 +0000 (20:37 +0000)
git-svn-id: https://svn.dealii.org/trunk@23950 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/include/deal.II/base/thread_local_storage.h [new file with mode: 0644]
tests/base/thread_local_storage_01.cc [new file with mode: 0644]
tests/base/thread_local_storage_01/cmp/generic [new file with mode: 0644]
tests/base/thread_local_storage_02.cc [new file with mode: 0644]
tests/base/thread_local_storage_02/cmp/generic [new file with mode: 0644]

index 049f03b4f4a30424cdece02e340f8c547abc0d23..9eaec5290e4aa6887355246cf33d13957f374409 100644 (file)
@@ -216,6 +216,12 @@ should be fixed now.
 <h3>Specific improvements</h3>
 
 <ol>
+<li> New: There is now a class Threads::ThreadLocalStorage that allows threads
+to have their own copy of an object without having to fear interference from
+other threads in accessing this object.
+<br>
+(Wolfgang Bangerth 2011/07/07)
+
 <li> Fixed: The 2d grid reordering algorithm that is used by all grid readers had
 a component that was quadratic in its complexity, sometimes leading to cases
 where reading in a mesh in debug mode could take minutes for just a few tens
diff --git a/deal.II/include/deal.II/base/thread_local_storage.h b/deal.II/include/deal.II/base/thread_local_storage.h
new file mode 100644 (file)
index 0000000..9f7caab
--- /dev/null
@@ -0,0 +1,153 @@
+//---------------------------------------------------------------------------
+//    $Id$
+//
+//    Copyright (C) 2011 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//---------------------------------------------------------------------------
+#ifndef __deal2__thread_local_storage_h
+#define __deal2__thread_local_storage_h
+
+
+#include <deal.II/base/config.h>
+
+#if DEAL_II_USE_MT == 1
+#  include <tbb/enumerable_thread_specific.h>
+#endif
+
+
+
+DEAL_II_NAMESPACE_OPEN
+
+/*!@addtogroup threads */
+/*@{*/
+
+
+namespace Threads
+{
+  template <typename T>
+  /**
+   * @brief A class that provides a separate storage location on each thread that accesses the object.
+   * 
+   * This class offers ways so that every thread that accesses it has its
+   * own copy of an object of type T. In essence, accessing this object
+   * can never result in race conditions in multithreaded programs since
+   * no other thread than the current one can ever access it.
+   * 
+   * The class builds on the Threading Building Blocks's tbb::enumerable_thread_specific
+   * class but wraps it in such a way that this class can also be used when deal.II
+   * is configured not to use threads at all -- in that case, this class simply
+   * stores a single copy of an object of type T.
+   * 
+   * <h3>Construction and destruction</h3>
+   * 
+   * Objects of this class can either be default constructed or by providing an
+   * "exemplar", i.e. an object of type T so that everytime we need to create
+   * a T on a thread that doesn't already have such an object, it is copied from
+   * the exemplar.
+   * 
+   * Upon destruction of objects of this class, all T objects that correspond
+   * to threads that have accessed this object are destroyed. Note that this may be
+   * before the time when a thread is terminated.
+   * 
+   * <h3>Access</h3>
+   * 
+   * The T object stored by this object can be accessed using the get() function. It
+   * provides a reference to a unique object when accessed from different threads.
+   * Objects of type T are created lazily, i.e. they are only created whenever a
+   * thread actually calls get().
+   **/
+  class ThreadLocalStorage 
+  {
+  public:
+    /**
+     * Default constructor. Initialize each thread local object
+     * using its default constructor.
+     **/
+    ThreadLocalStorage ();
+    
+    /**
+     * A kind of copy constructor. Initialize each thread local object
+     * by copying the given object.
+     **/
+    ThreadLocalStorage (const T &t);
+
+    /**
+     * Copy constructor. Initialize each thread local object
+     * with the corresponding object of the given object.
+     **/
+    ThreadLocalStorage (const ThreadLocalStorage<T> &t);
+
+    /**
+     * Return a reference to the data stored by this object for the current
+     * thread this function is called on.
+     */
+    T & get ();
+    
+  private:
+#if DEAL_II_USE_MT == 1
+    /**
+     * The data element we store. If we support threads, then this
+     * object will be of a type that provides a separate object
+     * for each thread. Otherwise, it is simply a single object 
+     * of type T.
+     **/
+    tbb::enumerable_thread_specific<T> data;
+#else
+    T data;
+#endif    
+  };
+  
+// ----------------- inline and template functions ----------------------------
+  
+  template <typename T>
+  inline
+  ThreadLocalStorage<T>::ThreadLocalStorage()
+  {}
+
+
+  template <typename T>
+  inline
+  ThreadLocalStorage<T>::ThreadLocalStorage(const T &t)
+  :
+  data (t)
+  {}
+
+
+  template <typename T>
+  inline
+  ThreadLocalStorage<T>::ThreadLocalStorage(const ThreadLocalStorage<T> &t)
+  :
+  data (t)
+  {}
+
+
+  template <typename T>
+  inline
+  T &
+  ThreadLocalStorage<T>::get ()
+  {
+#if DEAL_II_USE_MT == 1
+    return data.local();
+#else
+    return data;
+#endif        
+  }
+
+  
+}   // end of implementation of namespace Threads
+
+/**
+ * @}
+ */
+
+
+//---------------------------------------------------------------------------
+DEAL_II_NAMESPACE_CLOSE
+// end of #ifndef __deal2__thread_local_storage_h
+#endif
+//---------------------------------------------------------------------------
diff --git a/tests/base/thread_local_storage_01.cc b/tests/base/thread_local_storage_01.cc
new file mode 100644 (file)
index 0000000..1db83c1
--- /dev/null
@@ -0,0 +1,102 @@
+//-----------------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$ 
+//
+//    Copyright (C) 2008, 2011 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//-----------------------------------------------------------------------------
+
+// verify that thread local storage works as advertised
+
+#include "../tests.h"
+#include <iomanip>
+#include <fstream>
+
+#include <deal.II/base/thread_management.h>
+#include <deal.II/base/thread_local_storage.h>
+
+
+struct X
+{
+    X ()  { deallog << "Creating" << std::endl; };
+    ~X () { deallog << "Destroying " << std::endl; };
+    int i;
+};
+
+Threads::ThreadLocalStorage<X> *tls_data;
+
+volatile int counter = 0;
+
+void execute (int i)
+{
+  tls_data->get().i = i;
+
+                                  // indicate that the TLS object has been
+                                  // accessed
+  static Threads::Mutex m;
+  {
+    Threads::Mutex::ScopedLock l(m);
+    ++counter;
+  }
+
+                                  // wait in order to make sure that the
+                                  // thread lives longer than the TLS object
+  sleep (5);
+}
+
+
+void test () 
+{
+                                  // create a thread local storage object
+  tls_data = new Threads::ThreadLocalStorage<X>();
+  
+                                  // create 5 threads and wait for their
+                                  // return. the OS will create 5 individual
+                                  // thread ids, which means that we will
+                                  // create 5 individual thread specific
+                                  // storage locations
+  Threads::ThreadGroup<> tg;
+  for (unsigned int i=10; i<15; ++i)
+    tg += Threads::new_thread (execute, i);
+
+                                  // spin lock until all threads have created
+                                  // their objects
+  while (counter != 5);
+  
+                                  // delete the TLS object. this should also
+                                  // destroy all the objects created so far,
+                                  // namely 5 of them; while that happens, we
+                                  // should get output into logstream
+  delete tls_data;
+
+                                  // write something into the output
+                                  // file. this also makes sure that the
+                                  // output file records the order in which
+                                  // this output is generated and in which
+                                  // the TLS objects were destroyed. we want
+                                  // that these objects are destroyed when
+                                  // the TLS object is destroyed, not when
+                                  // the threads exit
+  deallog << "Done." << std::endl;
+
+                                  // now make sure the threads all finish
+  tg.join_all ();
+}
+
+  
+  
+
+int main()
+{
+  std::ofstream logfile("thread_local_storage_01/output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  test ();
+}
diff --git a/tests/base/thread_local_storage_01/cmp/generic b/tests/base/thread_local_storage_01/cmp/generic
new file mode 100644 (file)
index 0000000..6353809
--- /dev/null
@@ -0,0 +1,12 @@
+
+DEAL::Creating
+DEAL::Creating
+DEAL::Creating
+DEAL::Creating
+DEAL::Creating
+DEAL::Destroying 
+DEAL::Destroying 
+DEAL::Destroying 
+DEAL::Destroying 
+DEAL::Destroying 
+DEAL::Done.
diff --git a/tests/base/thread_local_storage_02.cc b/tests/base/thread_local_storage_02.cc
new file mode 100644 (file)
index 0000000..f8ab901
--- /dev/null
@@ -0,0 +1,112 @@
+//-----------------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$ 
+//
+//    Copyright (C) 2008, 2011 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//-----------------------------------------------------------------------------
+
+// verify that thread local storage works as advertised. like _01 but using
+// the initialization with an exemplar
+
+#include "../tests.h"
+#include <iomanip>
+#include <fstream>
+
+#include <deal.II/base/thread_management.h>
+#include <deal.II/base/thread_local_storage.h>
+
+
+struct X
+{
+    X ()  { Assert (false, ExcInternalError()); };
+    X (int n)  { deallog << "Creating" << std::endl;
+      Assert (n==42, ExcInternalError()); };
+    X (const X&)  { deallog << "Copying" << std::endl; };
+    ~X () { deallog << "Destroying " << std::endl; };
+    int i;
+};
+
+Threads::ThreadLocalStorage<X> *tls_data;
+
+volatile int counter = 0;
+
+void execute (int i)
+{
+  tls_data->get().i = i;
+
+                                  // indicate that the TLS object has been
+                                  // accessed
+  static Threads::Mutex m;
+  {
+    Threads::Mutex::ScopedLock l(m);
+    ++counter;
+  }
+
+                                  // wait in order to make sure that the
+                                  // thread lives longer than the TLS object
+  sleep (5);
+}
+
+
+void test () 
+{
+                                  // create a thread local storage object
+  X exemplar(42);
+  tls_data = new Threads::ThreadLocalStorage<X>(exemplar);
+  
+                                  // create 5 threads and wait for their
+                                  // return. the OS will create 5 individual
+                                  // thread ids, which means that we will
+                                  // create 5 individual thread specific
+                                  // storage locations
+  Threads::ThreadGroup<> tg;
+  for (unsigned int i=10; i<15; ++i)
+    tg += Threads::new_thread (execute, i);
+
+                                  // spin lock until all threads have created
+                                  // their objects
+  while (counter != 5);
+  
+                                  // delete the TLS object. this should also
+                                  // destroy all the objects created so far,
+                                  // namely 5 of them, plus the one copy the
+                                  // TLS object stores; while that happens,
+                                  // we should get output into logstream
+  delete tls_data;
+
+                                  // write something into the output
+                                  // file. this also makes sure that the
+                                  // output file records the order in which
+                                  // this output is generated and in which
+                                  // the TLS objects were destroyed. we want
+                                  // that these objects are destroyed when
+                                  // the TLS object is destroyed, not when
+                                  // the threads exit
+  deallog << "Done." << std::endl;
+
+                                  // now make sure the threads all finish
+  tg.join_all ();
+
+                                  // at this point, the seventh object will
+                                  // be destroyed, which is the exemplar
+                                  // object local to this function
+}
+
+  
+  
+
+int main()
+{
+  std::ofstream logfile("thread_local_storage_02/output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  test ();
+}
diff --git a/tests/base/thread_local_storage_02/cmp/generic b/tests/base/thread_local_storage_02/cmp/generic
new file mode 100644 (file)
index 0000000..3b23f38
--- /dev/null
@@ -0,0 +1,16 @@
+
+DEAL::Creating
+DEAL::Copying
+DEAL::Copying
+DEAL::Copying
+DEAL::Copying
+DEAL::Copying
+DEAL::Copying
+DEAL::Destroying 
+DEAL::Destroying 
+DEAL::Destroying 
+DEAL::Destroying 
+DEAL::Destroying 
+DEAL::Destroying 
+DEAL::Done.
+DEAL::Destroying 

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.