]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Move some collective add/set functions into the cc file to not make the header file...
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Mon, 6 Apr 2009 15:47:02 +0000 (15:47 +0000)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Mon, 6 Apr 2009 15:47:02 +0000 (15:47 +0000)
git-svn-id: https://svn.dealii.org/trunk@18558 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/compressed_simple_sparsity_pattern.h
deal.II/lac/include/lac/sparse_matrix.h
deal.II/lac/include/lac/sparse_matrix.templates.h
deal.II/lac/source/compressed_simple_sparsity_pattern.cc
deal.II/lac/source/sparse_matrix.inst.in

index 0d25f1dddc71c7a1bd4b77c7907aa074cf5687a1..d588f72f7be49eea30bbe67081c0a30fedb95784 100644 (file)
@@ -197,7 +197,7 @@ class CompressedSimpleSparsityPattern : public Subscriptor
     void add_entries (const unsigned int row, 
                      ForwardIterator    begin,
                      ForwardIterator    end,
-                     const bool         indices_are_sorted = false);
+                     const bool         indices_are_unique_and_sorted = false);
 
                                     /**
                                      * Check if a value at a certain
@@ -415,133 +415,6 @@ CompressedSimpleSparsityPattern::Line::add (const unsigned int j)
 
 
 
-template <typename ForwardIterator>
-inline
-void
-CompressedSimpleSparsityPattern::Line::add_entries (ForwardIterator begin,
-                                                   ForwardIterator end,
-                                                   const bool      indices_are_sorted)
-{
-  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;
-
-  // 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.
-  if (stop_size > entries.capacity())
-    entries.reserve (stop_size);
-
-  unsigned int col = *my_it;
-  std::vector<unsigned int>::iterator it, it2;
-  // insert the first element as for one entry only first check the last
-  // element (or if line is still empty)
-  if ( (entries.size()==0) || ( entries.back() < col) ) {
-    entries.push_back(col);
-    it = entries.end()-1;
-  }
-  else { 
-    // do a binary search to find the place where to insert:
-    it2 = std::lower_bound(entries.begin(), entries.end(), col); 
-
-    // If this entry is a duplicate, continue immediately Insert at the
-    // right place in the vector. Vector grows automatically to fit
-    // elements. Always doubles its size.
-    if (*it2 != col)
-      it = entries.insert(it2, col);
-    else
-      it = it2;
-  }
-
-  ++my_it;
-  // Now try to be smart and insert with bias in the direction we are
-  // walking. This has the advantage that for sorted lists, we always search
-  // in the right direction, what should decrease the work needed in here.
-  for ( ; my_it != end; ++my_it)
-    {
-      col = *my_it;
-      // need a special insertion command when we're at the end of the list
-      if (col > entries.back()) {
-       entries.push_back(col);
-       it = entries.end()-1;
-      }
-      // search to the right (preferred search direction)
-      else if (col > *it) {
-       it2 = std::lower_bound(it++, entries.end(), col);
-       if (*it2 != col)
-         it = entries.insert(it2, col);
-      }
-      // search to the left
-      else if (col < *it) {
-       it2 = std::lower_bound(entries.begin(), it, col);
-       if (*it2 != col)
-         it = entries.insert(it2, col);
-      }
-      // if we're neither larger nor smaller, then this was a duplicate and
-      // we can just continue.
-    }
-}
-
-
-
 inline
 unsigned int
 CompressedSimpleSparsityPattern::n_rows () const
index fdf0d34e96f2a2524fcd7ffff143a3d2f66eb6d2..76ddb7edae2ba6228ac0e8f2299e1591ad95715e 100644 (file)
@@ -2019,11 +2019,7 @@ unsigned int SparseMatrix<number>::n () const
 
                                        // Inline the set() and add()
                                        // functions, since they will be 
-                                        // called frequently, and the 
-                                       // compiler can optimize away 
-                                       // some unnecessary loops when
-                                       // the sizes are given at 
-                                       // compile time.
+                                        // called frequently.
 template <typename number>  
 inline
 void
@@ -2031,12 +2027,25 @@ SparseMatrix<number>::set (const unsigned int i,
                           const unsigned int j,
                           const number       value)
 {
-
   Assert (numbers::is_finite(value),
          ExcMessage("The given value is not finite but either "
                     "infinite or Not A Number (NaN)"));
 
-  set (i, 1, &j, &value, false);
+  const unsigned int index = cols->operator()(i, j);
+
+                                  // it is allowed to set elements of
+                                  // the matrix that are not part of
+                                  // the sparsity pattern, if the
+                                  // value to which we set it is zero
+  if (index == SparsityPattern::invalid_entry)
+    {
+      Assert ((index != SparsityPattern::invalid_entry) ||
+             (value == 0.),
+             ExcInvalidIndex(i, j));
+      return;
+    }
+
+  val[index] = value;
 }
 
 
@@ -2099,71 +2108,6 @@ SparseMatrix<number>::set (const unsigned int               row,
 
 
 
-template <typename number>
-template <typename number2>
-inline
-void
-SparseMatrix<number>::set (const unsigned int  row,
-                          const unsigned int  n_cols,
-                          const unsigned int *col_indices,
-                          const number2       *values,
-                          const bool          elide_zero_values)
-{
-  Assert (cols != 0, ExcNotInitialized());
-
-  unsigned int n_columns = 0;
-
-                                  // Otherwise, extract nonzero values in
-                                  // each row and pass on to the other 
-                                  // function.
-  global_indices.resize(n_cols);
-  column_values.resize(n_cols);
-
-                                  // First, search all the indices to find
-                                  // out which values we actually need to
-                                  // set.
-                                  // TODO: Could probably made more
-                                  // efficient.
-  for (unsigned int j=0; j<n_cols; ++j)
-    {
-      const number value = values[j];
-      Assert (numbers::is_finite(value),
-             ExcMessage("The given value is not finite but either "
-                        "infinite or Not A Number (NaN)"));
-
-      if (value == 0 && elide_zero_values == true)
-       continue;
-
-      const unsigned int index = cols->operator()(row, col_indices[j]);
-
-                                  // it is allowed to set elements in
-                                  // the matrix that are not part of
-                                  // the sparsity pattern, if the
-                                  // value to which we set it is zero
-      if (index == SparsityPattern::invalid_entry)
-       {
-         Assert ((index != SparsityPattern::invalid_entry) ||
-                 (value == 0.),
-                 ExcInvalidIndex(row,col_indices[j]));
-         continue;
-       }
-
-      global_indices[n_columns] = index;
-      column_values[n_columns] = value;
-      n_columns++;
-    }
-
-  const unsigned int * index_ptr = &global_indices[0];
-  const number * value_ptr = &column_values[0];
-  
-                                   // Finally, go through the index list 
-                                   // and set the elements one by one.
-  for (unsigned int j=0; j<n_columns; ++j)
-    val[*index_ptr++] = *value_ptr++;
-}
-
-
-
 template <typename number>  
 inline
 void
@@ -2171,12 +2115,28 @@ SparseMatrix<number>::add (const unsigned int i,
                           const unsigned int j,
                           const number       value)
 {
-
   Assert (numbers::is_finite(value),
          ExcMessage("The given value is not finite but either "
                     "infinite or Not A Number (NaN)"));
 
-  add (i, 1, &j, &value, false);
+  if (value == 0)
+    return;
+
+  const unsigned int index = cols->operator()(i, j);
+
+                                  // it is allowed to add elements to
+                                  // the matrix that are not part of
+                                  // the sparsity pattern, if the
+                                  // value to which we set it is zero
+  if (index == SparsityPattern::invalid_entry)
+    {
+      Assert ((index != SparsityPattern::invalid_entry) ||
+             (value == 0.),
+             ExcInvalidIndex(i, j));
+      return;
+    }
+
+  val[index] += value;
 }
 
 
@@ -2239,141 +2199,6 @@ SparseMatrix<number>::add (const unsigned int               row,
 
 
 
-template <typename number>
-template <typename number2>
-inline
-void
-SparseMatrix<number>::add (const unsigned int  row,
-                          const unsigned int  n_cols,
-                          const unsigned int *col_indices,
-                          const number2      *values,
-                          const bool          elide_zero_values,
-                          const bool          col_indices_are_sorted)
-{
-  Assert (cols != 0, ExcNotInitialized());
-
-                                  // if we have sufficiently many columns
-                                  // and sorted indices it is faster to
-                                  // just go through the column indices and
-                                  // look whether we found one, rather than
-                                  // doing a binary search for every column
-  if (col_indices_are_sorted == true)
-   if (n_cols * 8 > cols->row_length(row))
-    {
-                                  // check whether the given indices are
-                                  // really sorted
-#ifdef DEBUG
-      for (unsigned int i=1; i<n_cols; ++i)
-       Assert (col_indices[i] > col_indices[i-1], 
-               ExcMessage("Indices are not sorted."));
-#endif
-      if (cols->optimize_diagonal() == true)
-       {
-         const unsigned int * this_cols = 
-           &cols->get_column_numbers()[cols->get_rowstart_indices()[row]];
-         number * val_ptr = &val[cols->get_rowstart_indices()[row]];
-         Assert (this_cols[0] == row, ExcInternalError());
-         unsigned int counter = 1;
-         for (unsigned int i=0; i<n_cols; ++i)
-           {
-                                  // diagonal
-             if (col_indices[i] == row)
-               {
-                 if (values[i] != 0)
-                   val_ptr[0] += values[i];
-                 continue;
-               }
-             
-             Assert (col_indices[i] >= this_cols[counter], ExcInternalError());
-
-             while (this_cols[counter] < col_indices[i])
-               ++counter;
-
-             Assert (this_cols[counter] == col_indices[i] || values[i] == 0,
-                     ExcInvalidIndex(row,col_indices[i]));
-
-             if (values[i] != 0)
-               val_ptr[counter] += values[i];
-           }
-         Assert (counter < cols->row_length(row), ExcInternalError());
-       }
-      else
-       {
-         const unsigned int * this_cols = 
-           &cols->get_column_numbers()[cols->get_rowstart_indices()[row]];
-         number * val_ptr = &val[cols->get_rowstart_indices()[row]];
-         unsigned int counter = 0;
-         for (unsigned int i=0; i<n_cols; ++i)
-           {
-             Assert (col_indices[i] >= this_cols[counter], ExcInternalError());
-             
-             while (this_cols[counter] < col_indices[i])
-               ++counter;
-
-             Assert (this_cols[counter] == col_indices[i] || values[i] == 0,
-                     ExcInvalidIndex(row,col_indices[i]));
-
-             if (values[i] != 0)
-               val_ptr[counter] += values[i];
-           }
-         Assert (counter < cols->row_length(row), ExcInternalError());
-       }
-      return;
-    }
-
-  unsigned int n_columns = 0;
-
-  if (global_indices.size() < n_cols)
-    {
-      global_indices.resize(n_cols);
-      column_values.resize(n_cols);
-    }
-
-                                  // First, search all the indices to find
-                                  // out which values we actually need to
-                                  // add.
-                                  // TODO: Could probably be made more
-                                  // efficient.
-  for (unsigned int j=0; j<n_cols; ++j)
-    {
-      const number value = values[j];
-      Assert (numbers::is_finite(value),
-             ExcMessage("The given value is not finite but either "
-                        "infinite or Not A Number (NaN)"));
-
-      if (value == 0 && elide_zero_values == true)
-       continue;
-
-      const unsigned int index = cols->operator()(row, col_indices[j]);
-
-                                  // it is allowed to add elements to
-                                  // the matrix that are not part of
-                                  // the sparsity pattern, if the
-                                  // value to which we set it is zero
-      if (index == SparsityPattern::invalid_entry)
-       {
-         Assert ((index != SparsityPattern::invalid_entry) ||
-                 (value == 0.),
-                 ExcInvalidIndex(row,col_indices[j]));
-         continue;
-       }
-
-      global_indices[n_columns] = index;
-      column_values[n_columns] = value;
-      n_columns++;
-    }
-
-  const unsigned int * index_ptr = &global_indices[0];
-  const number * value_ptr = &column_values[0];
-
-                                   // Finally, go through the index list 
-                                   // and add the elements one by one.
-  for (unsigned int j=0; j<n_columns; ++j)
-    val[*index_ptr++] += *value_ptr++;
-}
-
-
-
 template <typename number>
 inline
 SparseMatrix<number> &
index 0816f00ec31bc2bda080a2ddc9f7a7b110a7b7aa..21639c9c12d018fafe4afd17966c32a3b673a038 100644 (file)
@@ -325,6 +325,201 @@ SparseMatrix<number>::add (const number factor,
 }
 
 
+
+template <typename number>
+template <typename number2>
+void
+SparseMatrix<number>::add (const unsigned int  row,
+                          const unsigned int  n_cols,
+                          const unsigned int *col_indices,
+                          const number2      *values,
+                          const bool          elide_zero_values,
+                          const bool          col_indices_are_sorted)
+{
+  Assert (cols != 0, ExcNotInitialized());
+
+                                  // if we have sufficiently many columns
+                                  // and sorted indices it is faster to
+                                  // just go through the column indices and
+                                  // look whether we found one, rather than
+                                  // doing a binary search for every column
+  if (col_indices_are_sorted == true)
+   if (n_cols * 8 > cols->row_length(row))
+    {
+                                  // check whether the given indices are
+                                  // really sorted
+#ifdef DEBUG
+      for (unsigned int i=1; i<n_cols; ++i)
+       Assert (col_indices[i] > col_indices[i-1], 
+               ExcMessage("Indices are not sorted."));
+#endif
+      if (cols->optimize_diagonal() == true)
+       {
+         const unsigned int * this_cols = 
+           &cols->get_column_numbers()[cols->get_rowstart_indices()[row]];
+         number * val_ptr = &val[cols->get_rowstart_indices()[row]];
+         Assert (this_cols[0] == row, ExcInternalError());
+         unsigned int counter = 1;
+         for (unsigned int i=0; i<n_cols; ++i)
+           {
+                                  // diagonal
+             if (col_indices[i] == row)
+               {
+                 if (values[i] != 0)
+                   val_ptr[0] += values[i];
+                 continue;
+               }
+             
+             Assert (col_indices[i] >= this_cols[counter], ExcInternalError());
+
+             while (this_cols[counter] < col_indices[i])
+               ++counter;
+
+             Assert (this_cols[counter] == col_indices[i] || values[i] == 0,
+                     ExcInvalidIndex(row,col_indices[i]));
+
+             if (values[i] != 0)
+               val_ptr[counter] += values[i];
+           }
+         Assert (counter < cols->row_length(row), ExcInternalError());
+       }
+      else
+       {
+         const unsigned int * this_cols = 
+           &cols->get_column_numbers()[cols->get_rowstart_indices()[row]];
+         number * val_ptr = &val[cols->get_rowstart_indices()[row]];
+         unsigned int counter = 0;
+         for (unsigned int i=0; i<n_cols; ++i)
+           {
+             Assert (col_indices[i] >= this_cols[counter], ExcInternalError());
+             
+             while (this_cols[counter] < col_indices[i])
+               ++counter;
+
+             Assert (this_cols[counter] == col_indices[i] || values[i] == 0,
+                     ExcInvalidIndex(row,col_indices[i]));
+
+             if (values[i] != 0)
+               val_ptr[counter] += values[i];
+           }
+         Assert (counter < cols->row_length(row), ExcInternalError());
+       }
+      return;
+    }
+
+  unsigned int n_columns = 0;
+
+  if (global_indices.size() < n_cols)
+    {
+      global_indices.resize(n_cols);
+      column_values.resize(n_cols);
+    }
+
+                                  // First, search all the indices to find
+                                  // out which values we actually need to
+                                  // add.
+  for (unsigned int j=0; j<n_cols; ++j)
+    {
+      const number value = values[j];
+      Assert (numbers::is_finite(value),
+             ExcMessage("The given value is not finite but either "
+                        "infinite or Not A Number (NaN)"));
+
+      if (value == 0 && elide_zero_values == true)
+       continue;
+
+      const unsigned int index = cols->operator()(row, col_indices[j]);
+
+                                  // it is allowed to add elements to
+                                  // the matrix that are not part of
+                                  // the sparsity pattern, if the
+                                  // value to which we set it is zero
+      if (index == SparsityPattern::invalid_entry)
+       {
+         Assert ((index != SparsityPattern::invalid_entry) ||
+                 (value == 0.),
+                 ExcInvalidIndex(row,col_indices[j]));
+         continue;
+       }
+
+      global_indices[n_columns] = index;
+      column_values[n_columns] = value;
+      n_columns++;
+    }
+
+  const unsigned int * index_ptr = &global_indices[0];
+  const number * value_ptr = &column_values[0];
+
+                                   // Finally, go through the index list 
+                                   // and add the elements one by one.
+  for (unsigned int j=0; j<n_columns; ++j)
+    val[*index_ptr++] += *value_ptr++;
+}
+
+
+
+template <typename number>
+template <typename number2>
+void
+SparseMatrix<number>::set (const unsigned int  row,
+                          const unsigned int  n_cols,
+                          const unsigned int *col_indices,
+                          const number2      *values,
+                          const bool          elide_zero_values)
+{
+  Assert (cols != 0, ExcNotInitialized());
+
+  unsigned int n_columns = 0;
+
+                                  // Otherwise, extract nonzero values in
+                                  // each row and pass on to the other 
+                                  // function.
+  global_indices.resize(n_cols);
+  column_values.resize(n_cols);
+
+                                  // First, search all the indices to find
+                                  // out which values we actually need to
+                                  // set.
+  for (unsigned int j=0; j<n_cols; ++j)
+    {
+      const number value = values[j];
+      Assert (numbers::is_finite(value),
+             ExcMessage("The given value is not finite but either "
+                        "infinite or Not A Number (NaN)"));
+
+      if (value == 0 && elide_zero_values == true)
+       continue;
+
+      const unsigned int index = cols->operator()(row, col_indices[j]);
+
+                                  // it is allowed to set elements in
+                                  // the matrix that are not part of
+                                  // the sparsity pattern, if the
+                                  // value to which we set it is zero
+      if (index == SparsityPattern::invalid_entry)
+       {
+         Assert ((index != SparsityPattern::invalid_entry) ||
+                 (value == 0.),
+                 ExcInvalidIndex(row,col_indices[j]));
+         continue;
+       }
+
+      global_indices[n_columns] = index;
+      column_values[n_columns] = value;
+      n_columns++;
+    }
+
+  const unsigned int * index_ptr = &global_indices[0];
+  const number * value_ptr = &column_values[0];
+  
+                                   // Finally, go through the index list 
+                                   // and set the elements one by one.
+  for (unsigned int j=0; j<n_columns; ++j)
+    val[*index_ptr++] = *value_ptr++;
+}
+
+
+
 template <typename number>
 template <class OutVector, class InVector>
 void
index 4a3c7024349ac8bf84532775302791a263c347c1..c3ae666ded0bb895c35693e7893270dd62b8bf8e 100644 (file)
 DEAL_II_NAMESPACE_OPEN
 
 
+
+template <typename ForwardIterator>
+void
+CompressedSimpleSparsityPattern::Line::add_entries (ForwardIterator begin,
+                                                   ForwardIterator end,
+                                                   const bool      indices_are_sorted)
+{
+  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)
+    {
+                                  // in debug mode, check whether the
+                                  // indices really are sorted.
+#ifdef DEBUG
+      {
+       ForwardIterator test = begin, test1 = begin;
+       ++test1;
+       for ( ; test1 != end; ++test, ++test1)
+         Assert (*test1 > *test, ExcInternalError());
+      }
+#endif
+
+      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 and return
+      const unsigned int new_size = it - entries.begin();
+      Assert (new_size <= stop_size, ExcInternalError());
+      entries.resize (new_size);
+      return;
+    }
+
+                                  // unsorted case or case with too few
+                                  // elements
+  ForwardIterator my_it = begin;
+
+                                  // If necessary, increase the size of the
+                                  // array.
+  if (stop_size > entries.capacity())
+    entries.reserve (stop_size);
+
+  unsigned int col = *my_it;
+  std::vector<unsigned int>::iterator it, it2;
+                                  // insert the first element as for one
+                                  // entry only first check the last
+                                  // element (or if line is still empty)
+  if ( (entries.size()==0) || ( entries.back() < col) ) {
+    entries.push_back(col);
+    it = entries.end()-1;
+  }
+  else { 
+                                  // do a binary search to find the place
+                                  // where to insert:
+    it2 = std::lower_bound(entries.begin(), entries.end(), col); 
+
+                                  // If this entry is a duplicate, continue
+                                  // immediately Insert at the right place
+                                  // in the vector. Vector grows
+                                  // automatically to fit elements. Always
+                                  // doubles its size.
+    if (*it2 != col)
+      it = entries.insert(it2, col);
+    else
+      it = it2;
+  }
+
+  ++my_it;
+                                  // Now try to be smart and insert with
+                                  // bias in the direction we are
+                                  // walking. This has the advantage that
+                                  // for sorted lists, we always search in
+                                  // the right direction, what should
+                                  // decrease the work needed in here.
+  for ( ; my_it != end; ++my_it)
+    {
+      col = *my_it;
+                                  // need a special insertion command when
+                                  // we're at the end of the list
+      if (col > entries.back()) {
+       entries.push_back(col);
+       it = entries.end()-1;
+      }
+                                  // search to the right (preferred search
+                                  // direction)
+      else if (col > *it) {
+       it2 = std::lower_bound(it++, entries.end(), col);
+       if (*it2 != col)
+         it = entries.insert(it2, col);
+      }
+                                  // search to the left
+      else if (col < *it) {
+       it2 = std::lower_bound(entries.begin(), it, col);
+       if (*it2 != col)
+         it = entries.insert(it2, col);
+      }
+                                  // if we're neither larger nor smaller,
+                                  // then this was a duplicate and we can
+                                  // just continue.
+    }
+}
+
+
+
 CompressedSimpleSparsityPattern::CompressedSimpleSparsityPattern ()
                 :
                rows(0),
@@ -238,4 +387,15 @@ CompressedSimpleSparsityPattern::n_nonzero_elements () const
   return n;
 }
 
+
+// explicit instantiations
+template void CompressedSimpleSparsityPattern::Line::add_entries(unsigned int *,
+                                                                unsigned int *,
+                                                                const bool);
+template void CompressedSimpleSparsityPattern::Line::
+add_entries(std::vector<unsigned int>::iterator,
+           std::vector<unsigned int>::iterator,
+           const bool);
+
+
 DEAL_II_NAMESPACE_CLOSE
index 9db1b9faee2c2af499e73c6fb20cbc763c3e0fe9..7766b8c87c6590c3b5853751bdfd8bd8ead49dfd 100644 (file)
@@ -29,6 +29,19 @@ for (S1, S2 : REAL_SCALARS)
 
     template void SparseMatrix<S1>::add<S2> (const S1,
                                             const SparseMatrix<S2> &);
+
+    template void SparseMatrix<S1>::add<S2> (const unsigned int,
+                                            const unsigned int,    
+                                            const unsigned int *,    
+                                            const S2 *,
+                                            const bool,    
+                                            const bool);
+
+    template void SparseMatrix<S1>::set<S2> (const unsigned int,
+                                            const unsigned int,    
+                                            const unsigned int *,    
+                                            const S2 *,
+                                            const bool);
   }
 
 

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.