From 0d2e59d13a3d59daf379e652b0a99512e80cb147 Mon Sep 17 00:00:00 2001 From: Denis Davydov Date: Tue, 26 Sep 2017 16:00:32 +0200 Subject: [PATCH] add thread mutex to const functions which use mutable variables --- include/deal.II/lac/lapack_full_matrix.h | 11 +++++++++++ source/lac/lapack_full_matrix.cc | 24 +++++++++++++++--------- 2 files changed, 26 insertions(+), 9 deletions(-) diff --git a/include/deal.II/lac/lapack_full_matrix.h b/include/deal.II/lac/lapack_full_matrix.h index 593794810b..6bf20a3a8f 100644 --- a/include/deal.II/lac/lapack_full_matrix.h +++ b/include/deal.II/lac/lapack_full_matrix.h @@ -22,6 +22,7 @@ #include #include #include +#include #include #include @@ -672,6 +673,11 @@ private: */ mutable std::vector work; + /** + * Integer working array used for some LAPACK functions. + */ + mutable std::vector iwork; + /** * The vector storing the permutations applied for pivoting in the LU- * factorization. @@ -717,6 +723,11 @@ private: * USVT. */ std::shared_ptr > svd_vt; + + /** + * Thread mutex. + */ + mutable Threads::Mutex mutex; }; diff --git a/source/lac/lapack_full_matrix.cc b/source/lac/lapack_full_matrix.cc index 763daac138..28961749ad 100644 --- a/source/lac/lapack_full_matrix.cc +++ b/source/lac/lapack_full_matrix.cc @@ -178,6 +178,8 @@ LAPACKFullMatrix::vmult ( const Vector &v, const bool adding) const { + Threads::Mutex::ScopedLock lock (mutex); + const int mm = this->n_rows(); const int nn = this->n_cols(); const number alpha = 1.; @@ -236,6 +238,8 @@ LAPACKFullMatrix::Tvmult ( const Vector &v, const bool adding) const { + Threads::Mutex::ScopedLock lock (mutex); + const int mm = this->n_rows(); const int nn = this->n_cols(); const number alpha = 1.; @@ -560,6 +564,8 @@ number LAPACKFullMatrix::frobenius_norm() const template number LAPACKFullMatrix::norm(const char type) const { + Threads::Mutex::ScopedLock lock (mutex); + Assert (state == LAPACKSupport::matrix || state == LAPACKSupport::inverse_matrix, ExcMessage("norms can be called in matrix state only.")); @@ -567,15 +573,14 @@ number LAPACKFullMatrix::norm(const char type) const const int N = this->n_cols(); const int M = this->n_rows(); const number *values = &this->values[0]; - std::vector w; if (property == symmetric) { const int lda = std::max(1,N); const int lwork = (type == 'I' || type == 'O') ? std::max(1,N) : 0; - w.resize(lwork); - return lansy (&type, &LAPACKSupport::L, &N, values, &lda, &w[0]); + work.resize(lwork); + return lansy (&type, &LAPACKSupport::L, &N, values, &lda, &work[0]); } else { @@ -583,8 +588,8 @@ number LAPACKFullMatrix::norm(const char type) const const int lwork = (type == 'I') ? std::max(1,M) : 0; - w.resize(lwork); - return lange (&type, &M, &N, values, &lda, &w[0]); + work.resize(lwork); + return lange (&type, &M, &N, values, &lda, &work[0]); } } @@ -619,6 +624,7 @@ template number LAPACKFullMatrix::reciprocal_condition_number(const number a_norm) const { + Threads::Mutex::ScopedLock lock (mutex); Assert(state == cholesky, ExcState(state)); number rcond = 0.; @@ -626,13 +632,13 @@ LAPACKFullMatrix::reciprocal_condition_number(const number a_norm) const const number *values = &this->values[0]; int info = 0; const int lda = std::max(1,N); - std::vector w(3*N); - std::vector iw(N); + work.resize(3*N); + iwork.resize(N); // use the same uplo as in Cholesky pocon (&LAPACKSupport::L, &N, values, &lda, &a_norm, &rcond, - &w[0], &iw[0], &info); + &work[0], &iwork[0], &info); Assert(info >= 0, ExcInternalError()); @@ -1002,7 +1008,7 @@ LAPACKFullMatrix::compute_generalized_eigenvalues_symmetric( const char *const uplo(&U); const char *const range(&V); const int *const dummy(&one); - std::vector iwork(static_cast (5*nn)); + iwork.resize(static_cast (5*nn)); std::vector ifail(static_cast (nn)); -- 2.39.5