From: wolf Date: Tue, 28 Jan 2003 21:34:05 +0000 (+0000) Subject: Implement scoped locking. Use it throughout the library. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=c69ababeba7d865d6ec72c1acdfc4de76868a602;p=dealii-svn.git Implement scoped locking. Use it throughout the library. git-svn-id: https://svn.dealii.org/trunk@6977 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/base/include/base/thread_management.h b/deal.II/base/include/base/thread_management.h index d7b4779dac..5d9eb31fd6 100644 --- a/deal.II/base/include/base/thread_management.h +++ b/deal.II/base/include/base/thread_management.h @@ -46,11 +46,64 @@ namespace Threads * reason to provide such a function is that the program can be * compiled both in MT and non-MT mode without difference. * - * @author Wolfgang Bangerth, 2000 + * @author Wolfgang Bangerth, 2000, 2003 */ class DummyThreadMutex { public: + /** + * Scoped lock class. When you + * declare an object of this + * type, you have to pass it a + * mutex, which is locked in + * the constructor of this + * class and unlocked in the + * destructor. The lock is thus + * held during the entire + * lifetime of this object, + * i.e. until the end of the + * present scope, which + * explains where the name + * comes from. This pattern of + * using locks with mutexes + * follows the + * resource-acquisition-is-initialization + * pattern, and was used first + * for mutexes by Doug + * Schmidt. It has the + * advantage that locking a + * mutex this way is + * thread-safe, i.e. when an + * exception is thrown between + * the locking and unlocking + * point, the destructor makes + * sure that the mutex is + * unlocked; this would not + * automatically be the case + * when you lock and unlock the + * mutex "by hand", i.e. using + * @p{acquire} and @p{release}. + */ + class ScopedLock + { + public: + /** + * Constructor. Lock the + * mutex. Since this is a + * dummy mutex class, this + * of course does nothing. + */ + ScopedLock (DummyThreadMutex &) {}; + + /** + * Destructor. Unlock the + * mutex. Since this is a + * dummy mutex class, this + * of course does nothing. + */ + ~ScopedLock () {}; + }; + /** * Simulate acquisition of the * mutex. As this class does @@ -247,11 +300,69 @@ namespace Threads * Class implementing a Mutex with * the help of POSIX functions. * - * @author Wolfgang Bangerth, 2002 + * @author Wolfgang Bangerth, 2002, 2003 */ class PosixThreadMutex { public: + /** + * Scoped lock class. When you + * declare an object of this + * type, you have to pass it a + * mutex, which is locked in + * the constructor of this + * class and unlocked in the + * destructor. The lock is thus + * held during the entire + * lifetime of this object, + * i.e. until the end of the + * present scope, which + * explains where the name + * comes from. This pattern of + * using locks with mutexes + * follows the + * resource-acquisition-is-initialization + * pattern, and was used first + * for mutexes by Doug + * Schmidt. It has the + * advantage that locking a + * mutex this way is + * thread-safe, i.e. when an + * exception is thrown between + * the locking and unlocking + * point, the destructor makes + * sure that the mutex is + * unlocked; this would not + * automatically be the case + * when you lock and unlock the + * mutex "by hand", i.e. using + * @p{acquire} and @p{release}. + */ + class ScopedLock + { + public: + /** + * Constructor. Lock the + * mutex. + */ + ScopedLock (PosixThreadMutex &m) : mutex(m) { mutex.acquire(); }; + + /** + * Destructor. Unlock the + * mutex. Since this is a + * dummy mutex class, this + * of course does nothing. + */ + ~ScopedLock () { mutex.release (); }; + + private: + /** + * Store the address of the + * mutex object. + */ + PosixThreadMutex &mutex; + }; + /** * Constructor. Initialize the * underlying POSIX mutex data diff --git a/deal.II/base/source/polynomial.cc b/deal.II/base/source/polynomial.cc index 1796d40e4f..2598d7addf 100644 --- a/deal.II/base/source/polynomial.cc +++ b/deal.II/base/source/polynomial.cc @@ -540,9 +540,11 @@ namespace Polynomials unsigned int k = k_; // first make sure that no other - // thread intercepts the operation - // of this function - coefficients_lock.acquire (); + // thread intercepts the + // operation of this function; + // for this, acquire the lock + // until we quit this function + Threads::ThreadMutex::ScopedLock lock(coefficients_lock); // The first 2 coefficients are hard-coded if (k==0) @@ -633,10 +635,6 @@ namespace Polynomials shifted_coefficients[k] = ck; }; }; - - // now, everything is done, so - // release the lock again - coefficients_lock.release (); } @@ -652,15 +650,8 @@ namespace Polynomials // then get a pointer to the array // of coefficients. do that in a MT // safe way - coefficients_lock.acquire (); - const std::vector *p = shifted_coefficients[k]; - coefficients_lock.release (); - - // return the object pointed - // to. since this object does not - // change any more once computed, - // this is MT safe - return *p; + Threads::ThreadMutex::ScopedLock lock (coefficients_lock); + return *shifted_coefficients[k]; } @@ -705,34 +696,36 @@ namespace Polynomials { unsigned int k = k_; - // first make sure that no other - // thread intercepts the operation - // of this function - coefficients_lock.acquire (); + // first make sure that no other + // thread intercepts the operation + // of this function + // for this, acquire the lock + // until we quit this function + Threads::ThreadMutex::ScopedLock lock(coefficients_lock); - // The first 2 coefficients - // are hard-coded + // The first 2 coefficients + // are hard-coded if (k==0) k=1; - // check: does the information - // already exist? + // check: does the information + // already exist? if ( (recursive_coefficients.size() < k+1) || ((recursive_coefficients.size() >= k+1) && (recursive_coefficients[k] == 0)) ) - // no, then generate the - // respective coefficients + // no, then generate the + // respective coefficients { recursive_coefficients.resize (k+1, 0); if (k<=1) { - // create coefficients - // vectors for k=0 and k=1 - // - // allocate the respective - // amount of memory and - // later assign it to the - // coefficients array to - // make it const + // create coefficients + // vectors for k=0 and k=1 + // + // allocate the respective + // amount of memory and + // later assign it to the + // coefficients array to + // make it const std::vector *c0 = new std::vector(2); (*c0)[0] = 1.0; (*c0)[1] = -1.0; @@ -741,8 +734,8 @@ namespace Polynomials (*c1)[0] = 0.0; (*c1)[1] = 1.0; - // now make these arrays - // const + // now make these arrays + // const recursive_coefficients[0] = c0; recursive_coefficients[1] = c1; } @@ -760,14 +753,14 @@ namespace Polynomials } else { - // for larger numbers, - // compute the coefficients - // recursively. to do so, - // we have to release the - // lock temporarily to - // allow the called - // function to acquire it - // itself + // for larger numbers, + // compute the coefficients + // recursively. to do so, + // we have to release the + // lock temporarily to + // allow the called + // function to acquire it + // itself coefficients_lock.release (); compute_coefficients(k-1); coefficients_lock.acquire (); @@ -781,25 +774,21 @@ namespace Polynomials - (*recursive_coefficients[k-1])[i] ); (*ck)[k] = 2.*(*recursive_coefficients[k-1])[k-1]; - // for even degrees, we need - // to add a multiple of - // basis fcn phi_2 + // for even degrees, we need + // to add a multiple of + // basis fcn phi_2 if ( (k%2) == 0 ) { (*ck)[1] += (*recursive_coefficients[2])[1]; (*ck)[2] += (*recursive_coefficients[2])[2]; } - // finally assign the newly - // created vector to the - // const pointer in the - // coefficients array + // finally assign the newly + // created vector to the + // const pointer in the + // coefficients array recursive_coefficients[k] = ck; }; }; - - // now, everything is done, so - // release the lock again - coefficients_lock.release (); } diff --git a/deal.II/base/source/subscriptor.cc b/deal.II/base/source/subscriptor.cc index c63b85d565..5216286251 100644 --- a/deal.II/base/source/subscriptor.cc +++ b/deal.II/base/source/subscriptor.cc @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 1998, 1999, 2000, 2001, 2002 by the deal authors +// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003 by the deal authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -129,17 +129,15 @@ void Subscriptor::subscribe () const object_info = &typeid(*this); #endif - subscription_lock.acquire(); + Threads::ThreadMutex::ScopedLock lock (subscription_lock); ++counter; - subscription_lock.release(); } void Subscriptor::unsubscribe () const { Assert (counter>0, ExcNotUsed()); - subscription_lock.acquire(); + Threads::ThreadMutex::ScopedLock lock (subscription_lock); --counter; - subscription_lock.release(); } diff --git a/deal.II/base/source/thread_management.cc b/deal.II/base/source/thread_management.cc index 8ac3e55497..718234dda7 100644 --- a/deal.II/base/source/thread_management.cc +++ b/deal.II/base/source/thread_management.cc @@ -37,20 +37,18 @@ namespace Threads void register_new_thread () { - n_existing_threads_mutex.acquire (); + ThreadMutex::ScopedLock lock (n_existing_threads_mutex); ++n_existing_threads_counter; - n_existing_threads_mutex.release (); } void deregister_new_thread () { - n_existing_threads_mutex.acquire (); + ThreadMutex::ScopedLock lock (n_existing_threads_mutex); --n_existing_threads_counter; Assert (n_existing_threads_counter >= 1, ExcInternalError()); - n_existing_threads_mutex.release (); } @@ -100,9 +98,8 @@ namespace Threads unsigned int n_existing_threads () { - n_existing_threads_mutex.acquire (); + ThreadMutex::ScopedLock lock (n_existing_threads_mutex); const unsigned int n = n_existing_threads_counter; - n_existing_threads_mutex.release (); return n; } @@ -233,10 +230,9 @@ namespace Threads // wait for all threads, and // release memory wait (); - list_mutex.acquire (); + ThreadMutex::ScopedLock lock (list_mutex); if (thread_id_list != 0) delete reinterpret_cast*>(thread_id_list); - list_mutex.release (); } @@ -248,12 +244,13 @@ namespace Threads { std::list &tid_list = *reinterpret_cast*>(thread_id_list); - - list_mutex.acquire (); - tid_list.push_back (pthread_t()); - pthread_t *tid = &tid_list.back(); - list_mutex.release (); - + + { + ThreadMutex::ScopedLock lock (list_mutex); + tid_list.push_back (pthread_t()); + pthread_t *tid = &tid_list.back(); + } + // start new thread. retry until // we either succeed or get an // error other than EAGAIN @@ -272,7 +269,7 @@ namespace Threads void PosixThreadManager::wait () const { - list_mutex.acquire (); + ThreadMutex::ScopedLock lock (list_mutex); std::list &tid_list = *reinterpret_cast*>(thread_id_list); @@ -290,8 +287,6 @@ namespace Threads // for expired threads with their // invalid handles again tid_list.clear (); - - list_mutex.release (); } # endif diff --git a/deal.II/doc/news/2002/c-3-4.html b/deal.II/doc/news/2002/c-3-4.html index 6c750d226c..9eef20c7cd 100644 --- a/deal.II/doc/news/2002/c-3-4.html +++ b/deal.II/doc/news/2002/c-3-4.html @@ -615,6 +615,15 @@ contributor's names are abbreviated by WB (Wolfgang Bangerth), GK

deal.II

    +
  1. + New: The ThreadMutex classes now have a + member class ScopedLock that implements the + scoped thread-safe locking pattern of Doug Schmidt. It is also used in + various places of the code now. +
    + (WB 2003/01/28) +

    +
  2. Fixed: The GridReordering tried to be thread-safe in the initialization of some data, but was not due to a diff --git a/deal.II/examples/step-13/step-13.cc b/deal.II/examples/step-13/step-13.cc index 4d03c9f8d3..d550b40139 100644 --- a/deal.II/examples/step-13/step-13.cc +++ b/deal.II/examples/step-13/step-13.cc @@ -1268,7 +1268,8 @@ namespace LaplaceSolver // the domain. Since it's actions // have all been explained in // previous programs, we do not - // comment on it any more. + // comment on it any more, except + // for one pointe below. template void Solver::assemble_matrix (LinearSystem &linear_system, @@ -1303,13 +1304,99 @@ namespace LaplaceSolver cell->get_dof_indices (local_dof_indices); - mutex.acquire (); + + // In the step-9 program, we + // have shown that you have + // to use the mutex to lock + // the matrix when copying + // the elements from the + // local to the global + // matrix. This was necessary + // to avoid that two threads + // access it at the same + // time, eventually + // overwriting their + // respective + // work. Previously, we have + // used the `acquire'' and + // ``release'' functions of + // the mutex to lock and + // unlock the mutex, + // respectively. While this + // is valid, there is one + // possible catch: if between + // the locking operation and + // the unlocking operation an + // exception is thrown, the + // mutex remains in the + // locked state, and in some + // cases this might lead to + // deadlocks. A similar + // situation arises, when one + // changes the code to have a + // return statement somewhere + // in the middle of the + // locked block, and forgets + // that before we call + // ``return'', we also have + // to unlock the mutex. This + // all is not be a problem + // here, but we want to show + // the general technique to + // cope with these problems + // nevertheless: have an + // object that upon + // initialization (i.e. in + // its constructor) locks the + // mutex, and on running the + // destructor unlocks it + // again. This is called the + // ``scoped lock'' pattern + // (apparently invented by + // Doug Schmidt originally), + // and it works because + // destructors of local + // objects are also run when + // we exit the function + // either through a + // ``return'' statement, or + // when an exception is + // raised. Thus, it is + // guaranteed that the mutex + // will always be unlocked + // when we exit this part of + // the program, whether the + // operation completed + // successfully or not, + // whether the exit path was + // something we implemented + // willfully or whether the + // function was exited by an + // exception that we did not + // forsee. + // + // deal.II implements the + // scoped locking pattern in + // the + // ThreadMutex::ScopedLock + // class: it takes the mutex + // in the constructor and + // locks it; in its + // destructor, it unlocks it + // again. So here is how it + // is used: + Threads::ThreadMutex::ScopedLock lock (mutex); for (unsigned int i=0; iget_dof_indices (local_dof_indices); - mutex.acquire (); + Threads::ThreadMutex::ScopedLock lock (mutex); for (unsigned int i=0; i::operator= (const SwappableVector &v) // if in MT mode, block all other // operations. if not in MT mode, // this is a no-op - lock.acquire (); + Threads::ThreadMutex::ScopedLock lock(this->lock); Vector::operator = (v); data_is_preloaded = false; - lock.release (); - return *this; } @@ -103,7 +101,7 @@ void SwappableVector::swap_out (const std::string &name) // if in MT mode, block all other // operations. if not in MT mode, // this is a no-op - lock.acquire (); + Threads::ThreadMutex::ScopedLock lock(this->lock); // check that we have not called // @p{alert} without the respective @@ -115,8 +113,6 @@ void SwappableVector::swap_out (const std::string &name) tmp_out.close (); this->reinit (0); - - lock.release (); } @@ -218,7 +214,7 @@ void SwappableVector::kill_file () // (there should be none, but who // knows). if not in MT mode, // this is a no-op - lock.acquire(); + Threads::ThreadMutex::ScopedLock lock(this->lock); // this is too bad: someone // requested the vector in advance, @@ -233,8 +229,6 @@ void SwappableVector::kill_file () filename = ""; }; - - lock.release (); } diff --git a/deal.II/lac/include/lac/vector_memory.h b/deal.II/lac/include/lac/vector_memory.h index 1576ae8ded..4f117b7ed9 100644 --- a/deal.II/lac/include/lac/vector_memory.h +++ b/deal.II/lac/include/lac/vector_memory.h @@ -184,7 +184,7 @@ template GrowingVectorMemory::GrowingVectorMemory(const unsigned int initial_size) : pool(initial_size) { - mutex.acquire (); + Threads::ThreadMutex::ScopedLock lock(mutex); for (typename std::vector::iterator i=pool.begin(); i != pool.end(); ++i) @@ -195,7 +195,6 @@ GrowingVectorMemory::GrowingVectorMemory(const unsigned int initial_size // no vectors yet claimed n_alloc = 0; - mutex.release (); } @@ -203,7 +202,7 @@ GrowingVectorMemory::GrowingVectorMemory(const unsigned int initial_size template GrowingVectorMemory::~GrowingVectorMemory() { - mutex.acquire (); + Threads::ThreadMutex::ScopedLock lock(mutex); // deallocate all vectors and count // number of vectors that is still @@ -222,7 +221,6 @@ GrowingVectorMemory::~GrowingVectorMemory() deallog << "GrowingVectorMemory:Maximum allocated vectors: " << pool.size() << std::endl; pool.clear (); - mutex.release (); // write out warning if memory leak if (n!=0) @@ -233,10 +231,10 @@ GrowingVectorMemory::~GrowingVectorMemory() template -Vector* +Vector * GrowingVectorMemory::alloc() { - mutex.acquire (); + Threads::ThreadMutex::ScopedLock lock(mutex); ++n_alloc; for (typename std::vector::iterator i=pool.begin(); i != pool.end(); @@ -245,15 +243,12 @@ GrowingVectorMemory::alloc() if (i->first == false) { i->first = true; - mutex.release (); return (i->second); } } - entry_type t; - t.first = true; - t.second = new Vector; + + const entry_type t (true, new Vector); pool.push_back(t); - mutex.release (); return t.second; } @@ -264,18 +259,15 @@ template void GrowingVectorMemory::free(const Vector* const v) { - mutex.acquire (); + Threads::ThreadMutex::ScopedLock lock(mutex); for (typename std::vector::iterator i=pool.begin();i != pool.end() ;++i) { if (v == (i->second)) { i->first = false; - mutex.release (); return; } } - mutex.release (); - Assert(false, typename VectorMemory::ExcNotAllocatedHere()); } diff --git a/deal.II/lac/source/filtered_matrix.cc b/deal.II/lac/source/filtered_matrix.cc index c937618cc3..6923ccafec 100644 --- a/deal.II/lac/source/filtered_matrix.cc +++ b/deal.II/lac/source/filtered_matrix.cc @@ -93,9 +93,8 @@ void FilteredMatrix,Vector >:: allocate_tmp_vector () { - tmp_mutex.acquire (); + Threads::ThreadMutex::ScopedLock lock (tmp_mutex); tmp_vector.reinit (matrix->n(), true); - tmp_mutex.release (); } @@ -105,9 +104,8 @@ void FilteredMatrix,Vector >:: allocate_tmp_vector () { - tmp_mutex.acquire (); + Threads::ThreadMutex::ScopedLock lock (tmp_mutex); tmp_vector.reinit (matrix->n(), true); - tmp_mutex.release (); } @@ -121,9 +119,8 @@ allocate_tmp_vector () for (unsigned int i=0; iblock(i,i).n(); - tmp_mutex.acquire (); + Threads::ThreadMutex::ScopedLock lock (tmp_mutex); tmp_vector.reinit (block_sizes, true); - tmp_mutex.release (); } @@ -137,9 +134,8 @@ allocate_tmp_vector () for (unsigned int i=0; iblock(i,i).n(); - tmp_mutex.acquire (); + Threads::ThreadMutex::ScopedLock lock (tmp_mutex); tmp_vector.reinit (block_sizes, true); - tmp_mutex.release (); } diff --git a/deal.II/lac/source/sparse_direct.cc b/deal.II/lac/source/sparse_direct.cc index 82180008a2..be91d8d161 100644 --- a/deal.II/lac/source/sparse_direct.cc +++ b/deal.II/lac/source/sparse_direct.cc @@ -261,9 +261,8 @@ namespace CommunicationsLog { Record record = {child_pid, direction, &typeid(T), count, sizeof(T)*count, completed_bytes, descr}; - list_access_lock.acquire (); + Threads::ThreadMutex::ScopedLock lock (list_access_lock); communication_log.push_back (record); - list_access_lock.release (); } @@ -274,7 +273,7 @@ namespace CommunicationsLog */ void list_communication (const pid_t /*child_pid*/) { - list_access_lock.acquire (); + Threads::ThreadMutex::ScopedLock lock (list_access_lock); std::cerr << "++++++++++++++++++++++++++++++" << std::endl << "Communication log history:" << std::endl; @@ -292,8 +291,6 @@ namespace CommunicationsLog << i->description << std::endl; std::cerr << "++++++++++++++++++++++++++++++" << std::endl; - - list_access_lock.release (); } @@ -308,12 +305,11 @@ namespace CommunicationsLog */ void remove_history (const pid_t &child_pid) { - list_access_lock.acquire (); + Threads::ThreadMutex::ScopedLock lock (list_access_lock); for (std::list::iterator i=communication_log.begin(); i!=communication_log.end(); ++i) if (i->child_pid == child_pid) communication_log.erase (i); - list_access_lock.release (); } } @@ -494,9 +490,8 @@ SparseDirectMA27::~SparseDirectMA27() CommunicationsLog::remove_history (detached_mode_data->child_pid); // next close down client - detached_mode_data->mutex.acquire (); + Threads::ThreadMutex::ScopedLock lock (detached_mode_data->mutex); write (detached_mode_data->server_client_pipe[1], "7", 1); - detached_mode_data->mutex.release (); // then also delete data delete detached_mode_data; @@ -949,7 +944,7 @@ void SparseDirectMA27::call_ma27ad (const unsigned int *N, IKEEP, IW1, NSTEPS, IFLAG); else { - detached_mode_data->mutex.acquire (); + Threads::ThreadMutex::ScopedLock lock (detached_mode_data->mutex); // first write the data we have // to push over, i.e. first // function index, then array @@ -972,8 +967,6 @@ void SparseDirectMA27::call_ma27ad (const unsigned int *N, // next get back what we need // to know detached_mode_data->get (IFLAG, 1, "IFLAG"); - - detached_mode_data->mutex.release (); }; } @@ -1001,7 +994,7 @@ void SparseDirectMA27::call_ma27bd (const unsigned int *N, // basically, everything is // already over the line, // except for A and LA - detached_mode_data->mutex.acquire (); + Threads::ThreadMutex::ScopedLock lock (detached_mode_data->mutex); detached_mode_data->put ("2", 1, "ACTION 2"); detached_mode_data->put (LA, 1, "LA"); @@ -1010,8 +1003,6 @@ void SparseDirectMA27::call_ma27bd (const unsigned int *N, // next get back what we need // to know detached_mode_data->get (IFLAG, 1, "IFLAG"); - - detached_mode_data->mutex.release (); }; } @@ -1051,13 +1042,12 @@ void SparseDirectMA27::call_ma27x1 (unsigned int *NRLNEC) HSL::MA27::ma27x1_ (NRLNEC); else { - detached_mode_data->mutex.acquire (); + Threads::ThreadMutex::ScopedLock lock (detached_mode_data->mutex); // ma27x1 only reads data, so // don't send anything except // for the id detached_mode_data->put ("4", 1, "ACTION 4"); detached_mode_data->get (NRLNEC, 1, "NRLNEC"); - detached_mode_data->mutex.release (); }; } @@ -1069,13 +1059,12 @@ void SparseDirectMA27::call_ma27x2 (unsigned int *NIRNEC) HSL::MA27::ma27x2_ (NIRNEC); else { - detached_mode_data->mutex.acquire (); + Threads::ThreadMutex::ScopedLock lock (detached_mode_data->mutex); // ma27x2 only reads data, so // don't send anything except // for the id detached_mode_data->put ("5", 1, "ACTION 5"); detached_mode_data->get (NIRNEC, 1, "NIRNEC"); - detached_mode_data->mutex.release (); }; } @@ -1087,13 +1076,12 @@ void SparseDirectMA27::call_ma27x3 (const unsigned int *LP) HSL::MA27::ma27x3_ (LP); else { - detached_mode_data->mutex.acquire (); + Threads::ThreadMutex::ScopedLock lock (detached_mode_data->mutex); // ma27x2 only reads data, so // don't send anything except // for the id detached_mode_data->put ("6", 1, "ACTION 6"); detached_mode_data->put (LP, 1, "LP"); - detached_mode_data->mutex.release (); }; }