]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Have a number of reinit functions, and merge code with the constructors.
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 15 Mar 2004 16:16:58 +0000 (16:16 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 15 Mar 2004 16:16:58 +0000 (16:16 +0000)
git-svn-id: https://svn.dealii.org/trunk@8741 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/petsc_sparse_matrix.h
deal.II/lac/source/petsc_sparse_matrix.cc

index 0c4a797eebb8de201b5bc28a5a9ac32a071bf992..bdfeb53421c4e5db6062aec8859fb3ff372ef152 100644 (file)
@@ -40,6 +40,12 @@ namespace PETScWrappers
   class SparseMatrix : public MatrixBase
   {
     public:
+                                       /**
+                                        * Default constructor. Create an empty
+                                        * matrix.
+                                        */
+      SparseMatrix ();
+      
                                        /**
                                         * Create a sparse matrix of dimensions
                                         * <tt>m</tt> times <tt>n</tt>, with an
@@ -98,6 +104,62 @@ namespace PETScWrappers
                     const unsigned int               n,
                     const std::vector<unsigned int> &row_lengths,
                     const bool                       is_symmetric = false);
+
+                                       /**
+                                        * Set all matrix entries to zero, but
+                                        * retain the sparsity pattern. This
+                                        * function simply calls the respective
+                                        * function of the base class.
+                                        */
+      void reinit ();
+
+                                       /**
+                                        * Throw away the present matrix and
+                                        * generate one that has the same
+                                        * properties as if it were created by
+                                        * the constructor of this class with
+                                        * the same argument list as the
+                                        * present function.
+                                        */
+      void reinit (const unsigned int m,
+                   const unsigned int n,
+                   const unsigned int n_nonzero_per_row,
+                   const bool         is_symmetric = false);
+
+                                       /**
+                                        * Throw away the present matrix and
+                                        * generate one that has the same
+                                        * properties as if it were created by
+                                        * the constructor of this class with
+                                        * the same argument list as the
+                                        * present function.
+                                        */
+      void reinit (const unsigned int               m,
+                   const unsigned int               n,
+                   const std::vector<unsigned int> &row_lengths,
+                   const bool                       is_symmetric = false);
+
+    private:
+
+                                       /**
+                                        * Do the actual work for the
+                                        * respective reinit() function and the
+                                        * matching constructor, i.e. create a
+                                        * matrix. Getting rid of the previous
+                                        * matrix is left to the caller.
+                                        */
+      void do_reinit (const unsigned int m,
+                      const unsigned int n,
+                      const unsigned int n_nonzero_per_row,
+                      const bool         is_symmetric = false);
+
+                                       /**
+                                        * Same as previous function.
+                                        */
+      void do_reinit (const unsigned int               m,
+                      const unsigned int               n,
+                      const std::vector<unsigned int> &row_lengths,
+                      const bool                       is_symmetric = false);
   };
 }
 
index b690583d2579a3fbfc8fcc908facc0a137198994..5f13542903e72c7e8f10d355ad622c9b8ee88bd7 100644 (file)
 
 namespace PETScWrappers
 {
+  SparseMatrix::SparseMatrix ()
+  {
+    const int m=0, n=0, n_nonzero_per_row=0;
+    const int ierr
+      = MatCreateSeqAIJ(PETSC_COMM_SELF, m, n, n_nonzero_per_row,
+                        0, &matrix);
+    AssertThrow (ierr == 0, ExcPETScError(ierr));
+  }
+
+
+
   SparseMatrix::SparseMatrix (const unsigned int m,
                               const unsigned int n,
                               const unsigned int n_nonzero_per_row,
                               const bool         is_symmetric)
+  {
+    do_reinit (m, n, n_nonzero_per_row, is_symmetric);
+  }
+
+
+
+  SparseMatrix::SparseMatrix (const unsigned int m,
+                              const unsigned int n,
+                              const std::vector<unsigned int> &row_lengths,
+                              const bool         is_symmetric)
+  {
+    do_reinit (m, n, row_lengths, is_symmetric);
+  }
+
+
+
+  void
+  SparseMatrix::reinit (const unsigned int m,
+                        const unsigned int n,
+                        const unsigned int n_nonzero_per_row,
+                        const bool         is_symmetric)
+  {
+                                     // get rid of old matrix and generate a
+                                     // new one
+    const int ierr = MatDestroy (matrix);
+    AssertThrow (ierr == 0, ExcPETScError(ierr));    
+    
+    do_reinit (m, n, n_nonzero_per_row, is_symmetric);
+  }
+
+
+
+  void
+  SparseMatrix::reinit (const unsigned int m,
+                        const unsigned int n,
+                        const std::vector<unsigned int> &row_lengths,
+                        const bool         is_symmetric)
+  {
+                                     // get rid of old matrix and generate a
+                                     // new one
+    const int ierr = MatDestroy (matrix);
+    AssertThrow (ierr == 0, ExcPETScError(ierr));    
+
+    do_reinit (m, n, row_lengths, is_symmetric);
+  }  
+
+
+
+  void
+  SparseMatrix::do_reinit (const unsigned int m,
+                           const unsigned int n,
+                           const unsigned int n_nonzero_per_row,
+                           const bool         is_symmetric)
   {
                                      // use the call sequence indicating only
                                      // a maximal number of elements per row
@@ -44,10 +108,11 @@ namespace PETScWrappers
 
 
 
-  SparseMatrix::SparseMatrix (const unsigned int m,
-                              const unsigned int n,
-                              const std::vector<unsigned int> &row_lengths,
-                              const bool         is_symmetric)
+  void
+  SparseMatrix::do_reinit (const unsigned int m,
+                           const unsigned int n,
+                           const std::vector<unsigned int> &row_lengths,
+                           const bool         is_symmetric)
   {
     Assert (row_lengths.size() == m,
             ExcDimensionMismatch (row_lengths.size(), m));

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.