From: Wolfgang Bangerth Date: Tue, 3 May 2011 22:35:11 +0000 (+0000) Subject: Patch by Martin Genet: Make SparseDirectUMFPACK work with SparseMatrixEZ as well. X-Git-Tag: v8.0.0~4086 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=cf62734abd94dc187307de4be7e6fbfd19cffc98;p=dealii.git Patch by Martin Genet: Make SparseDirectUMFPACK work with SparseMatrixEZ as well. git-svn-id: https://svn.dealii.org/trunk@23677 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index d45761a19b..315aa47ef0 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -30,6 +30,11 @@ inconvenience this causes.
    +
  1. New: The SparseDirectUMFPACK class can now also deal with matrices +provided in SparseMatrixEZ format. +
    +(Martin Genet, 2011/05/04) +
  2. New: The new tutorial program step-46 shows how to couple different models defined on subsets of the domain, in this case Stokes flow around an elastic solid. The diff --git a/deal.II/include/deal.II/lac/block_sparse_matrix.h b/deal.II/include/deal.II/lac/block_sparse_matrix.h index a51df16d76..699557df9a 100644 --- a/deal.II/include/deal.II/lac/block_sparse_matrix.h +++ b/deal.II/include/deal.II/lac/block_sparse_matrix.h @@ -192,6 +192,12 @@ class BlockSparseMatrix : public BlockMatrixBase > */ bool empty () const; + /** + * Return the number of entries + * in a specific row. + */ + unsigned int get_row_length (const unsigned int row) const; + /** * Return the number of nonzero * elements of this diff --git a/deal.II/include/deal.II/lac/block_sparse_matrix.templates.h b/deal.II/include/deal.II/lac/block_sparse_matrix.templates.h index 4da7ebb4d0..6453d3215c 100644 --- a/deal.II/include/deal.II/lac/block_sparse_matrix.templates.h +++ b/deal.II/include/deal.II/lac/block_sparse_matrix.templates.h @@ -127,6 +127,15 @@ BlockSparseMatrix::empty () const +template +unsigned int +BlockSparseMatrix::get_row_length (const unsigned int row) const +{ + return sparsity_pattern->row_length(row); +} + + + template unsigned int BlockSparseMatrix::n_nonzero_elements () const diff --git a/deal.II/include/deal.II/lac/sparse_direct.h b/deal.II/include/deal.II/lac/sparse_direct.h index 22f5060d62..34bc1a9a14 100644 --- a/deal.II/include/deal.II/lac/sparse_direct.h +++ b/deal.II/include/deal.II/lac/sparse_direct.h @@ -21,6 +21,7 @@ #include #include #include +#include #include #ifdef DEAL_II_USE_MUMPS @@ -1241,12 +1242,16 @@ class SparseDirectUMFPACK : public Subscriptor * Make sure that the arrays Ai * and Ap are sorted in each * row. UMFPACK wants it this - * way. We need to have two + * way. We need to have three * versions of this function, one - * for the usual SparseMatrix and + * for the usual SparseMatrix, one + * for the SparseMatrixEZ, and * one for the BlockSparseMatrix * classes */ + template + void sort_arrays (const SparseMatrixEZ &); + template void sort_arrays (const SparseMatrix &); diff --git a/deal.II/include/deal.II/lac/sparse_matrix.h b/deal.II/include/deal.II/lac/sparse_matrix.h index e68b183c6f..1912965049 100644 --- a/deal.II/include/deal.II/lac/sparse_matrix.h +++ b/deal.II/include/deal.II/lac/sparse_matrix.h @@ -749,15 +749,21 @@ class SparseMatrix : public virtual Subscriptor */ unsigned int n () const; - /** - * Return the number of nonzero - * elements of this - * matrix. Actually, it returns - * the number of entries in the - * sparsity pattern; if any of - * the entries should happen to - * be zero, it is counted anyway. - */ + /** + * Return the number of entries + * in a specific row. + */ + unsigned int get_row_length (const unsigned int row) const; + + /** + * Return the number of nonzero + * elements of this + * matrix. Actually, it returns + * the number of entries in the + * sparsity pattern; if any of + * the entries should happen to + * be zero, it is counted anyway. + */ unsigned int n_nonzero_elements () const; /** diff --git a/deal.II/include/deal.II/lac/sparse_matrix.templates.h b/deal.II/include/deal.II/lac/sparse_matrix.templates.h index c4d1711b57..65d1f2d707 100644 --- a/deal.II/include/deal.II/lac/sparse_matrix.templates.h +++ b/deal.II/include/deal.II/lac/sparse_matrix.templates.h @@ -214,6 +214,16 @@ SparseMatrix::empty () const +template +unsigned int +SparseMatrix::get_row_length (const unsigned int row) const +{ + Assert (cols != 0, ExcNotInitialized()); + return cols->row_length(row); +} + + + template unsigned int SparseMatrix::n_nonzero_elements () const diff --git a/deal.II/include/deal.II/lac/sparse_matrix_ez.h b/deal.II/include/deal.II/lac/sparse_matrix_ez.h index aec09a947a..25232fa98f 100644 --- a/deal.II/include/deal.II/lac/sparse_matrix_ez.h +++ b/deal.II/include/deal.II/lac/sparse_matrix_ez.h @@ -451,6 +451,18 @@ class SparseMatrixEZ : public Subscriptor */ unsigned int n () const; + /** + * Return the number of entries + * in a specific row. + */ + unsigned int get_row_length (const unsigned int row) const; + + /** + * Return the number of nonzero + * elements of this matrix. + */ + unsigned int n_nonzero_elements () const; + /** * Set the element (i,j) to * @p value. Allocates the entry, diff --git a/deal.II/include/deal.II/lac/sparse_matrix_ez.templates.h b/deal.II/include/deal.II/lac/sparse_matrix_ez.templates.h index 7f03da83fd..698e195d37 100644 --- a/deal.II/include/deal.II/lac/sparse_matrix_ez.templates.h +++ b/deal.II/include/deal.II/lac/sparse_matrix_ez.templates.h @@ -385,6 +385,32 @@ SparseMatrixEZ::memory_consumption() const } + +template +unsigned int +SparseMatrixEZ::get_row_length (const unsigned int row) const +{ + return row_info[row].length; +} + + + +template +unsigned int +SparseMatrixEZ::n_nonzero_elements() const +{ + typename std::vector::const_iterator row = row_info.begin(); + const typename std::vector::const_iterator endrow = row_info.end(); + + // Add up entries actually used + unsigned int used = 0; + for (; row != endrow ; ++ row) + used += row->length; + return used; +} + + + template void SparseMatrixEZ::compute_statistics( diff --git a/deal.II/source/lac/sparse_direct.cc b/deal.II/source/lac/sparse_direct.cc index 30f9448ae6..55e8b13e93 100644 --- a/deal.II/source/lac/sparse_direct.cc +++ b/deal.II/source/lac/sparse_direct.cc @@ -1600,6 +1600,27 @@ sort_arrays (const SparseMatrix &matrix) +template +void +SparseDirectUMFPACK:: +sort_arrays (const SparseMatrixEZ &matrix) +{ + //same thing for SparseMatrixEZ + for (unsigned int row=0; row Ai[cursor+1])) + { + std::swap (Ai[cursor], Ai[cursor+1]); + std::swap (Ax[cursor], Ax[cursor+1]); + ++cursor; + } + } +} + + + template void SparseDirectUMFPACK:: @@ -1691,13 +1712,13 @@ factorize (const Matrix &matrix) // anyway. people are supposed to provide // accurate sparsity patterns. Ap.resize (N+1); - Ai.resize (matrix.get_sparsity_pattern().n_nonzero_elements()); - Ax.resize (matrix.get_sparsity_pattern().n_nonzero_elements()); + Ai.resize (matrix.n_nonzero_elements()); + Ax.resize (matrix.n_nonzero_elements()); // first fill row lengths array Ap[0] = 0; for (unsigned int row=1; row<=N; ++row) - Ap[row] = Ap[row-1] + matrix.get_sparsity_pattern().row_length(row-1); + Ap[row] = Ap[row-1] + matrix.get_row_length(row-1); Assert (static_cast(Ap.back()) == Ai.size(), ExcInternalError()); @@ -2028,6 +2049,8 @@ void SparseDirectMA27::solve (const SparseMatrix &matrix, InstantiateUMFPACK(SparseMatrix); InstantiateUMFPACK(SparseMatrix); +InstantiateUMFPACK(SparseMatrixEZ); +InstantiateUMFPACK(SparseMatrixEZ); InstantiateUMFPACK(BlockSparseMatrix); InstantiateUMFPACK(BlockSparseMatrix);