From 4e99a5294bda28e99ac7ba1acc9ea4519b68ab86 Mon Sep 17 00:00:00 2001 From: kronbichler Date: Thu, 27 Nov 2008 17:13:50 +0000 Subject: [PATCH] ConstraintMatrix::distribute_local_to_global() uses now collective add operations into sparse matrices. git-svn-id: https://svn.dealii.org/trunk@17765 0785d39b-7218-0410-832d-ea1e28bc413d --- .../deal.II/include/dofs/dof_constraints.h | 32 ++- .../include/dofs/dof_constraints.templates.h | 217 ++++++++++------ deal.II/doc/news/changes.h | 20 ++ deal.II/lac/include/lac/block_matrix_base.h | 243 ++++++++++-------- deal.II/lac/include/lac/sparse_matrix.h | 52 ++-- deal.II/lac/source/trilinos_sparse_matrix.cc | 5 +- 6 files changed, 356 insertions(+), 213 deletions(-) diff --git a/deal.II/deal.II/include/dofs/dof_constraints.h b/deal.II/deal.II/include/dofs/dof_constraints.h index 06ec28b472..197e7b1932 100644 --- a/deal.II/deal.II/include/dofs/dof_constraints.h +++ b/deal.II/deal.II/include/dofs/dof_constraints.h @@ -18,6 +18,7 @@ #include #include #include +#include #include #include @@ -269,7 +270,7 @@ class BlockIndices; * a second step. * * @ingroup dofs - * @author Wolfgang Bangerth, 1998, 2004 + * @author Wolfgang Bangerth, 1998, 2004, 2008 */ class ConstraintMatrix : public Subscriptor { @@ -1353,6 +1354,35 @@ class ConstraintMatrix : public Subscriptor * add_entries_local_to_global(). */ static const Table<2,bool> default_empty_table; + + /** + * A mutex that makes the function + * distribute_local_to_global() + * only accessible to one process at + * time. + */ + Threads::ThreadMutex mutex; + + /** + * A temporary array for the function + * distribute_local_to_global() + * that collects column indices for + * values that should be written into + * the same row (so that we can write + * the values collectively into the + * sparse matrix, which is faster). + */ + mutable std::vector column_indices; + + /** + * A temporary array for the function + * distribute_local_to_global() + * that collects values that should be + * written into the same row (so that + * we can write them collectively into + * the sparse matrix, which is faster). + */ + mutable std::vector column_values; }; diff --git a/deal.II/deal.II/include/dofs/dof_constraints.templates.h b/deal.II/deal.II/include/dofs/dof_constraints.templates.h index c9b41d4762..22f1ba11b0 100644 --- a/deal.II/deal.II/include/dofs/dof_constraints.templates.h +++ b/deal.II/deal.II/include/dofs/dof_constraints.templates.h @@ -101,8 +101,8 @@ ConstraintMatrix::condense (const SparseMatrix &uncondensed, uncondensed.global_entry(j)); else { - // let c point to the constraint - // of this column + // let c point to the + // constraint of this column std::vector::const_iterator c = lines.begin(); while (c->line != uncondensed_struct.get_column_numbers()[j]) ++c; @@ -649,35 +649,39 @@ distribute_local_to_global (const FullMatrix &local_matrix, Assert (sorted == true, ExcMatrixNotClosed()); const unsigned int n_local_dofs = local_dof_indices.size(); + + // A lock that allows only one thread at + // time to go on in this function. + mutex.acquire(); // have a special case where there are no // constraints at all, since then we can be // a lot faster if (lines.size() == 0) - { - for (unsigned int i=0; i constraint_lines (n_local_dofs, static_cast(0)); + unsigned int n_max_entries_per_row = 0; for (unsigned int i=0; i &local_matrix, lines.end(), index_comparison); - // if this dof is constrained, then - // set the respective entry in the - // array. otherwise leave it at the - // invalid position + // if this dof is constrained, + // then set the respective entry + // in the array. otherwise leave + // it at the invalid position if ((position != lines.end()) && (position->line == local_dof_indices[i])) - constraint_lines[i] = &*position; + { + constraint_lines[i] = &*position; + n_max_entries_per_row += position->entries.size(); + } } + // We need to add the number of + // entries in the local matrix in + // order to obtain a sufficient size + // for the scratch array. + n_max_entries_per_row += n_local_dofs; + column_indices.resize(n_max_entries_per_row); + column_values.resize(n_max_entries_per_row); - // now distribute entries + // now distribute entries row by row for (unsigned int i=0; i &local_matrix, (is_constrained_j == false)) { // neither row nor column - // is constrained - global_matrix.add (local_dof_indices[i], - local_dof_indices[j], - local_matrix(i,j)); + // is constrained, so + // write the value into + // the scratch array + column_indices[col_counter] = local_dof_indices[j]; + column_values[col_counter] = local_matrix(i,j); + col_counter++; } else if ((is_constrained_i == true) && (is_constrained_j == false)) { - // ok, row is constrained, - // but column is not + // ok, row is + // constrained, but + // column is not. This + // creates entries in + // several rows to the + // same column, which is + // not covered by the + // scratch array. Write + // the values directly + // into the matrix for (unsigned int q=0; qentries.size(); ++q) global_matrix.add (position_i->entries[q].first, local_dof_indices[j], @@ -733,19 +764,32 @@ distribute_local_to_global (const FullMatrix &local_matrix, (is_constrained_j == true)) { // simply the other way - // round: row ok, column is - // constrained + // round: row ok, column + // is constrained. This + // time, we can put + // everything into the + // scratch array, since + // we are in the correct + // row. for (unsigned int q=0; qentries.size(); ++q) - global_matrix.add (local_dof_indices[i], - position_j->entries[q].first, - local_matrix(i,j) * - position_j->entries[q].second); + { + column_indices[col_counter] = position_j->entries[q].first; + column_values[col_counter] = local_matrix(i,j) * + position_j->entries[q].second; + col_counter++; + } } else if ((is_constrained_i == true) && (is_constrained_j == true)) { - // last case: both row and - // column are constrained + // last case: both row + // and column are + // constrained. Again, + // this creates entries + // in other rows than the + // current one, so write + // the values again in + // the matrix directly for (unsigned int p=0; pentries.size(); ++p) for (unsigned int q=0; qentries.size(); ++q) global_matrix.add (position_i->entries[p].first, @@ -760,61 +804,68 @@ distribute_local_to_global (const FullMatrix &local_matrix, // do something with the // diagonal elements. add // the absolute value of - // the local matrix, so the - // resulting entry will - // always be positive and + // the local matrix, so + // the resulting entry + // will always be + // positive and // furthermore be in the - // same order of magnitude - // as the other elements of - // the matrix + // same order of + // magnitude as the other + // elements of the matrix // - // note that this - // also captures - // the special case - // that a dof is - // both constrained - // and fixed (this - // can happen for + // note that this also + // captures the special + // case that a dof is + // both constrained and + // fixed (this can happen + // for hanging nodes in + // 3d that also happen to + // be on the + // boundary). in that + // case, following the + // above program flow, it + // is realized that when + // distributing the row + // and column no elements + // of the matrix are + // actually touched if + // all the degrees of + // freedom to which this + // dof is constrained are + // also constrained (the + // usual case with // hanging nodes in - // 3d that also - // happen to be on - // the - // boundary). in - // that case, - // following the - // above program - // flow, it is - // realized that - // when - // distributing the - // row and column - // no elements of - // the matrix are - // actually touched - // if all the - // degrees of - // freedom to which - // this dof is - // constrained are - // also constrained - // (the usual case - // with hanging - // nodes in - // 3d). however, in - // the line below, - // we do actually - // do something + // 3d). however, in the + // line below, we do + // actually do something // with this dof if (i == j) - global_matrix.add (local_dof_indices[i], - local_dof_indices[i], - local_matrix(i,i)); + { + column_indices[col_counter] = local_dof_indices[j]; + column_values[col_counter] = local_matrix(i,j); + col_counter++; + } } else Assert (false, ExcInternalError()); } + + // Check whether we did remain within the + // arrays when adding elements into the + // scratch arrays. Moreover, there should + // be at least one element in the scratch + // array (the element diagonal). + Assert (col_counter <= n_max_entries_per_row, ExcInternalError()); + + // Finally, write the scratch array into + // the sparse matrix. + if (col_counter > 0) + global_matrix.add(local_dof_indices[i], col_counter, + &column_indices[0], &column_values[0], + false); } } + mutex.release(); } diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index 5764cfac01..c7f2f14c1d 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -374,6 +374,17 @@ inconvenience this causes.

lac

    +
  1. +

    + New: All SparseMatrix classes (SparseMatrix, PETSc sparse + matrices, Trilinos sparse matrices, block sparse matrices) can now + directly add and set a FullMatrix and some other arrays into their value + list. This is faster and more convenient than an element-by-element + addition/set. +
    + (Martin Kronbichler 2008/11/26) +

    +
  2. New: The class LAPACKFullMatrix can now invert full matrices using @@ -492,6 +503,15 @@ inconvenience this causes.

    deal.II

      +
    1. +

      + Updated: The function ConstraintMatrix::distribute_local_to_global() for + matrices does now use row-wise addition into sparse matrices, which + accelerates the transfer from local to global data. +
      + (Martin Kronbichler, WB 2008/11/27) +

      +
    2. Fixed: The VectorTools::interpolate_boundary_values function was implemented a bit diff --git a/deal.II/lac/include/lac/block_matrix_base.h b/deal.II/lac/include/lac/block_matrix_base.h index b408f40742..e31ee10ad0 100644 --- a/deal.II/lac/include/lac/block_matrix_base.h +++ b/deal.II/lac/include/lac/block_matrix_base.h @@ -506,8 +506,9 @@ class BlockMatrixBase : public Subscriptor * false, i.e., even zero * values are treated. */ + template void set (const std::vector &indices, - const FullMatrix &full_matrix, + const FullMatrix &full_matrix, const bool elide_zero_values = false); /** @@ -517,9 +518,10 @@ class BlockMatrixBase : public Subscriptor * different local-to-global indexing * on rows and columns, respectively. */ + template void set (const std::vector &row_indices, const std::vector &col_indices, - const FullMatrix &full_matrix, + const FullMatrix &full_matrix, const bool elide_zero_values = false); /** @@ -540,9 +542,10 @@ class BlockMatrixBase : public Subscriptor * false, i.e., even zero * values are treated. */ + template void set (const unsigned int row, const std::vector &col_indices, - const std::vector &values, + const std::vector &values, const bool elide_zero_values = false); /** @@ -561,10 +564,11 @@ class BlockMatrixBase : public Subscriptor * false, i.e., even zero * values are inserted/replaced. */ + template void set (const unsigned int row, const unsigned int n_cols, const unsigned int *col_indices, - const value_type *values, + const number *values, const bool elide_zero_values = false); /** @@ -607,8 +611,9 @@ class BlockMatrixBase : public Subscriptor * i.e., zero values won't be added * into the matrix. */ + template void add (const std::vector &indices, - const FullMatrix &full_matrix, + const FullMatrix &full_matrix, const bool elide_zero_values = true); /** @@ -618,9 +623,10 @@ class BlockMatrixBase : public Subscriptor * different local-to-global indexing * on rows and columns, respectively. */ + template void add (const std::vector &row_indices, const std::vector &col_indices, - const FullMatrix &full_matrix, + const FullMatrix &full_matrix, const bool elide_zero_values = true); /** @@ -640,9 +646,10 @@ class BlockMatrixBase : public Subscriptor * i.e., zero values won't be added * into the matrix. */ + template void add (const unsigned int row, const std::vector &col_indices, - const std::vector &values, + const std::vector &values, const bool elide_zero_values = true); /** @@ -662,10 +669,11 @@ class BlockMatrixBase : public Subscriptor * i.e., zero values won't be added * into the matrix. */ + template void add (const unsigned int row, const unsigned int n_cols, const unsigned int *col_indices, - const value_type *values, + const number *values, const bool elide_zero_values = true); /** @@ -1165,24 +1173,28 @@ class BlockMatrixBase : public Subscriptor private: /** - * Temporary vector for block indices - * when writing local to global data. + * Temporary vector for counting the + * elements written into the + * individual blocks when doing a + * collective add or set. */ - std::vector block_col_indices; + std::vector counter_within_block; /** * Temporary vector for column - * indices when writing local to - * global data and passing the call - * down to a SparseMatrix. + * indices on each block when writing + * local to global data on each + * sparse matrix. */ - std::vector local_col_indices; + std::vector > column_indices; /** * Temporary vector for storing the - * local lengths of the row data. + * local values (they need to be + * reordered when writing local to + * global). */ - std::vector local_row_length; + std::vector > column_values; }; @@ -1778,11 +1790,12 @@ BlockMatrixBase::set (const unsigned int i, template +template inline void BlockMatrixBase::set (const std::vector &row_indices, const std::vector &col_indices, - const FullMatrix &values, + const FullMatrix &values, const bool elide_zero_values) { Assert (row_indices.size() == values.m(), @@ -1798,10 +1811,11 @@ BlockMatrixBase::set (const std::vector &row_indices, template +template inline void BlockMatrixBase::set (const std::vector &indices, - const FullMatrix &values, + const FullMatrix &values, const bool elide_zero_values) { Assert (indices.size() == values.m(), @@ -1816,11 +1830,12 @@ BlockMatrixBase::set (const std::vector &indices, template +template inline void BlockMatrixBase::set (const unsigned int row, const std::vector &col_indices, - const std::vector &values, + const std::vector &values, const bool elide_zero_values) { Assert (col_indices.size() == values.size(), @@ -1836,61 +1851,64 @@ BlockMatrixBase::set (const unsigned int row, // we need to calculate to each position // the location in the global array. template +template inline void -BlockMatrixBase::set (const unsigned int row, - const unsigned int n_cols, - const unsigned int *col_indices, - const value_type *values, - const bool elide_zero_values) +BlockMatrixBase::set (const unsigned int row, + const unsigned int n_cols, + const unsigned int *col_indices, + const number *values, + const bool elide_zero_values) { // Resize scratch arrays - local_row_length.resize (this->n_block_cols()); - block_col_indices.resize (this->n_block_cols()); - local_col_indices.resize (n_cols); - - // Clear the content in local_row_length + column_indices.resize (this->n_block_cols()); + column_values.resize (this->n_block_cols()); + counter_within_block.resize (this->n_block_cols()); + + // Resize sub-arrays to n_cols. This is a + // bit wasteful, but we resize only a few + // times (then the maximum row length + // won't increase that much any more). for (unsigned int i=0; in_block_cols(); ++i) - local_row_length[i] = 0; + { + column_indices[i].resize(n_cols); + column_values[i].resize(n_cols); + counter_within_block[i] = 0; + } // Go through the column indices to find // out which portions of the values // should be written into which block - // matrix. - { - unsigned int current_block = 0, row_length = 0; - block_col_indices[0] = 0; - for (unsigned int j=0; j - col_index = this->column_block_indices.global_to_local(col_indices[j]); - local_col_indices[j] = col_index.second; - if (col_index.first > current_block) - { - local_row_length[current_block] = row_length; - row_length = 0; - while (col_index.first > current_block) - current_block++; - block_col_indices[current_block] = j; - } + // matrix. We need to restore all the + // data, since we can't be sure that the + // data of one block is stored + // contiguously (in fact, it will be + // intermixed when the data comes from an + // element matrix). + for (unsigned int j=0; jn_block_cols(), - ExcInternalError()); + if (value == 0 && elide_zero_values == true) + continue; + + const std::pair + col_index = this->column_block_indices.global_to_local(col_indices[j]); + + const unsigned int local_index = counter_within_block[col_index.first]++; + + column_indices[col_index.first][local_index] = col_index.second; + column_values[col_index.first][local_index] = value; + } #ifdef DEBUG // If in debug mode, do a check whether // the right length has been obtained. - unsigned int length = 0; - for (unsigned int i=0; in_block_cols(); ++i) - length += local_row_length[i]; - Assert (length == n_cols, - ExcDimensionMismatch(length, n_cols)); + unsigned int length = 0; + for (unsigned int i=0; in_block_cols(); ++i) + length += counter_within_block[i]; + Assert (length <= n_cols, ExcInternalError()); #endif - } // Now we found out about where the // individual columns should start and @@ -1901,15 +1919,15 @@ BlockMatrixBase::set (const unsigned int row, row_index = this->row_block_indices.global_to_local (row); for (unsigned int block_col=0; block_col::add (const unsigned int i, template +template inline void BlockMatrixBase::add (const std::vector &row_indices, const std::vector &col_indices, - const FullMatrix &values, + const FullMatrix &values, const bool elide_zero_values) { Assert (row_indices.size() == values.m(), @@ -1966,10 +1985,11 @@ BlockMatrixBase::add (const std::vector &row_indices, template +template inline void BlockMatrixBase::add (const std::vector &indices, - const FullMatrix &values, + const FullMatrix &values, const bool elide_zero_values) { Assert (indices.size() == values.m(), @@ -1984,11 +2004,12 @@ BlockMatrixBase::add (const std::vector &indices, template +template inline void BlockMatrixBase::add (const unsigned int row, const std::vector &col_indices, - const std::vector &values, + const std::vector &values, const bool elide_zero_values) { Assert (col_indices.size() == values.size(), @@ -2004,65 +2025,67 @@ BlockMatrixBase::add (const unsigned int row, // we need to calculate to each position // the location in the global array. template +template inline void BlockMatrixBase::add (const unsigned int row, const unsigned int n_cols, const unsigned int *col_indices, - const value_type *values, + const number *values, const bool elide_zero_values) { // TODO: Look over this to find out // whether we can do that more // efficiently. - // Resize scratch arrays - block_col_indices.resize (this->n_block_cols()); - local_row_length.resize (this->n_block_cols()); - local_col_indices.resize (n_cols); - - // Clear the content in local_row_length + column_indices.resize (this->n_block_cols()); + column_values.resize (this->n_block_cols()); + counter_within_block.resize (this->n_block_cols()); + + // Resize sub-arrays to n_cols. This is a + // bit wasteful, but we resize only a few + // times (then the maximum row length + // won't increase that much any more). for (unsigned int i=0; in_block_cols(); ++i) - local_row_length[i] = 0; + { + column_indices[i].resize(n_cols); + column_values[i].resize(n_cols); + counter_within_block[i] = 0; + } // Go through the column indices to find // out which portions of the values // should be written into which block - // matrix. - { - unsigned int current_block = 0, row_length = 0; - block_col_indices[0] = 0; - for (unsigned int j=0; j - col_index = this->column_block_indices.global_to_local(col_indices[j]); - local_col_indices[j] = col_index.second; - if (col_index.first > current_block) - { - local_row_length[current_block] = row_length; - row_length = 0; - while (col_index.first > current_block) - current_block++; - block_col_indices[current_block] = j; - } + // matrix. We need to restore all the + // data, since we can't be sure that the + // data of one block is stored + // contiguously (in fact, it will be + // intermixed when the data comes from an + // element matrix). + for (unsigned int j=0; jn_block_cols(), - ExcInternalError()); + if (value == 0 && elide_zero_values == true) + continue; + + const std::pair + col_index = this->column_block_indices.global_to_local(col_indices[j]); + + const unsigned int local_index = counter_within_block[col_index.first]++; + + column_indices[col_index.first][local_index] = col_index.second; + column_values[col_index.first][local_index] = value; + } #ifdef DEBUG // If in debug mode, do a check whether // the right length has been obtained. - unsigned int length = 0; - for (unsigned int i=0; in_block_cols(); ++i) - length += local_row_length[i]; - Assert (length == n_cols, - ExcDimensionMismatch(length, n_cols)); + unsigned int length = 0; + for (unsigned int i=0; in_block_cols(); ++i) + length += counter_within_block[i]; + Assert (length <= n_cols, ExcInternalError()); #endif - } // Now we found out about where the // individual columns should start and @@ -2073,15 +2096,15 @@ BlockMatrixBase::add (const unsigned int row, row_index = this->row_block_indices.global_to_local (row); for (unsigned int block_col=0; block_colfalse, i.e., even zero * values are treated. */ + template void set (const std::vector &indices, - const FullMatrix &full_matrix, + const FullMatrix &full_matrix, const bool elide_zero_values = false); /** @@ -804,9 +805,10 @@ class SparseMatrix : public virtual Subscriptor * different local-to-global indexing * on rows and columns, respectively. */ + template void set (const std::vector &row_indices, const std::vector &col_indices, - const FullMatrix &full_matrix, + const FullMatrix &full_matrix, const bool elide_zero_values = false); /** @@ -827,9 +829,10 @@ class SparseMatrix : public virtual Subscriptor * false, i.e., even zero * values are treated. */ + template void set (const unsigned int row, const std::vector &col_indices, - const std::vector &values, + const std::vector &values, const bool elide_zero_values = false); /** @@ -848,10 +851,11 @@ class SparseMatrix : public virtual Subscriptor * false, i.e., even zero * values are inserted/replaced. */ + template void set (const unsigned int row, const unsigned int n_cols, const unsigned int *col_indices, - const number *values, + const number2 *values, const bool elide_zero_values = false); /** @@ -894,8 +898,9 @@ class SparseMatrix : public virtual Subscriptor * i.e., zero values won't be added * into the matrix. */ + template void add (const std::vector &indices, - const FullMatrix &full_matrix, + const FullMatrix &full_matrix, const bool elide_zero_values = true); /** @@ -905,9 +910,10 @@ class SparseMatrix : public virtual Subscriptor * different local-to-global indexing * on rows and columns, respectively. */ + template void add (const std::vector &row_indices, const std::vector &col_indices, - const FullMatrix &full_matrix, + const FullMatrix &full_matrix, const bool elide_zero_values = true); /** @@ -927,9 +933,10 @@ class SparseMatrix : public virtual Subscriptor * i.e., zero values won't be added * into the matrix. */ + template void add (const unsigned int row, const std::vector &col_indices, - const std::vector &values, + const std::vector &values, const bool elide_zero_values = true); /** @@ -949,10 +956,11 @@ class SparseMatrix : public virtual Subscriptor * i.e., zero values won't be added * into the matrix. */ + template void add (const unsigned int row, const unsigned int n_cols, const unsigned int *col_indices, - const number *values, + const number2 *values, const bool elide_zero_values = true); /** @@ -2033,10 +2041,11 @@ SparseMatrix::set (const unsigned int i, template +template inline void SparseMatrix::set (const std::vector &indices, - const FullMatrix &values, + const FullMatrix &values, const bool elide_zero_values) { Assert (indices.size() == values.m(), @@ -2051,11 +2060,12 @@ SparseMatrix::set (const std::vector &indices, template +template inline void SparseMatrix::set (const std::vector &row_indices, const std::vector &col_indices, - const FullMatrix &values, + const FullMatrix &values, const bool elide_zero_values) { Assert (row_indices.size() == values.m(), @@ -2071,11 +2081,12 @@ SparseMatrix::set (const std::vector &row_indices, template +template inline void SparseMatrix::set (const unsigned int row, const std::vector &col_indices, - const std::vector &values, + const std::vector &values, const bool elide_zero_values) { Assert (col_indices.size() == values.size(), @@ -2087,13 +2098,14 @@ SparseMatrix::set (const unsigned int row, -template +template +template inline void SparseMatrix::set (const unsigned int row, const unsigned int n_cols, const unsigned int *col_indices, - const number *values, + const number2 *values, const bool elide_zero_values) { Assert (cols != 0, ExcNotInitialized()); @@ -2169,10 +2181,11 @@ SparseMatrix::add (const unsigned int i, template +template inline void SparseMatrix::add (const std::vector &indices, - const FullMatrix &values, + const FullMatrix &values, const bool elide_zero_values) { Assert (indices.size() == values.m(), @@ -2187,11 +2200,12 @@ SparseMatrix::add (const std::vector &indices, template +template inline void SparseMatrix::add (const std::vector &row_indices, const std::vector &col_indices, - const FullMatrix &values, + const FullMatrix &values, const bool elide_zero_values) { Assert (row_indices.size() == values.m(), @@ -2207,11 +2221,12 @@ SparseMatrix::add (const std::vector &row_indices, template +template inline void SparseMatrix::add (const unsigned int row, const std::vector &col_indices, - const std::vector &values, + const std::vector &values, const bool elide_zero_values) { Assert (col_indices.size() == values.size(), @@ -2223,13 +2238,14 @@ SparseMatrix::add (const unsigned int row, -template +template +template inline void SparseMatrix::add (const unsigned int row, const unsigned int n_cols, const unsigned int *col_indices, - const number *values, + const number2 *values, const bool elide_zero_values) { Assert (cols != 0, ExcNotInitialized()); diff --git a/deal.II/lac/source/trilinos_sparse_matrix.cc b/deal.II/lac/source/trilinos_sparse_matrix.cc index d8fed0969d..a71a854cf9 100755 --- a/deal.II/lac/source/trilinos_sparse_matrix.cc +++ b/deal.II/lac/source/trilinos_sparse_matrix.cc @@ -283,7 +283,10 @@ namespace TrilinosWrappers // into it, and then build up the matrix // from the graph. This is considerable // faster than directly filling elements - // into the matrix. + // into the matrix. Moreover, it consumes + // less memory, since the internal + // reordering does not need to be done on + // all the double values. // TODO: There seems to be problem in // Epetra when a Graph/matrix is -- 2.39.5