From: bangerth Date: Wed, 27 Aug 2008 12:20:47 +0000 (+0000) Subject: Add an add() function. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=fe66dca30ae85ec24497ccb99a0ef6df835494d2;p=dealii-svn.git Add an add() function. git-svn-id: https://svn.dealii.org/trunk@16673 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/lac/include/lac/trilinos_sparse_matrix.h b/deal.II/lac/include/lac/trilinos_sparse_matrix.h index bce08ab59c..98b7da3c3f 100755 --- a/deal.II/lac/include/lac/trilinos_sparse_matrix.h +++ b/deal.II/lac/include/lac/trilinos_sparse_matrix.h @@ -954,6 +954,15 @@ namespace TrilinosWrappers const Vector &x, const Vector &b) const; + /** + * Add matrix scaled by + * factor to this matrix, + * i.e. the matrix factor*matrix + * is added to this. + */ + void add (const TrilinosScalar factor, + const SparseMatrix &matrix); + /** * STL-like iterator with the * first entry. @@ -1021,14 +1030,14 @@ namespace TrilinosWrappers */ void write_ascii (); - /** - * Print the matrix to the given - * stream, using the format - * (line,col) value, - * i.e. one nonzero entry of the - * matrix per line. - */ - void print (std::ostream &out) const; + /** + * Print the matrix to the given + * stream, using the format + * (line,col) value, + * i.e. one nonzero entry of the + * matrix per line. + */ + void print (std::ostream &out) const; // TODO: Write an overloading // of the operator << for output. diff --git a/deal.II/lac/source/trilinos_sparse_matrix.cc b/deal.II/lac/source/trilinos_sparse_matrix.cc index ef7c0ef741..0dd2f390a8 100755 --- a/deal.II/lac/source/trilinos_sparse_matrix.cc +++ b/deal.II/lac/source/trilinos_sparse_matrix.cc @@ -792,6 +792,54 @@ namespace TrilinosWrappers + void + SparseMatrix::add (const TrilinosScalar factor, + const SparseMatrix &rhs) + { + Assert (rhs.m() == m(), ExcDimensionMismatch (rhs.m(), m())); + Assert (rhs.n() == n(), ExcDimensionMismatch (rhs.n(), n())); + + // I bet that there must be a better way + // to do this but it has not been found: + // currently, we simply go through each + // row of the argument matrix, copy it, + // scale it, and add it to the current + // matrix. that's probably not the most + // efficient way to do things. + const std::pair + local_range = rhs.local_range(); + + unsigned int max_row_length = 0; + for (unsigned int row=local_range.first; + row < local_range.second; ++row) + max_row_length + = std::max (max_row_length, + static_cast(rhs.matrix->NumGlobalEntries(row))); + + std::vector column_indices (max_row_length); + std::vector values (max_row_length); + + for (unsigned int row=local_range.first; + row < local_range.second; ++row) + { + int n_entries; + rhs.matrix->ExtractGlobalRowCopy (row, max_row_length, + n_entries, + &values[0], + &column_indices[0]); + for (int i=0; iSumIntoGlobalValues (row, n_entries, + &values[0], + &column_indices[0]); + } + + compress (); + } + + + // TODO: Currently this only flips // a flag that tells Trilinos that // any application should be done with