]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Use a replacement version of std::lower_bound that is optimized for small arrays...
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 12 Jul 2011 12:23:28 +0000 (12:23 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 12 Jul 2011 12:23:28 +0000 (12:23 +0000)
git-svn-id: https://svn.dealii.org/trunk@23944 0785d39b-7218-0410-832d-ea1e28bc413d

16 files changed:
deal.II/include/deal.II/base/index_set.h
deal.II/include/deal.II/base/utilities.h
deal.II/include/deal.II/lac/block_matrix_base.h
deal.II/include/deal.II/lac/compressed_simple_sparsity_pattern.h
deal.II/include/deal.II/lac/constraint_matrix.templates.h
deal.II/include/deal.II/lac/precondition.h
deal.II/include/deal.II/lac/sparse_decomposition.templates.h
deal.II/include/deal.II/lac/sparse_matrix.templates.h
deal.II/include/deal.II/lac/sparsity_pattern.h
deal.II/source/grid/grid_reordering.cc
deal.II/source/lac/compressed_simple_sparsity_pattern.cc
deal.II/source/lac/constraint_matrix.cc
deal.II/source/lac/sparsity_pattern.cc
deal.II/source/multigrid/mg_tools.cc
deal.II/source/numerics/matrices.cc
deal.II/source/numerics/time_dependent.cc

index b0ac759f757d064bccebd483643945e0587ce605..47b7bddc7bff82801dfbe604b5b3aa93b92b65dc 100644 (file)
@@ -1,7 +1,7 @@
 //---------------------------------------------------------------------------
 //    $Id$
 //
-//    Copyright (C) 2009, 2010 by the deal.II authors
+//    Copyright (C) 2009, 2010, 2011 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -13,6 +13,7 @@
 #define __deal2__index_set_h
 
 #include <deal.II/base/config.h>
+#include <deal.II/base/utilities.h>
 #include <deal.II/base/exceptions.h>
 
 #include <vector>
@@ -618,7 +619,7 @@ IndexSet::add_range (const unsigned int begin,
       if (ranges.size() == 0 || begin > ranges.back().end)
        ranges.push_back(new_range);
       else
-       ranges.insert (std::lower_bound (ranges.begin(),
+       ranges.insert (Utilities::lower_bound (ranges.begin(),
                                         ranges.end(),
                                         new_range),
                       new_range);
@@ -641,7 +642,7 @@ IndexSet::add_index (const unsigned int index)
   else if (index == ranges.back().end)
     ranges.back().end++;
   else
-    ranges.insert (std::lower_bound (ranges.begin(),
+    ranges.insert (Utilities::lower_bound (ranges.begin(),
                                     ranges.end(),
                                     new_range),
                   new_range);
@@ -824,7 +825,7 @@ IndexSet::nth_index_in_set (const unsigned int n) const
     }
 
   std::vector<Range>::const_iterator 
-    p = std::lower_bound(range_begin, range_end, r,
+    p = Utilities::lower_bound(range_begin, range_end, r,
                         Range::nth_index_compare);
 
   if (p != ranges.end())
@@ -870,7 +871,7 @@ IndexSet::index_within_set (const unsigned int n) const
     }
 
   std::vector<Range>::const_iterator 
-    p = std::lower_bound(range_begin, range_end, r,
+    p = Utilities::lower_bound(range_begin, range_end, r,
                         Range::end_compare);
 
   Assert(p!=ranges.end(), ExcInternalError());
index 55940ed21539f11ec0dedfe97e4476af1ba7a2ee..0fc6374cfd12a4ca2e03b67c19d3601ad49a010d 100644 (file)
@@ -1,7 +1,7 @@
 //---------------------------------------------------------------------------
 //    $Id$
 //
-//    Copyright (C) 2005, 2006, 2007, 2008, 2009, 2010 by the deal.II authors
+//    Copyright (C) 2005, 2006, 2007, 2008, 2009, 2010, 2011 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -17,6 +17,7 @@
 
 #include <vector>
 #include <utility>
+#include <functional>
 #include <string>
 
 #if defined(DEAL_II_COMPILER_SUPPORTS_MPI) || defined(DEAL_II_USE_PETSC)
@@ -172,6 +173,68 @@ namespace Utilities
   T
   fixed_power (const T t);
 
+                                  /**
+                                   * Optimized replacement for
+                                   * <tt>std::lower_bound</tt> for
+                                   * searching within the range of
+                                   * column indices. Slashes
+                                   * execution time by
+                                   * approximately one half for the
+                                   * present application, partly
+                                   * because because the
+                                   * binary search is replaced by a
+                                   * linear search for small loop
+                                   * lengths.
+                                   *
+                                   * Another reason for this function is
+                                   * rather obscure:
+                                   * when using the GCC libstdc++
+                                   * function std::lower_bound,
+                                   * complexity is O(log(N)) as
+                                   * required. However, when using the
+                                   * debug version of the GCC libstdc++
+                                   * as we do when running the testsuite,
+                                   * then std::lower_bound tests whether
+                                   * the sequence is in fact partitioned
+                                   * with respect to the pivot 'value'
+                                   * (i.e. in essence that the sequence
+                                   * is sorted as required for binary
+                                   * search to work). However, verifying
+                                   * this means that the complexity of
+                                   * std::lower_bound jumps to O(N); we
+                                   * call this function O(N) times below,
+                                   * making the overall complexity
+                                   * O(N**2). The consequence is that a
+                                   * few tests with big meshes completely
+                                   * run off the wall time limit for
+                                   * tests and fail with the libstdc++
+                                   * debug mode
+                                   *
+                                   * This function simply makes the
+                                   * assumption that the sequence is
+                                   * sorted, and we simply don't do the
+                                   * additional check.
+                                   */
+  template<typename Iterator, typename T>
+  Iterator
+  lower_bound (Iterator  first,
+              Iterator  last,
+              const T  &val);
+  
+
+                                  /**
+                                   * The same function as above, but taking
+                                   * an argument that is used to compare
+                                   * individual elements of the sequence of
+                                   * objects pointed to by the iterators.
+                                   */
+  template<typename Iterator, typename T, typename Comp>
+  Iterator
+  lower_bound (Iterator   first,
+              Iterator   last,
+              const T   &val,
+              const Comp comp);
+    
                                   /**
                                    * Given a permutation vector (i.e. a
                                    * vector $p_0\ldots p_{N-1}$ where each
@@ -675,8 +738,168 @@ namespace Utilities
       }
   }
 
+
+
+  template<typename Iterator, typename T>
+  inline
+  Iterator
+  lower_bound (Iterator  first,
+              Iterator  last,
+              const T  &val)
+  {
+    return Utilities::lower_bound (first, last, val,
+                                  std::less<T>());
+  }  
+
+
+
+  template<typename Iterator, typename T, typename Comp>
+  inline
+  Iterator
+  lower_bound (Iterator    first,
+              Iterator    last,
+              const T    &val,
+              const Comp  comp)
+  {
+    unsigned int len = last-first;
+
+    if (len==0)
+      return first;
+
+    while (true)
+      {
+                                        // if length equals 8 or less,
+                                        // then do a rolled out
+                                        // search. use a switch without
+                                        // breaks for that and roll-out
+                                        // the loop somehow
+       if (len < 8)
+         {
+           switch (len)
+             {
+               case 7:
+                     if (!comp(*first, val))
+                       return first;
+                     ++first;
+               case 6:
+                     if (!comp(*first, val))
+                       return first;
+                     ++first;
+               case 5:
+                     if (!comp(*first, val))
+                       return first;
+                     ++first;
+               case 4:
+                     if (!comp(*first, val))
+                       return first;
+                     ++first;
+               case 3:
+                     if (!comp(*first, val))
+                       return first;
+                     ++first;
+               case 2:
+                     if (!comp(*first, val))
+                       return first;
+                     ++first;
+               case 1:
+                     if (!comp(*first, val))
+                       return first;
+                     return first+1;
+               default:
+                                                      // indices seem
+                                                      // to not be
+                                                      // sorted
+                                                      // correctly!? or
+                                                      // did len
+                                                      // become==0
+                                                      // somehow? that
+                                                      // shouln't have
+                                                      // happened
+                     Assert (false, ExcInternalError());
+             }
+         }
+
+
+
+       const unsigned int half   = len >> 1;
+       const Iterator     middle = first + half;
+
+                                        // if the value is larger than
+                                        // that pointed to by the
+                                        // middle pointer, then the
+                                        // insertion point must be
+                                        // right of it
+       if (comp(*middle, val))
+         {
+           first = middle + 1;
+           len  -= half + 1;
+         }
+       else
+         len = half;
+      }
+  }  
 }
 
+                                      /**
+                                       * A replacement for std::lower_bound
+                                       * which does a binary search for a
+                                       * position in a sorted array. The
+                                       * reason for this function is obscure:
+                                       * when using the GCC libstdc++
+                                       * function std::lower_bound,
+                                       * complexity is O(log(N)) as
+                                       * required. However, when using the
+                                       * debug version of the GCC libstdc++
+                                       * as we do when running the testsuite,
+                                       * then std::lower_bound tests whether
+                                       * the sequence is in fact partitioned
+                                       * with respect to the pivot 'value'
+                                       * (i.e. in essence that the sequence
+                                       * is sorted as required for binary
+                                       * search to work). However, verifying
+                                       * this means that the complexity of
+                                       * std::lower_bound jumps to O(N); we
+                                       * call this function O(N) times below,
+                                       * making the overall complexity
+                                       * O(N**2). The consequence is that a
+                                       * few tests with big meshes completely
+                                       * run off the wall time limit for
+                                       * tests and fail with the libstdc++
+                                       * debug mode
+                                       *
+                                       * Here, we know that the sequence is
+                                       * sorted, so we don't need the
+                                       * additional check
+                                       */
+      template<typename Iterator, typename T, typename Comp>
+      Iterator
+      my_lower_bound (Iterator first,
+                     Iterator last,
+                     const T &value,
+                     Comp    &comp)
+      {
+       unsigned int length = std::distance(first, last);
+       unsigned int half;
+       Iterator middle;
+
+       while (length > 0)
+         {
+           half = length >> 1;
+           middle = first;
+           std::advance(middle, half);
+           if (comp(*middle, value))
+             {
+               first = middle;
+               ++first;
+               length = length - half - 1;
+             }
+           else
+             length = half;
+         }
+       return first;
+      }
+    
+
 
 DEAL_II_NAMESPACE_CLOSE
 
index 9be10c4c4b0f2b4ec57c924d7c281eb882974a31..f21f5651b4a3341e540971245c0511fb20a0b72d 100644 (file)
@@ -15,6 +15,7 @@
 
 #include <deal.II/base/config.h>
 #include <deal.II/base/table.h>
+#include <deal.II/base/utilities.h>
 #include <deal.II/base/smartpointer.h>
 #include <deal.II/base/memory_consumption.h>
 #include <deal.II/lac/block_indices.h>
@@ -2147,7 +2148,7 @@ BlockMatrixBase<MatrixType>::add (const unsigned int   row,
 
       if (this->n_block_cols() > 1)
        {
-         const unsigned int * first_block = std::lower_bound (col_indices,
+         const unsigned int * first_block = Utilities::lower_bound (col_indices,
                                                               col_indices+n_cols,
                                                               this->column_block_indices.block_start(1));
 
index cf4f5ed828e84d9200ff150e65ec2f597483f220..eb6f9543e892e4d2f73c92d6b852f6f850d281e9 100644 (file)
@@ -1,7 +1,7 @@
 //---------------------------------------------------------------------------
 //    $Id$
 //
-//    Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 by the deal.II authors
+//    Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2011 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -15,6 +15,7 @@
 
 #include <deal.II/base/config.h>
 #include <deal.II/base/subscriptor.h>
+#include <deal.II/base/utilities.h>
 #include <deal.II/lac/exceptions.h>
 #include <deal.II/base/index_set.h>
 
@@ -463,9 +464,10 @@ CompressedSimpleSparsityPattern::Line::add (const unsigned int j)
 
                                   // do a binary search to find the place
                                   // where to insert:
-  std::vector<unsigned int>::iterator it = std::lower_bound(entries.begin(),
-                                                           entries.end(),
-                                                           j);
+  std::vector<unsigned int>::iterator
+    it = Utilities::lower_bound(entries.begin(),
+                                         entries.end(),
+                                         j);
 
                                   // If this entry is a duplicate, exit
                                   // immediately
index ee924b9a392bc5b804194832ee96984dfdc3c50d..f4dd01fbe37e654629210366beb9dfe50c334997 100644 (file)
@@ -1325,7 +1325,7 @@ namespace internals
     for (unsigned int i=1; i<n_active_rows; ++i)
       Assert (total_row_indices[i-1] < total_row_indices[i], ExcInternalError());
 
-    pos = std::lower_bound (total_row_indices.begin(),
+    pos = Utilities::lower_bound (total_row_indices.begin(),
                            total_row_indices.begin()+n_active_rows,
                            row_value);
     if (pos->global_row == global_row)
@@ -1429,7 +1429,7 @@ namespace internals
     for (unsigned int i=1;i<num_blocks;++i)
       {
        row_iterator first_block =
-         std::lower_bound (block_indices,
+         Utilities::lower_bound (block_indices,
                            global_rows.total_row_indices.begin()+n_active_rows,
                            Distributing(block_object.get_row_indices().block_start(i)));
        block_starts[i] = first_block - global_rows.total_row_indices.begin();
@@ -1468,7 +1468,7 @@ namespace internals
     for (unsigned int i=1;i<num_blocks;++i)
       {
        row_iterator first_block =
-         std::lower_bound (col_indices,
+         Utilities::lower_bound (col_indices,
                            row_indices.end(),
                            block_object.get_row_indices().block_start(i));
        block_starts[i] = first_block - row_indices.begin();
@@ -2122,7 +2122,7 @@ make_sorted_row_list (const std::vector<unsigned int> &local_dof_indices,
          else
            {
              std::vector<unsigned int>::iterator it =
-               std::lower_bound(active_dofs.begin(),
+               Utilities::lower_bound(active_dofs.begin(),
                                 active_dofs.end()-i+1,
                                 new_index);
              if (*it != new_index)
index d399dc58009a0ab12ead4a8811eafa0588b3bed3..914f06aca966e2449720e4426e552165ddb1d68b 100644 (file)
@@ -1,7 +1,7 @@
 //---------------------------------------------------------------------------
 //    $Id$
 //
-//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2009, 2010 by the deal.II authors
+//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2009, 2010, 2011 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -16,6 +16,7 @@
 
 #include <deal.II/base/config.h>
 #include <deal.II/base/smartpointer.h>
+#include <deal.II/base/utilities.h>
 #include <deal.II/base/template_constraints.h>
 #include <deal.II/lac/tridiagonal_matrix.h>
 #include <deal.II/lac/vector_memory.h>
@@ -1336,7 +1337,7 @@ PreconditionSSOR<MATRIX>::initialize (const MATRIX &rA,
                                       // line denotes the diagonal element,
                                       // which we need not check.
          pos_right_of_diagonal[row] =
-           std::lower_bound (&colnums[*rowstart_ptr+1],
+           Utilities::lower_bound (&colnums[*rowstart_ptr+1],
                              &colnums[*(rowstart_ptr+1)],
                              row)
            - colnums;
index 37ad5836752cb0c83b66a26670090b89ebfee92a..bb07c301236aa1b83775abf096014d2afc877064 100644 (file)
@@ -15,6 +15,7 @@
 
 #include <deal.II/base/memory_consumption.h>
 #include <deal.II/base/template_constraints.h>
+#include <deal.II/base/utilities.h>
 #include <deal.II/lac/sparse_decomposition.h>
 #include <algorithm>
 #include <cstring>
@@ -170,7 +171,7 @@ SparseLUDecomposition<number>::prebuild_lower_bound()
 
   for(unsigned int row=0; row<N; row++) {
     prebuilt_lower_bound[row]
-      = std::lower_bound (&column_numbers[rowstart_indices[row]+1],
+      = Utilities::lower_bound (&column_numbers[rowstart_indices[row]+1],
                           &column_numbers[rowstart_indices[row+1]],
                           row);
   }
index 75db988f2998f4f0a6dd7102340be939697ce9aa..e73595b042a90cde6827de1ea52b7f4f09ae5630 100644 (file)
@@ -1,7 +1,7 @@
 //---------------------------------------------------------------------------
 //    $Id$
 //
-//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010 by the deal.II authors
+//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010, 2011 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -19,6 +19,7 @@
 #include <deal.II/base/parallel.h>
 #include <deal.II/base/thread_management.h>
 #include <deal.II/base/multithread_info.h>
+#include <deal.II/base/utilities.h>
 #include <deal.II/lac/sparse_matrix.h>
 #include <deal.II/lac/vector.h>
 #include <deal.II/lac/full_matrix.h>
@@ -472,9 +473,9 @@ SparseMatrix<number>::add (const unsigned int  row,
                                   // find diagonal and add it if found
          Assert (this_cols[0] == row, ExcInternalError());
          const unsigned int * diag_pos =
-           internals::SparsityPatternTools::optimized_lower_bound (col_indices,
-                                                                   col_indices+n_cols,
-                                                                   row);
+           Utilities::lower_bound (col_indices,
+                                             col_indices+n_cols,
+                                             row);
          const unsigned int diag = diag_pos - col_indices;
          unsigned int post_diag = diag;
          if (diag != n_cols && *diag_pos == row)
@@ -1433,9 +1434,9 @@ SparseMatrix<number>::precondition_SSOR (Vector<somenumber>              &dst,
                                       // line denotes the diagonal element,
                                       // which we need not check.
       const unsigned int first_right_of_diagonal_index
-       = (internals::SparsityPatternTools::optimized_lower_bound (&cols->colnums[*rowstart_ptr+1],
-                                                                  &cols->colnums[*(rowstart_ptr+1)],
-                                                                  row)
+       = (Utilities::lower_bound (&cols->colnums[*rowstart_ptr+1],
+                                            &cols->colnums[*(rowstart_ptr+1)],
+                                            row)
           -
           &cols->colnums[0]);
 
@@ -1461,9 +1462,9 @@ SparseMatrix<number>::precondition_SSOR (Vector<somenumber>              &dst,
     {
       const unsigned int end_row = *(rowstart_ptr+1);
       const unsigned int first_right_of_diagonal_index
-       = (internals::SparsityPatternTools::optimized_lower_bound (&cols->colnums[*rowstart_ptr+1],
-                                                                  &cols->colnums[end_row],
-                                                                  static_cast<unsigned int>(row)) -
+       = (Utilities::lower_bound (&cols->colnums[*rowstart_ptr+1],
+                                            &cols->colnums[end_row],
+                                            static_cast<unsigned int>(row)) -
           &cols->colnums[0]);
       number s = 0;
       for (unsigned int j=first_right_of_diagonal_index; j<end_row; ++j)
index 0d944dd717f8e47f16c4b944e359d8d690904084..b5defef65103e6f548d812faefcf63d18f3d2493 100644 (file)
@@ -44,115 +44,6 @@ namespace internals
 {
   namespace SparsityPatternTools
   {
-                                    /**
-                                     * Optimized replacement for
-                                     * <tt>std::lower_bound</tt> for
-                                     * searching within the range of
-                                     * column indices. Slashes
-                                     * execution time by
-                                     * approximately one half for the
-                                     * present application, partly
-                                     * because we have eliminated
-                                     * templates and the compiler
-                                     * seems to be able to optimize
-                                     * better, and partly because the
-                                     * binary search is replaced by a
-                                     * linear search for small loop
-                                     * lengths.
-                                     */
-    static
-    inline
-    const unsigned int *
-    optimized_lower_bound (const unsigned int *first,
-                          const unsigned int *last,
-                          const unsigned int &val)
-    {
-                                      // this function is mostly copied
-                                      // over from the STL __lower_bound
-                                      // function, but with template args
-                                      // replaced by the actual data
-                                      // types needed here, and above all
-                                      // with a rolled out search on the
-                                      // last four elements
-      unsigned int len = last-first;
-
-      if (len==0)
-       return first;
-
-      while (true)
-       {
-                                          // if length equals 8 or less,
-                                          // then do a rolled out
-                                          // search. use a switch without
-                                          // breaks for that and roll-out
-                                          // the loop somehow
-         if (len < 8)
-           {
-             switch (len)
-               {
-                 case 7:
-                       if (*first >= val)
-                         return first;
-                       ++first;
-                 case 6:
-                       if (*first >= val)
-                         return first;
-                       ++first;
-                 case 5:
-                       if (*first >= val)
-                         return first;
-                       ++first;
-                 case 4:
-                       if (*first >= val)
-                         return first;
-                       ++first;
-                 case 3:
-                       if (*first >= val)
-                         return first;
-                       ++first;
-                 case 2:
-                       if (*first >= val)
-                         return first;
-                       ++first;
-                 case 1:
-                       if (*first >= val)
-                         return first;
-                       return first+1;
-                 default:
-                                                        // indices seem
-                                                        // to not be
-                                                        // sorted
-                                                        // correctly!? or
-                                                        // did len
-                                                        // become==0
-                                                        // somehow? that
-                                                        // shouln't have
-                                                        // happened
-                       Assert (false, ExcInternalError());
-               }
-           }
-
-
-
-         const unsigned int         half   = len >> 1;
-         const unsigned int * const middle = first + half;
-
-                                          // if the value is larger than
-                                          // that pointed to by the
-                                          // middle pointer, then the
-                                          // insertion point must be
-                                          // right of it
-         if (*middle < val)
-           {
-             first = middle + 1;
-             len  -= half + 1;
-           }
-         else
-           len = half;
-       }
-    }
-
-
                                     /**
                                      * Helper function to get the
                                      * column index from a
index 2b885f513d1a5997e1af0f831468a1bf25ab0a2f..c5e29b7a6c087c70e57d9dca007a2a9946676427 100644 (file)
@@ -14,6 +14,7 @@
 #include <deal.II/grid/grid_reordering.h>
 #include <deal.II/grid/grid_reordering_internal.h>
 #include <deal.II/grid/grid_tools.h>
+#include <deal.II/base/utilities.h>
 
 #include <algorithm>
 #include <set>
@@ -268,66 +269,6 @@ namespace internal
 
     namespace
     {
-                                      /**
-                                       * A replacement for std::lower_bound
-                                       * which does a binary search for a
-                                       * position in a sorted array. The
-                                       * reason for this function is obscure:
-                                       * when using the GCC libstdc++
-                                       * function std::lower_bound,
-                                       * complexity is O(log(N)) as
-                                       * required. However, when using the
-                                       * debug version of the GCC libstdc++
-                                       * as we do when running the testsuite,
-                                       * then std::lower_bound tests whether
-                                       * the sequence is in fact partitioned
-                                       * with respect to the pivot 'value'
-                                       * (i.e. in essence that the sequence
-                                       * is sorted as required for binary
-                                       * search to work). However, verifying
-                                       * this means that the complexity of
-                                       * std::lower_bound jumps to O(N); we
-                                       * call this function O(N) times below,
-                                       * making the overall complexity
-                                       * O(N**2). The consequence is that a
-                                       * few tests with big meshes completely
-                                       * run off the wall time limit for
-                                       * tests and fail with the libstdc++
-                                       * debug mode
-                                       *
-                                       * Here, we know that the sequence is
-                                       * sorted, so we don't need the
-                                       * additional check
-                                       */
-      template<typename Iterator, typename T, typename Comp>
-      Iterator
-      my_lower_bound (Iterator first,
-                     Iterator last,
-                     const T &value,
-                     Comp comp)
-      {
-       unsigned int length = std::distance(first, last);
-       unsigned int half;
-       Iterator middle;
-
-       while (length > 0)
-         {
-           half = length >> 1;
-           middle = first;
-           std::advance(middle, half);
-           if (comp(*middle, value))
-             {
-               first = middle;
-               ++first;
-               length = length - half - 1;
-             }
-           else
-             length = half;
-         }
-       return first;
-      }
-    
-
                                       /**
                                        * Create an MQuad object from the
                                        * indices of the four vertices by
@@ -348,10 +289,10 @@ namespace internal
                                  numbers::invalid_unsigned_int };
 
        for (unsigned int i=0; i<4; ++i)
-         edges[i] = (my_lower_bound (elist.begin(),
-                                     elist.end(),
-                                     quadside(q,i),
-                                     MSide::SideSortLess())
+         edges[i] = (Utilities::lower_bound (elist.begin(),
+                                                       elist.end(),
+                                                       quadside(q,i),
+                                                       MSide::SideSortLess())
                      -
                      elist.begin());
 
index 876724206c052603ec4e07a4b24bdd7c971d5e95..4611beea82308e81932155c4312793221b73ef49 100644 (file)
@@ -62,7 +62,7 @@ CompressedSimpleSparsityPattern::Line::add_entries (ForwardIterator begin,
       ForwardIterator my_it = begin;
       unsigned int col = *my_it;
       std::vector<unsigned int>::iterator it =
-       std::lower_bound(entries.begin(), entries.end(), col);
+       Utilities::lower_bound(entries.begin(), entries.end(), col);
       while (*it == col)
        {
          ++my_it;
@@ -80,7 +80,7 @@ CompressedSimpleSparsityPattern::Line::add_entries (ForwardIterator begin,
            continue;
                                   // ok, it wasn't the very next one, do a
                                   // binary search to find the insert point
-         it = std::lower_bound(it, entries.end(), col);
+         it = Utilities::lower_bound(it, entries.end(), col);
          if (it == entries.end())
            break;
        }
@@ -152,7 +152,7 @@ CompressedSimpleSparsityPattern::Line::add_entries (ForwardIterator begin,
   else {
                                   // do a binary search to find the place
                                   // where to insert:
-    it2 = std::lower_bound(entries.begin(), entries.end(), col);
+    it2 = Utilities::lower_bound(entries.begin(), entries.end(), col);
 
                                   // If this entry is a duplicate, continue
                                   // immediately Insert at the right place
@@ -184,13 +184,13 @@ CompressedSimpleSparsityPattern::Line::add_entries (ForwardIterator begin,
                                   // search to the right (preferred search
                                   // direction)
       else if (col > *it) {
-       it2 = std::lower_bound(it++, entries.end(), col);
+       it2 = Utilities::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);
+       it2 = Utilities::lower_bound(entries.begin(), it, col);
        if (*it2 != col)
          it = entries.insert(it2, col);
       }
index 09c327f9ea5eda4e6b45b30ef3c0442b851634b9..5f9eb71580494051eb74b476075982c32081403f 100644 (file)
@@ -1858,11 +1858,11 @@ ConstraintMatrix::distribute (TrilinosWrappers::MPI::Vector &vec) const
   ConstraintLine index_comparison;
   index_comparison.line = vec.local_range().first;
   const constraint_iterator begin_my_constraints =
-    std::lower_bound (lines.begin(),lines.end(),index_comparison);
+    Utilities::lower_bound (lines.begin(),lines.end(),index_comparison);
 
   index_comparison.line = vec.local_range().second;
   const constraint_iterator end_my_constraints
-    = std::lower_bound(lines.begin(),lines.end(),index_comparison);
+    = Utilities::lower_bound(lines.begin(),lines.end(),index_comparison);
 
                                   // Here we search all the indices that we
                                   // need to have read-access to - the
@@ -1940,13 +1940,13 @@ ConstraintMatrix::distribute (TrilinosWrappers::MPI::BlockVector &vec) const
       index_comparison.line = vec.block(block).local_range().first
        +vec.get_block_indices().block_start(block);
       const constraint_iterator begin_my_constraints =
-       std::lower_bound (lines.begin(),lines.end(),index_comparison);
+       Utilities::lower_bound (lines.begin(),lines.end(),index_comparison);
 
       index_comparison.line = vec.block(block).local_range().second
        +vec.get_block_indices().block_start(block);
 
       const constraint_iterator end_my_constraints
-       = std::lower_bound(lines.begin(),lines.end(),index_comparison);
+       = Utilities::lower_bound(lines.begin(),lines.end(),index_comparison);
 
                                   // Here we search all the indices that we
                                   // need to have read-access to - the local
@@ -1995,13 +1995,13 @@ ConstraintMatrix::distribute (TrilinosWrappers::MPI::BlockVector &vec) const
       index_comparison.line = vec.block(block).local_range().first
        +vec.get_block_indices().block_start(block);
       const constraint_iterator begin_my_constraints =
-       std::lower_bound (lines.begin(),lines.end(),index_comparison);
+       Utilities::lower_bound (lines.begin(),lines.end(),index_comparison);
 
       index_comparison.line = vec.block(block).local_range().second
        +vec.get_block_indices().block_start(block);
 
       const constraint_iterator end_my_constraints
-       = std::lower_bound(lines.begin(),lines.end(),index_comparison);
+       = Utilities::lower_bound(lines.begin(),lines.end(),index_comparison);
 
       for (constraint_iterator it = begin_my_constraints;
           it != end_my_constraints; ++it)
@@ -2040,11 +2040,11 @@ ConstraintMatrix::distribute (PETScWrappers::MPI::Vector &vec) const
   ConstraintLine index_comparison;
   index_comparison.line = vec.local_range().first;
   const constraint_iterator begin_my_constraints =
-    std::lower_bound (lines.begin(),lines.end(),index_comparison);
+    Utilities::lower_bound (lines.begin(),lines.end(),index_comparison);
 
   index_comparison.line = vec.local_range().second;
   const constraint_iterator end_my_constraints
-    = std::lower_bound(lines.begin(),lines.end(),index_comparison);
+    = Utilities::lower_bound(lines.begin(),lines.end(),index_comparison);
 
                                   // all indices we need to read from
   IndexSet my_indices (vec.size());
index a0d34c0b5eb746ae50c81cd3128099968c0a9c6c..d1d09ed4a9c14fbf18968012a27116156c0a2def 100644 (file)
@@ -12,6 +12,7 @@
 //---------------------------------------------------------------------------
 
 #include <deal.II/base/vector_slice.h>
+#include <deal.II/base/utilities.h>
 #include <deal.II/lac/sparsity_pattern.h>
 #include <deal.II/lac/sparsity_tools.h>
 #include <deal.II/lac/full_matrix.h>
@@ -183,7 +184,7 @@ SparsityPattern::SparsityPattern (const SparsityPattern &original,
       const unsigned int * const
        original_last_before_side_diagonals
        = (row > extra_off_diagonals ?
-          std::lower_bound (original_row_start,
+          Utilities::lower_bound (original_row_start,
                             original_row_end,
                             row-extra_off_diagonals) :
           original_row_start);
@@ -734,13 +735,13 @@ SparsityPattern::operator () (const unsigned int i,
                                   // top of this function, so it may
                                   // not be called for noncompressed
                                   // structures.
-  const unsigned int * const sorted_region_start = (diagonal_optimized ?
-                                                   &colnums[rowstart[i]+1] :
-                                                   &colnums[rowstart[i]]);
+  const unsigned int * sorted_region_start = (diagonal_optimized ?
+                                             &colnums[rowstart[i]+1] :
+                                             &colnums[rowstart[i]]);
   const unsigned int * const p
-    = internals::SparsityPatternTools::optimized_lower_bound (sorted_region_start,
-                                                             &colnums[rowstart[i+1]],
-                                                             j);
+    = Utilities::lower_bound<const unsigned int *> (sorted_region_start,
+                                       &colnums[rowstart[i+1]],
+                                       j);
   if ((p != &colnums[rowstart[i+1]])  &&  (*p == j))
     return (p - &colnums[0]);
   else
index 99892d46b9cb1db774e3f2d5fc5f88d29ecc04c0..1f440b920179262586a4ff12fcfc0ba0b786dbd3 100644 (file)
@@ -1710,7 +1710,7 @@ MGTools::apply_boundary_values (
                                               // element
                                               // (row,dof_number)
              const unsigned int *
-               p = std::lower_bound(&sparsity_colnums[sparsity_rowstart[row]+1],
+               p = Utilities::lower_bound(&sparsity_colnums[sparsity_rowstart[row]+1],
                                     &sparsity_colnums[sparsity_rowstart[row+1]],
                                     dof_number);
 
@@ -1957,14 +1957,14 @@ MGTools::apply_boundary_values (
                        p = &this_sparsity.get_column_numbers()
                            [this_sparsity.get_rowstart_indices()[row]];
                      else
-                       p = std::lower_bound(&this_sparsity.get_column_numbers()
+                       p = Utilities::lower_bound(&this_sparsity.get_column_numbers()
                                             [this_sparsity.get_rowstart_indices()[row]+1],
                                             &this_sparsity.get_column_numbers()
                                             [this_sparsity.get_rowstart_indices()[row+1]],
                                             block_index.second);
                    }
                  else
-                   p = std::lower_bound(&this_sparsity.get_column_numbers()
+                   p = Utilities::lower_bound(&this_sparsity.get_column_numbers()
                                         [this_sparsity.get_rowstart_indices()[row]],
                                         &this_sparsity.get_column_numbers()
                                         [this_sparsity.get_rowstart_indices()[row+1]],
index 553be095af31fe45ae119b2e3eca4be18dd04dad..3d7a1e2dcaeb58b22dda0da23a5329cb2434fe2c 100644 (file)
@@ -1993,7 +1993,7 @@ MatrixTools::apply_boundary_values (const std::map<unsigned int,double> &boundar
                                               // element
                                               // (row,dof_number)
              const unsigned int *
-               p = Utilities::optimized_lower_bound(&sparsity_colnums[sparsity_rowstart[row]+1],
+               p = Utilities::lower_bound(&sparsity_colnums[sparsity_rowstart[row]+1],
                                     &sparsity_colnums[sparsity_rowstart[row+1]],
                                     dof_number);
 
@@ -2291,14 +2291,14 @@ MatrixTools::apply_boundary_values (const std::map<unsigned int,double> &boundar
                        p = &this_sparsity.get_column_numbers()
                            [this_sparsity.get_rowstart_indices()[row]];
                      else
-                       p = Utilities::optimized_lower_bound(&this_sparsity.get_column_numbers()
+                       p = Utilities::lower_bound(&this_sparsity.get_column_numbers()
                                             [this_sparsity.get_rowstart_indices()[row]+1],
                                             &this_sparsity.get_column_numbers()
                                             [this_sparsity.get_rowstart_indices()[row+1]],
                                             block_index.second);
                    }
                  else
-                   p = Utilities::optimized_lower_bound(&this_sparsity.get_column_numbers()
+                   p = Utilities::lower_bound(&this_sparsity.get_column_numbers()
                                         [this_sparsity.get_rowstart_indices()[row]],
                                         &this_sparsity.get_column_numbers()
                                         [this_sparsity.get_rowstart_indices()[row+1]],
index 7d4e4f49b8a3bc9449285560d217e2717ed792e1..f4e5ccfcf0eb64c14844a9db2fd30b9bdc3e3e27 100644 (file)
@@ -15,6 +15,7 @@
 #include <deal.II/numerics/time_dependent.h>
 #include <deal.II/base/memory_consumption.h>
 #include <deal.II/base/thread_management.h>
+#include <deal.II/base/utilities.h>
 #include <deal.II/grid/tria.h>
 #include <deal.II/grid/tria_accessor.h>
 #include <deal.II/grid/tria_iterator.h>
@@ -766,7 +767,7 @@ void TimeStepBase_Tria<dim>::refine_grid (const RefinementData refinement_data)
       sorted_criteria = criteria;
       std::sort (sorted_criteria.begin(),
                 sorted_criteria.end());
-      p_refinement_threshold = std::lower_bound (sorted_criteria.begin(),
+      p_refinement_threshold = Utilities::lower_bound (sorted_criteria.begin(),
                                                 sorted_criteria.end(),
                                                 static_cast<float>(refinement_threshold));
       p_coarsening_threshold = std::upper_bound (sorted_criteria.begin(),

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.