Use a smarter algorithm for inserting a list of sorted elements into a CompressedSimp...

git-svn-id: https://svn.dealii.org/trunk@18547 0785d39b-7218-0410-832d-ea1e28bc413d

diff --git a/deal.II/lac/include/lac/block_vector.h b/deal.II/lac/include/lac/block_vector.h
index 77bd88e14d3d8af49b9c2a05b76b4e734694eb32..29eb14f25f2c093b40cb3d72fa179a473905c1d7 100644 (file)
--- a/deal.II/lac/include/lac/block_vector.h
+++ b/deal.II/lac/include/lac/block_vector.h
@@ -19,7 +19,6 @@
 #include <lac/vector.h>
 #include <lac/block_indices.h>
 #include <lac/block_vector_base.h>
-#include <lac/trilinos_block_vector.h>
 
 #include <cstdio>
 #include <vector>
@@ -503,16 +502,6 @@ BlockVector<Number>::operator = (const BlockVector<Number2> &v)
 }
 
 
-#ifdef DEAL_II_USE_TRILINOS
-template <typename Number>
-inline
-BlockVector<Number> &
-BlockVector<Number>::operator = (const TrilinosWrappers::BlockVector &v)
-{
-  BaseClass::operator = (v);
-  return *this;
-}
-#endif
 
 template <typename Number>
 void BlockVector<Number>::scale (const value_type factor)

diff --git a/deal.II/lac/include/lac/block_vector.templates.h b/deal.II/lac/include/lac/block_vector.templates.h
index 0a6e1e8b78317d779b78067c0fa312f3044c356f..c5c7ae89c68a87a8ed2cef5ea0fab8cb488053b0 100644 (file)
--- a/deal.II/lac/include/lac/block_vector.templates.h
+++ b/deal.II/lac/include/lac/block_vector.templates.h
@@ -16,6 +16,7 @@
 
 #include <base/config.h>
 #include <lac/block_vector.h>
+#include <lac/trilinos_block_vector.h>
 #include <cmath>
 #include <algorithm>
 
@@ -122,6 +123,17 @@ BlockVector<Number>::~BlockVector ()
 {}
 
 
+#ifdef DEAL_II_USE_TRILINOS
+template <typename Number>
+inline
+BlockVector<Number> &
+BlockVector<Number>::operator = (const TrilinosWrappers::BlockVector &v)
+{
+  BaseClass::operator = (v);
+  return *this;
+}
+#endif
+
 
 template <typename Number>
 void BlockVector<Number>::swap (BlockVector<Number> &v)

diff --git a/deal.II/lac/include/lac/compressed_simple_sparsity_pattern.h b/deal.II/lac/include/lac/compressed_simple_sparsity_pattern.h
index fb8f716127ed5a106d5cc98f2dce46bf5e889750..0d25f1dddc71c7a1bd4b77c7907aa074cf5687a1 100644 (file)
--- a/deal.II/lac/include/lac/compressed_simple_sparsity_pattern.h
+++ b/deal.II/lac/include/lac/compressed_simple_sparsity_pattern.h
@@ -20,6 +20,7 @@
 
 #include <vector>
 #include <algorithm>
+#include <iostream>
 
 DEAL_II_NAMESPACE_OPEN
 
@@ -366,7 +367,8 @@ class CompressedSimpleSparsityPattern : public Subscriptor
                                           */
         template <typename ForwardIterator>
        void add_entries (ForwardIterator begin,
-                         ForwardIterator end);
+                         ForwardIterator end,
+                         const bool indices_are_sorted);
     };
 
 
@@ -417,18 +419,74 @@ template <typename ForwardIterator>
 inline
 void
 CompressedSimpleSparsityPattern::Line::add_entries (ForwardIterator begin,
-                                                   ForwardIterator end)
+                                                   ForwardIterator end,
+                                                   const bool      indices_are_sorted)
 {
-  if (end - begin <= 0)
+  const int n_elements = end - begin;
+  if (n_elements <= 0)
     return;
 
+  const unsigned int n_cols = static_cast<unsigned int>(n_elements);
+  const unsigned int stop_size = entries.size() + n_cols;
+
+  if (indices_are_sorted == true && n_elements > 3)
+    {
+      if (entries.size() == 0 || entries.back() < *begin)
+       {
+         entries.insert(entries.end(), begin, end);
+         return;
+       }
+
+                                  // resize vector by just inserting the
+                                  // list
+      const unsigned int col = *begin;
+      std::vector<unsigned int>::iterator it = 
+       std::lower_bound(entries.begin(), entries.end(), col);
+      const unsigned int pos1 = it - entries.begin();
+      entries.insert (it, begin, end);
+      it = entries.begin() + pos1;
+
+                                  // now merge the two lists.
+      ForwardIterator my_it = begin;
+      std::vector<unsigned int>::iterator it2 = it + n_cols;
+
+                                  // as long as there are indices both in
+                                  // the end of the entries list and in the
+                                  // input list
+      while (my_it != end && it2 != entries.end())
+       {
+         if (*my_it < *it2)
+           *it++ = *my_it++;
+         else if (*my_it == *it2)
+           {
+             *it++ = *it2++;
+             ++my_it;
+           }
+         else
+           *it++ = *it2++;
+       }
+                                  // in case there are indices left in the
+                                  // input list
+      while (my_it != end)
+       *it++ = *my_it++;
+
+                                  // in case there are indices left in the
+                                  // end of entries
+      while (it2 != entries.end())
+       *it++ = *it2++;
+
+                                  // resize
+      const unsigned int new_size = it - entries.begin();
+      Assert (new_size <= stop_size, ExcInternalError());
+      entries.resize (new_size);
+      return;
+    }
+
   ForwardIterator my_it = begin;
-  const unsigned int n_cols = static_cast<unsigned int>(end - begin);
 
   // If necessary, increase the size of the array. In order to avoid
   // allocating just a few entries every time, use five elements at a time
   // at least.
-  const unsigned int stop_size = entries.size() + (n_cols > 5 ? n_cols : 5);
   if (stop_size > entries.capacity())
     entries.reserve (stop_size);
 
@@ -521,11 +579,11 @@ void
 CompressedSimpleSparsityPattern::add_entries (const unsigned int row,
                                              ForwardIterator begin,
                                              ForwardIterator end,
-                                             const bool         /*indices_are_sorted*/)
+                                             const bool      indices_are_sorted)
 {
   Assert (row < rows, ExcIndexRange (row, 0, rows));
 
-  lines[row].add_entries (begin, end);
+  lines[row].add_entries (begin, end, indices_are_sorted);
 }
 
 

diff --git a/deal.II/lac/include/lac/constraint_matrix.templates.h b/deal.II/lac/include/lac/constraint_matrix.templates.h
index 56de0b588bf9d0397901b199f860de228ada1348..307417bac8627f19ad8adae4908d7426988164dd 100644 (file)
--- a/deal.II/lac/include/lac/constraint_matrix.templates.h
+++ b/deal.II/lac/include/lac/constraint_matrix.templates.h
@@ -2006,7 +2006,8 @@ add_entries_local_to_global (const std::vector<unsigned int> &local_dof_indices,
                                   // into the vector
       Threads::ThreadMutex::ScopedLock lock(mutex);
       if (col_counter > 0)
-       sparsity_pattern.add_entries(row, cols.begin(), cols.begin()+col_counter);
+       sparsity_pattern.add_entries(row, cols.begin(), cols.begin()+col_counter,
+                                    true);
     }
 }
 
@@ -2337,7 +2338,8 @@ add_entries_local_to_global (const std::vector<unsigned int> &local_dof_indices,
            if (current_length > 0)
              sparsity_pattern.block(row_block.first, i).add_entries(row_block.second, 
                                                                     cols.begin()+col_begins, 
-                                                                    cols.begin()+col_begins+current_length);
+                                                                    cols.begin()+col_begins+current_length,
+                                                                    true);
            col_begins = actual_block_ends[i];
          }
       }

diff --git a/deal.II/lac/include/lac/trilinos_block_sparse_matrix.h b/deal.II/lac/include/lac/trilinos_block_sparse_matrix.h
index 973a204fea7574864c89c05da25cece975766304..c03cf2ab77f79031a228a405a8670efb1d05230d 100644 (file)
--- a/deal.II/lac/include/lac/trilinos_block_sparse_matrix.h
+++ b/deal.II/lac/include/lac/trilinos_block_sparse_matrix.h
@@ -17,7 +17,6 @@
 #include <base/config.h>
 #include <base/table.h>
 #include <lac/block_matrix_base.h>
-#include <lac/block_sparse_matrix.h>
 #include <lac/trilinos_sparse_matrix.h>
 #include <lac/trilinos_block_vector.h>
 #include <lac/full_matrix.h>
@@ -37,6 +36,7 @@ class BlockSparsityPattern;
 class BlockCompressedSparsityPattern;
 class BlockCompressedSetSparsityPattern;
 class BlockCompressedSimpleSparsityPattern;
+template <typename number> class BlockSparseMatrix;
 
 
 namespace TrilinosWrappers

diff --git a/deal.II/lac/include/lac/trilinos_precondition.h b/deal.II/lac/include/lac/trilinos_precondition.h
index a9f2cc93cb0073436b6392b1b17b5249748c6b4f..791f59ce6664eeb5d9dcf00406ffb1bf01644660 100755 (executable)
--- a/deal.II/lac/include/lac/trilinos_precondition.h
+++ b/deal.II/lac/include/lac/trilinos_precondition.h
@@ -16,8 +16,6 @@
 
 #include <base/config.h>
 #include <base/subscriptor.h>
-#include <lac/vector.h>
-#include <lac/sparse_matrix.h>
 
 #include <base/std_cxx0x/shared_ptr.h>
 
@@ -45,6 +43,10 @@ namespace ML_Epetra
 
 DEAL_II_NAMESPACE_OPEN
 
+// forward declarations
+template <typename number> class SparseMatrix;
+template <typename number> class Vector;
+
 /*! @addtogroup TrilinosWrappers
  *@{
  */

diff --git a/deal.II/lac/include/lac/trilinos_precondition_block.h b/deal.II/lac/include/lac/trilinos_precondition_block.h
index 7615f9c5942aa490899380b10b7e1e2b4cceb51b..69d40311c129167a626ddf78e2abb7853d0b306a 100644 (file)
--- a/deal.II/lac/include/lac/trilinos_precondition_block.h
+++ b/deal.II/lac/include/lac/trilinos_precondition_block.h
@@ -16,9 +16,8 @@
 
 #include <base/config.h>
 #include <lac/exceptions.h>
-#include <lac/trilinos_sparse_matrix.h>
-#include <lac/trilinos_block_sparse_matrix.h>
 #include <lac/trilinos_precondition.h>
+#include <lac/trilinos_sparse_matrix.h>
 
 
 #ifdef DEAL_II_USE_TRILINOS
@@ -39,6 +38,8 @@ DEAL_II_NAMESPACE_OPEN
 
 namespace TrilinosWrappers
 {
+                                  // forward declarations
+  class BlockSparseMatrix;
 
 /**
  * This implements an interface class of preconditioners for block

diff --git a/deal.II/lac/include/lac/trilinos_sparse_matrix.h b/deal.II/lac/include/lac/trilinos_sparse_matrix.h
index 452f6f6aedb4767c9d8d65f7d242ba874e977942..3873c70b4d4671d750b1880ad4f80da888acd04c 100755 (executable)
--- a/deal.II/lac/include/lac/trilinos_sparse_matrix.h
+++ b/deal.II/lac/include/lac/trilinos_sparse_matrix.h
@@ -16,7 +16,6 @@
 
 #include <base/config.h>
 #include <base/subscriptor.h>
-#include <lac/sparse_matrix.h>
 #include <lac/full_matrix.h>
 #include <lac/exceptions.h>
 #include <lac/trilinos_vector_base.h>
@@ -41,6 +40,8 @@
 
 DEAL_II_NAMESPACE_OPEN
 
+// forward declaration
+template <typename number> class SparseMatrix;
 
 namespace TrilinosWrappers
 {

diff --git a/deal.II/lac/source/trilinos_block_sparse_matrix.cc b/deal.II/lac/source/trilinos_block_sparse_matrix.cc
index 7778b00db5f26ca596bc805d7af7e64c3ffefc57..46469bd6fec358d4de0e887dba14e99794335594 100644 (file)
--- a/deal.II/lac/source/trilinos_block_sparse_matrix.cc
+++ b/deal.II/lac/source/trilinos_block_sparse_matrix.cc
@@ -13,6 +13,7 @@
 
 #include <lac/trilinos_block_sparse_matrix.h>
 
+#include <lac/block_sparse_matrix.h>
 #include <lac/block_sparsity_pattern.h>
 
 #ifdef DEAL_II_USE_TRILINOS

diff --git a/deal.II/lac/source/trilinos_precondition.cc b/deal.II/lac/source/trilinos_precondition.cc
index 5019b3fc8c3b8a21e300f84ef05bca2f90a3de83..7db151a005af7e1476fd4bde24970d0840e0273b 100755 (executable)
--- a/deal.II/lac/source/trilinos_precondition.cc
+++ b/deal.II/lac/source/trilinos_precondition.cc
@@ -12,6 +12,8 @@
 //---------------------------------------------------------------------------
 
 
+#include <lac/vector.h>
+#include <lac/sparse_matrix.h>
 #include <lac/trilinos_precondition.h>
 #include <lac/trilinos_vector_base.h>
 #include <lac/trilinos_sparse_matrix.h>

diff --git a/deal.II/lac/source/trilinos_precondition_block.cc b/deal.II/lac/source/trilinos_precondition_block.cc
index e313aa660daccad77d453ecabd6df122b977d07a..c3700b545eb7f00555199e96479ecffeb8ccefb9 100755 (executable)
--- a/deal.II/lac/source/trilinos_precondition_block.cc
+++ b/deal.II/lac/source/trilinos_precondition_block.cc
@@ -13,6 +13,9 @@
 
 #include <lac/trilinos_precondition_block.h>
 
+#include <lac/trilinos_sparse_matrix.h>
+#include <lac/trilinos_block_sparse_matrix.h>
+
 #ifdef DEAL_II_USE_TRILINOS
 
 #include <Thyra_EpetraLinearOp.hpp>

diff --git a/deal.II/lac/source/trilinos_sparse_matrix.cc b/deal.II/lac/source/trilinos_sparse_matrix.cc
index 7cde3d2772b105c76c57cf122e6375231d0126f2..4aa9ea3dce6bdcacf6a3a553ad3385238b923cb2 100755 (executable)
--- a/deal.II/lac/source/trilinos_sparse_matrix.cc
+++ b/deal.II/lac/source/trilinos_sparse_matrix.cc
@@ -13,6 +13,7 @@
 
 #include <lac/trilinos_sparse_matrix.h>
 
+#include <lac/sparse_matrix.h>
 #include <lac/trilinos_sparsity_pattern.h>
 #include <lac/sparsity_pattern.h>
 #include <lac/compressed_sparsity_pattern.h>