]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Patch by Martin Genet: Make SparseDirectUMFPACK work with SparseMatrixEZ as well.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Tue, 3 May 2011 22:35:11 +0000 (22:35 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Tue, 3 May 2011 22:35:11 +0000 (22:35 +0000)
git-svn-id: https://svn.dealii.org/trunk@23677 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/include/deal.II/lac/block_sparse_matrix.h
deal.II/include/deal.II/lac/block_sparse_matrix.templates.h
deal.II/include/deal.II/lac/sparse_direct.h
deal.II/include/deal.II/lac/sparse_matrix.h
deal.II/include/deal.II/lac/sparse_matrix.templates.h
deal.II/include/deal.II/lac/sparse_matrix_ez.h
deal.II/include/deal.II/lac/sparse_matrix_ez.templates.h
deal.II/source/lac/sparse_direct.cc

index d45761a19b379b9771611d0dffe9edf854e3f18d..315aa47ef0390065802f5876df5b8a9695056ae6 100644 (file)
@@ -30,6 +30,11 @@ inconvenience this causes.
 
 <ol>
 
+<li> New: The SparseDirectUMFPACK class can now also deal with matrices
+provided in SparseMatrixEZ format.
+<br>
+(Martin Genet, 2011/05/04)
+
 <li> 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
index a51df16d76f1a3ace27784299d54430c8f17c6c3..699557df9a89dc47e6bc98eb6a529e46124ca51d 100644 (file)
@@ -192,6 +192,12 @@ class BlockSparseMatrix : public BlockMatrixBase<SparseMatrix<number> >
                                      */
     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
index 4da7ebb4d0c870073c2a882e7d00a9cad01d9a14..6453d3215cfbef17b9e304c7922981a93f99a64c 100644 (file)
@@ -127,6 +127,15 @@ BlockSparseMatrix<number>::empty () const
 
 
 
+template <typename number>
+unsigned int
+BlockSparseMatrix<number>::get_row_length (const unsigned int row) const
+{
+  return sparsity_pattern->row_length(row);
+}
+
+
+
 template <typename number>
 unsigned int
 BlockSparseMatrix<number>::n_nonzero_elements () const
index 22f5060d62be4e96aa997412247f93529540fbcd..34bc1a9a149ab604deed466c854ff216107f0f8f 100644 (file)
@@ -21,6 +21,7 @@
 #include <base/thread_management.h>
 #include <lac/vector.h>
 #include <lac/sparse_matrix.h>
+#include <lac/sparse_matrix_ez.h>
 #include <lac/block_sparse_matrix.h>
 
 #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 <typename number>
+    void sort_arrays (const SparseMatrixEZ<number> &);
+
     template <typename number>
     void sort_arrays (const SparseMatrix<number> &);
     
index e68b183c6fea16ede753b38c60dc747822796f81..19129650492fec835812d0c34c330b0904eddb88 100644 (file)
@@ -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;
 
                                     /**
index c4d1711b57f4b73ffbc3d5b83b8c8e9e42ecdef0..65d1f2d7072847cf618369f45d0d9989013d410c 100644 (file)
@@ -214,6 +214,16 @@ SparseMatrix<number>::empty () const
 
 
 
+template <typename number>
+unsigned int
+SparseMatrix<number>::get_row_length (const unsigned int row) const
+{
+  Assert (cols != 0, ExcNotInitialized());
+  return cols->row_length(row);
+}
+
+
+
 template <typename number>
 unsigned int
 SparseMatrix<number>::n_nonzero_elements () const
index aec09a947a1caf09f0e592caad4c7746836d8b06..25232fa98f5c0b847c44c49fcfe1198359c3d9dd 100644 (file)
@@ -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 <tt>(i,j)</tt> to
                                      * @p value. Allocates the entry,
index 7f03da83fd94c9f6f6f55dc4090039c1f828a425..698e195d37eaf250140c4d1617a25d9861081d93 100644 (file)
@@ -385,6 +385,32 @@ SparseMatrixEZ<number>::memory_consumption() const
 }
 
 
+
+template <typename number>
+unsigned int
+SparseMatrixEZ<number>::get_row_length (const unsigned int row) const
+{
+    return row_info[row].length;
+}
+
+
+
+template <typename number>
+unsigned int
+SparseMatrixEZ<number>::n_nonzero_elements() const
+{
+  typename std::vector<RowInfo>::const_iterator row = row_info.begin();
+  const typename std::vector<RowInfo>::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 <typename number>
 void
 SparseMatrixEZ<number>::compute_statistics(
index 30f9448ae696dcee498e5d2fdb66ccf619595f81..55e8b13e934effae3661724b4e398632f7543f6f 100644 (file)
@@ -1600,6 +1600,27 @@ sort_arrays (const SparseMatrix<number> &matrix)
 
 
 
+template <typename number>
+void
+SparseDirectUMFPACK::
+sort_arrays (const SparseMatrixEZ<number> &matrix)
+{
+                                   //same thing for SparseMatrixEZ
+  for (unsigned int row=0; row<matrix.m(); ++row)
+    {
+      long int cursor = Ap[row];
+      while ((cursor < Ap[row+1]-1) &&
+             (Ai[cursor] > Ai[cursor+1]))
+        {
+          std::swap (Ai[cursor], Ai[cursor+1]);
+          std::swap (Ax[cursor], Ax[cursor+1]);
+          ++cursor;
+        }
+    }
+}
+
+
+
 template <typename number>
 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<unsigned int>(Ap.back()) == Ai.size(),
           ExcInternalError());
   
@@ -2028,6 +2049,8 @@ void SparseDirectMA27::solve (const SparseMatrix<float>  &matrix,
 
 InstantiateUMFPACK(SparseMatrix<double>);
 InstantiateUMFPACK(SparseMatrix<float>);
+InstantiateUMFPACK(SparseMatrixEZ<double>);
+InstantiateUMFPACK(SparseMatrixEZ<float>);
 InstantiateUMFPACK(BlockSparseMatrix<double>);
 InstantiateUMFPACK(BlockSparseMatrix<float>);
 

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.