<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
--- /dev/null
+//---------------------------------------------------------------------------
+// $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
+//---------------------------------------------------------------------------
--- /dev/null
+//-----------------------------------------------------------------------------
+// $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 ();
+}
--- /dev/null
+
+DEAL::Creating
+DEAL::Creating
+DEAL::Creating
+DEAL::Creating
+DEAL::Creating
+DEAL::Destroying
+DEAL::Destroying
+DEAL::Destroying
+DEAL::Destroying
+DEAL::Destroying
+DEAL::Done.
--- /dev/null
+//-----------------------------------------------------------------------------
+// $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 ();
+}
--- /dev/null
+
+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