From ca1792c1cd65c3b4d77e67a49967bf714112677f Mon Sep 17 00:00:00 2001 From: wolf Date: Wed, 13 Nov 2002 19:15:58 +0000 Subject: [PATCH] Optimize SparseMatrix in multithreaded mode, when only one thread is required. git-svn-id: https://svn.dealii.org/trunk@6760 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/doc/news/2002/c-3-4.html | 14 ++ .../lac/include/lac/sparse_matrix.templates.h | 169 ++++++++++-------- 2 files changed, 104 insertions(+), 79 deletions(-) diff --git a/deal.II/doc/news/2002/c-3-4.html b/deal.II/doc/news/2002/c-3-4.html index 08517d935a..eb1b40b48f 100644 --- a/deal.II/doc/news/2002/c-3-4.html +++ b/deal.II/doc/news/2002/c-3-4.html @@ -332,6 +332,20 @@ contributor's names are abbreviated by WB (Wolfgang Bangerth), GK

lac

    +
  1. + Changed: In multithread mode, the SparseMatrix would spawn + multithread_info.n_default_threads threads to + perform matrix-vector multiplications and similar + operations. It would even do so if + multithread_info.n_default_threads was equal to + one. In that case, we now do the operation on the thread we are + presently on, eliminating the overhead of spawning a single + thread, and later waiting and terminating it. +
    + (WB 2002/11/13) +

    +
  2. Fixed: In the SparseDirectMA27 class, wrapping the MA27 solver written in Fortran77 into some structure amenable to C++, diff --git a/deal.II/lac/include/lac/sparse_matrix.templates.h b/deal.II/lac/include/lac/sparse_matrix.templates.h index 9a1cca1697..689d935d5a 100644 --- a/deal.II/lac/include/lac/sparse_matrix.templates.h +++ b/deal.II/lac/include/lac/sparse_matrix.templates.h @@ -32,14 +32,11 @@ #include #include -#ifdef DEAL_II_USE_MT -# include -# include - -# include -# include -#endif +#include +#include +#include +#include @@ -288,11 +285,13 @@ SparseMatrix::vmult (Vector& dst, const unsigned int n_rows = m(); -#ifdef DEAL_II_USE_MT - // in MT mode: start new threads only - // if the matrix is sufficiently large. - // the limit is mostly artificial - if (n_rows/multithread_info.n_default_threads > 2000) + // in MT mode: start new threads + // only if the matrix is + // sufficiently large. the limit + // is mostly artificial + if (DEAL_II_USE_MT && + (multithread_info.n_default_threads > 1) && + (n_rows/multithread_info.n_default_threads > 2000)) { const unsigned int n_threads = multithread_info.n_default_threads; @@ -326,21 +325,22 @@ SparseMatrix::vmult (Vector& dst, thread_manager.wait (); return; - }; -#endif - - // if not in MT mode or size<2000 - // do it in an oldfashioned way - const number *val_ptr = &val[cols->rowstart[0]]; - const unsigned int *colnum_ptr = &cols->colnums[cols->rowstart[0]]; - somenumber *dst_ptr = &dst(0); - for (unsigned int row=0; rowrowstart[row+1]]; - while (val_ptr != val_end_of_row) - s += *val_ptr++ * src(*colnum_ptr++); - *dst_ptr++ = s; + // if not in MT mode or size<2000 + // do it in an oldfashioned way + const number *val_ptr = &val[cols->rowstart[0]]; + const unsigned int *colnum_ptr = &cols->colnums[cols->rowstart[0]]; + somenumber *dst_ptr = &dst(0); + for (unsigned int row=0; rowrowstart[row+1]]; + while (val_ptr != val_end_of_row) + s += *val_ptr++ * src(*colnum_ptr++); + *dst_ptr++ = s; + }; }; }; @@ -453,11 +453,13 @@ SparseMatrix::matrix_norm_square (const Vector& v) const Assert(n() == v.size(), ExcDimensionMismatch(n(),v.size())); const unsigned int n_rows = m(); -#ifdef DEAL_II_USE_MT + // if in MT mode and size sufficiently // large: do it in parallel; the limit // is mostly artificial - if (n_rows/multithread_info.n_default_threads > 2000) + if (DEAL_II_USE_MT && + (multithread_info.n_default_threads > 1) && + (n_rows/multithread_info.n_default_threads > 2000)) { const unsigned int n_threads = multithread_info.n_default_threads; @@ -499,24 +501,26 @@ SparseMatrix::matrix_norm_square (const Vector& v) const return std::accumulate (partial_sums.begin(), partial_sums.end(), 0.); - }; -#endif - // if not in MT mode or the matrix is - // too small: do it one-by-one - somenumber sum = 0.; - const number *val_ptr = &val[cols->rowstart[0]]; - const unsigned int *colnum_ptr = &cols->colnums[cols->rowstart[0]]; - for (unsigned int row=0; rowrowstart[row+1]]; - while (val_ptr != val_end_of_row) - s += *val_ptr++ * v(*colnum_ptr++); - - sum += s* v(row); + // if not in MT mode or the matrix is + // too small: do it one-by-one + somenumber sum = 0.; + const number *val_ptr = &val[cols->rowstart[0]]; + const unsigned int *colnum_ptr = &cols->colnums[cols->rowstart[0]]; + for (unsigned int row=0; rowrowstart[row+1]]; + while (val_ptr != val_end_of_row) + s += *val_ptr++ * v(*colnum_ptr++); + + sum += s* v(row); + }; + + return sum; }; - - return sum; }; @@ -564,11 +568,13 @@ SparseMatrix::matrix_scalar_product (const Vector& u, Assert(n() == v.size(), ExcDimensionMismatch(n(),v.size())); const unsigned int n_rows = m(); -#ifdef DEAL_II_USE_MT + // if in MT mode and size sufficiently // large: do it in parallel; the limit // is mostly artificial - if (n_rows/multithread_info.n_default_threads > 2000) + if (DEAL_II_USE_MT && + (multithread_info.n_default_threads != 1) && + (n_rows/multithread_info.n_default_threads > 2000)) { const unsigned int n_threads = multithread_info.n_default_threads; @@ -611,24 +617,26 @@ SparseMatrix::matrix_scalar_product (const Vector& u, return std::accumulate (partial_sums.begin(), partial_sums.end(), 0.); - }; -#endif - // if not in MT mode or the matrix is - // too small: do it one-by-one - somenumber sum = 0.; - const number *val_ptr = &val[cols->rowstart[0]]; - const unsigned int *colnum_ptr = &cols->colnums[cols->rowstart[0]]; - for (unsigned int row=0; rowrowstart[row+1]]; - while (val_ptr != val_end_of_row) - s += *val_ptr++ * v(*colnum_ptr++); - - sum += s* u(row); + // if not in MT mode or the matrix is + // too small: do it one-by-one + somenumber sum = 0.; + const number *val_ptr = &val[cols->rowstart[0]]; + const unsigned int *colnum_ptr = &cols->colnums[cols->rowstart[0]]; + for (unsigned int row=0; rowrowstart[row+1]]; + while (val_ptr != val_end_of_row) + s += *val_ptr++ * v(*colnum_ptr++); + + sum += s* u(row); + }; + + return sum; }; - - return sum; }; @@ -720,11 +728,13 @@ SparseMatrix::residual (Vector &dst, Assert(n() == u.size(), ExcDimensionMismatch(n(),u.size())); const unsigned int n_rows = m(); -#ifdef DEAL_II_USE_MT + // if in MT mode and size sufficiently // large: do it in parallel; the limit // is mostly artificial - if (n_rows/multithread_info.n_default_threads > 2000) + if (DEAL_II_USE_MT && + (multithread_info.n_default_threads > 1) && + (n_rows/multithread_info.n_default_threads > 2000)) { const unsigned int n_threads = multithread_info.n_default_threads; @@ -768,23 +778,24 @@ SparseMatrix::residual (Vector &dst, return std::sqrt(std::accumulate (partial_norms.begin(), partial_norms.end(), 0.)); - }; -#endif - - somenumber norm=0.; - - for (unsigned int i=0; irowstart[i]; jrowstart[i+1] ;j++) + somenumber norm=0.; + + for (unsigned int i=0; icolnums[j]; - s -= val[j] * u(p); + somenumber s = b(i); + for (unsigned int j=cols->rowstart[i]; jrowstart[i+1] ;j++) + { + const unsigned int p = cols->colnums[j]; + s -= val[j] * u(p); + } + dst(i) = s; + norm += dst(i)*dst(i); } - dst(i) = s; - norm += dst(i)*dst(i); - } - return std::sqrt(norm); + return std::sqrt(norm); + }; } @@ -797,10 +808,10 @@ SparseMatrix::threaded_residual (Vector &dst, const std::pair interval, somenumber *partial_norm) const { +#ifdef DEAL_II_USE_MT const unsigned int begin_row = interval.first, end_row = interval.second; -#ifdef DEAL_II_USE_MT somenumber norm=0.; for (unsigned int i=begin_row; i