]> https://gitweb.dealii.org/ - dealii.git/commitdiff
de-templatize SparsityPattern::copy_from. 4590/head
authorDavid Wells <wellsd2@rpi.edu>
Fri, 7 Jul 2017 23:52:28 +0000 (19:52 -0400)
committerDavid Wells <wellsd2@rpi.edu>
Sat, 8 Jul 2017 16:21:27 +0000 (12:21 -0400)
This template is instantiated for only SparsityPattern and is overloaded for
DynamicSparsityPattern, so it is a bit simpler to just overload a non-template
method.

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

index 369a58c554cde8d47d7e18a9898f401b1f9f98a3..f6f6f53b04e1a76434d7da226efdf38f0c8760e5 100644 (file)
@@ -38,6 +38,7 @@
 DEAL_II_NAMESPACE_OPEN
 
 class SparsityPattern;
+class DynamicSparsityPattern;
 class ChunkSparsityPattern;
 template <typename number> class FullMatrix;
 template <typename number> class SparseMatrix;
@@ -630,15 +631,16 @@ public:
                   const ForwardIterator end);
 
   /**
-   * Copy data from an object of type DynamicSparsityPattern. Although not a
-   * compressed sparsity pattern, this function is also instantiated if the
-   * argument is of type SparsityPattern (i.e., the current class). Previous
-   * content of this object is lost, and the sparsity pattern is in compressed
-   * mode afterwards.
+   * Copy data from a DynamicSparsityPattern. Previous content of this object
+   * is lost, and the sparsity pattern is in compressed mode afterwards.
    */
-  template <typename SparsityPatternType>
-  void copy_from (const SparsityPatternType &dsp);
+  void copy_from (const DynamicSparsityPattern &dsp);
 
+  /**
+   * Copy data from a SparsityPattern. Previous content of this object is
+   * lost, and the sparsity pattern is in compressed mode afterwards.
+   */
+  void copy_from (const SparsityPattern &sp);
 
   /**
    * Take a full matrix and use its nonzero entries to generate a sparse
index bc35387af8bb93b4af0e8168b2b68c8b89144625..3a51a46b9666f29a3cda20186054f7a8ccec22de 100644 (file)
@@ -482,33 +482,32 @@ SparsityPattern::compress ()
 
 
 
-template <typename SparsityPatternType>
 void
-SparsityPattern::copy_from (const SparsityPatternType &dsp)
+SparsityPattern::copy_from (const SparsityPattern &sp)
 {
   // first determine row lengths for each row. if the matrix is quadratic,
   // then we might have to add an additional entry for the diagonal, if that
   // is not yet present. as we have to call compress anyway later on, don't
   // bother to check whether that diagonal entry is in a certain row or not
-  const bool do_diag_optimize = (dsp.n_rows() == dsp.n_cols());
-  std::vector<unsigned int> row_lengths (dsp.n_rows());
-  for (size_type i=0; i<dsp.n_rows(); ++i)
+  const bool do_diag_optimize = (sp.n_rows() == sp.n_cols());
+  std::vector<unsigned int> row_lengths (sp.n_rows());
+  for (size_type i=0; i<sp.n_rows(); ++i)
     {
-      row_lengths[i] = dsp.row_length(i);
-      if (do_diag_optimize && !dsp.exists(i,i))
+      row_lengths[i] = sp.row_length(i);
+      if (do_diag_optimize && !sp.exists(i,i))
         ++row_lengths[i];
     }
-  reinit (dsp.n_rows(), dsp.n_cols(), row_lengths);
+  reinit (sp.n_rows(), sp.n_cols(), row_lengths);
 
   // now enter all the elements into the matrix, if there are any. note that
   // if the matrix is quadratic, then we already have the diagonal element
   // preallocated
   if (n_rows() != 0 && n_cols() != 0)
-    for (size_type row = 0; row<dsp.n_rows(); ++row)
+    for (size_type row = 0; row<sp.n_rows(); ++row)
       {
         size_type *cols = &colnums[rowstart[row]] + (do_diag_optimize ? 1 : 0);
-        typename SparsityPatternType::iterator col_num = dsp.begin (row),
-                                               end_row = dsp.end (row);
+        typename SparsityPattern::iterator col_num = sp.begin (row),
+                                           end_row = sp.end (row);
 
         for (; col_num != end_row; ++col_num)
           {
@@ -519,8 +518,8 @@ SparsityPattern::copy_from (const SparsityPatternType &dsp)
       }
 
   // do not need to compress the sparsity pattern since we already have
-  // allocated the right amount of data, and the SparsityPatternType data is sorted,
-  // too.
+  // allocated the right amount of data, and the SparsityPattern data is
+  // sorted, too.
   compressed = true;
 }
 
@@ -531,7 +530,6 @@ SparsityPattern::copy_from (const SparsityPatternType &dsp)
 // entries. DynamicSparsityPattern::iterator can show quadratic complexity in
 // case many rows are empty and the begin() method needs to jump to the next
 // free row. Otherwise, the code is exactly the same as above.
-template <>
 void
 SparsityPattern::copy_from (const DynamicSparsityPattern &dsp)
 {
@@ -1038,7 +1036,6 @@ SparsityPattern::memory_consumption () const
 
 
 // explicit instantiations
-template void SparsityPattern::copy_from<SparsityPattern> (const SparsityPattern &);
 template void SparsityPattern::copy_from<float> (const FullMatrix<float> &);
 template void SparsityPattern::copy_from<double> (const FullMatrix<double> &);
 

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.