]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Use boost's small_vector in a few places. 5024/head
authorDavid Wells <wellsd2@rpi.edu>
Mon, 4 Sep 2017 21:26:35 +0000 (17:26 -0400)
committerDavid Wells <wellsd2@rpi.edu>
Mon, 4 Sep 2017 21:26:35 +0000 (17:26 -0400)
This class allows us to simplify some logic where we previously switched
between stack-allocated built-in arrays and std::vectors; now we can
just use the object's internal buffer in most cases, but if
dofs_per_cell or indices.size() is large enough we will allocate a big
enough buffer.

source/fe/fe_values.cc
source/grid/manifold.cc
source/lac/trilinos_sparse_matrix.cc

index 9e4fc261c44ddfdb4f173c943f6688d5d70f54c6..4892d82e31973ea956969503582a6410d175bb40 100644 (file)
@@ -42,6 +42,8 @@
 
 #include <deal.II/differentiation/ad/sacado_product_types.h>
 
+#include <boost/container/small_vector.hpp>
+
 DEAL_II_NAMESPACE_OPEN
 
 
@@ -3096,22 +3098,10 @@ void FEValuesBase<dim,spacedim>::get_function_values (
   AssertDimension (fe->n_components(), 1);
   AssertDimension (indices.size(), dofs_per_cell);
 
-  // avoid allocation when the local size is small enough
-  if (dofs_per_cell <= 100)
-    {
-      Number dof_values[100];
-      for (unsigned int i=0; i<dofs_per_cell; ++i)
-        dof_values[i] = get_vector_element (fe_function, indices[i]);
-      internal::do_function_values(&dof_values[0], this->finite_element_output.shape_values, values);
-    }
-  else
-    {
-      Vector<Number> dof_values(dofs_per_cell);
-      for (unsigned int i=0; i<dofs_per_cell; ++i)
-        dof_values[i] = get_vector_element (fe_function, indices[i]);
-      internal::do_function_values(dof_values.begin(), this->finite_element_output.shape_values,
-                                   values);
-    }
+  boost::container::small_vector<Number, 200> dof_values(dofs_per_cell);
+  for (unsigned int i=0; i<dofs_per_cell; ++i)
+    dof_values[i] = get_vector_element (fe_function, indices[i]);
+  internal::do_function_values(dof_values.data(), this->finite_element_output.shape_values, values);
 }
 
 
@@ -3156,24 +3146,12 @@ void FEValuesBase<dim,spacedim>::get_function_values (
           ExcAccessToUninitializedField("update_values"));
 
   VectorSlice<std::vector<Vector<Number> > > val(values);
-  if (indices.size() <= 100)
-    {
-      Number dof_values[100];
-      for (unsigned int i=0; i<dofs_per_cell; ++i)
-        dof_values[i] = get_vector_element (fe_function, indices[i]);
-      internal::do_function_values(&dof_values[0], this->finite_element_output.shape_values, *fe,
-                                   this->finite_element_output.shape_function_to_row_table, val,
-                                   false, indices.size()/dofs_per_cell);
-    }
-  else
-    {
-      Vector<Number> dof_values(100);
-      for (unsigned int i=0; i<dofs_per_cell; ++i)
-        dof_values[i] = get_vector_element (fe_function, indices[i]);
-      internal::do_function_values(dof_values.begin(), this->finite_element_output.shape_values, *fe,
-                                   this->finite_element_output.shape_function_to_row_table, val,
-                                   false, indices.size()/dofs_per_cell);
-    }
+  boost::container::small_vector<Number, 200> dof_values(dofs_per_cell);
+  for (unsigned int i=0; i<dofs_per_cell; ++i)
+    dof_values[i] = get_vector_element (fe_function, indices[i]);
+  internal::do_function_values(dof_values.data(), this->finite_element_output.shape_values, *fe,
+                               this->finite_element_output.shape_function_to_row_table, val,
+                               false, indices.size()/dofs_per_cell);
 }
 
 
@@ -3195,26 +3173,13 @@ void FEValuesBase<dim,spacedim>::get_function_values (
   Assert (indices.size() % dofs_per_cell == 0,
           ExcNotMultiple(indices.size(), dofs_per_cell));
 
-  if (indices.size() <= 100)
-    {
-      Number dof_values[100];
-      for (unsigned int i=0; i<indices.size(); ++i)
-        dof_values[i] = get_vector_element (fe_function, indices[i]);
-      internal::do_function_values(&dof_values[0], this->finite_element_output.shape_values, *fe,
-                                   this->finite_element_output.shape_function_to_row_table, values,
-                                   quadrature_points_fastest,
-                                   indices.size()/dofs_per_cell);
-    }
-  else
-    {
-      Vector<Number> dof_values(indices.size());
-      for (unsigned int i=0; i<indices.size(); ++i)
-        dof_values[i] = get_vector_element (fe_function, indices[i]);
-      internal::do_function_values(dof_values.begin(), this->finite_element_output.shape_values, *fe,
-                                   this->finite_element_output.shape_function_to_row_table, values,
-                                   quadrature_points_fastest,
-                                   indices.size()/dofs_per_cell);
-    }
+  boost::container::small_vector<Number, 200> dof_values(indices.size());
+  for (unsigned int i=0; i<indices.size(); ++i)
+    dof_values[i] = get_vector_element (fe_function, indices[i]);
+  internal::do_function_values(dof_values.data(), this->finite_element_output.shape_values, *fe,
+                               this->finite_element_output.shape_function_to_row_table, values,
+                               quadrature_points_fastest,
+                               indices.size()/dofs_per_cell);
 }
 
 
@@ -3255,22 +3220,12 @@ void FEValuesBase<dim,spacedim>::get_function_gradients (
           ExcAccessToUninitializedField("update_gradients"));
   AssertDimension (fe->n_components(), 1);
   AssertDimension (indices.size(), dofs_per_cell);
-  if (dofs_per_cell <= 100)
-    {
-      Number dof_values[100];
-      for (unsigned int i=0; i<dofs_per_cell; ++i)
-        dof_values[i] = get_vector_element (fe_function, indices[i]);
-      internal::do_function_derivatives(&dof_values[0], this->finite_element_output.shape_gradients,
-                                        gradients);
-    }
-  else
-    {
-      Vector<Number> dof_values(dofs_per_cell);
-      for (unsigned int i=0; i<dofs_per_cell; ++i)
-        dof_values[i] = get_vector_element (fe_function, indices[i]);
-      internal::do_function_derivatives(dof_values.begin(), this->finite_element_output.shape_gradients,
-                                        gradients);
-    }
+
+  boost::container::small_vector<Number, 200> dof_values(dofs_per_cell);
+  for (unsigned int i=0; i<dofs_per_cell; ++i)
+    dof_values[i] = get_vector_element (fe_function, indices[i]);
+  internal::do_function_derivatives(dof_values.data(), this->finite_element_output.shape_gradients,
+                                    gradients);
 }
 
 
@@ -3317,26 +3272,13 @@ void FEValuesBase<dim,spacedim>::get_function_gradients (
   Assert (this->update_flags & update_gradients,
           ExcAccessToUninitializedField("update_gradients"));
 
-  if (indices.size() <= 100)
-    {
-      Number dof_values[100];
-      for (unsigned int i=0; i<indices.size(); ++i)
-        dof_values[i] = get_vector_element (fe_function, indices[i]);
-      internal::do_function_derivatives(&dof_values[0], this->finite_element_output.shape_gradients,
-                                        *fe, this->finite_element_output.shape_function_to_row_table,
-                                        gradients, quadrature_points_fastest,
-                                        indices.size()/dofs_per_cell);
-    }
-  else
-    {
-      Vector<Number> dof_values(indices.size());
-      for (unsigned int i=0; i<indices.size(); ++i)
-        dof_values[i] = get_vector_element (fe_function, indices[i]);
-      internal::do_function_derivatives(dof_values.begin(),this->finite_element_output.shape_gradients,
-                                        *fe, this->finite_element_output.shape_function_to_row_table,
-                                        gradients, quadrature_points_fastest,
-                                        indices.size()/dofs_per_cell);
-    }
+  boost::container::small_vector<Number, 200> dof_values(indices.size());
+  for (unsigned int i=0; i<indices.size(); ++i)
+    dof_values[i] = get_vector_element (fe_function, indices[i]);
+  internal::do_function_derivatives(dof_values.data(), this->finite_element_output.shape_gradients,
+                                    *fe, this->finite_element_output.shape_function_to_row_table,
+                                    gradients, quadrature_points_fastest,
+                                    indices.size()/dofs_per_cell);
 }
 
 
@@ -3377,22 +3319,12 @@ void FEValuesBase<dim,spacedim>::get_function_hessians (
           ExcAccessToUninitializedField("update_hessians"));
   AssertDimension (fe_function.size(), present_cell->n_dofs_for_dof_handler());
   AssertDimension (indices.size(), dofs_per_cell);
-  if (dofs_per_cell <= 100)
-    {
-      Number dof_values[100];
-      for (unsigned int i=0; i<dofs_per_cell; ++i)
-        dof_values[i] = get_vector_element (fe_function, indices[i]);
-      internal::do_function_derivatives(&dof_values[0], this->finite_element_output.shape_hessians,
-                                        hessians);
-    }
-  else
-    {
-      Vector<Number> dof_values(dofs_per_cell);
-      for (unsigned int i=0; i<dofs_per_cell; ++i)
-        dof_values[i] = get_vector_element (fe_function, indices[i]);
-      internal::do_function_derivatives(dof_values.begin(), this->finite_element_output.shape_hessians,
-                                        hessians);
-    }
+
+  boost::container::small_vector<Number, 200> dof_values(dofs_per_cell);
+  for (unsigned int i=0; i<dofs_per_cell; ++i)
+    dof_values[i] = get_vector_element (fe_function, indices[i]);
+  internal::do_function_derivatives(dof_values.data(), this->finite_element_output.shape_hessians,
+                                    hessians);
 }
 
 
@@ -3437,26 +3369,14 @@ void FEValuesBase<dim, spacedim>::get_function_hessians (
           ExcAccessToUninitializedField("update_hessians"));
   Assert (indices.size() % dofs_per_cell == 0,
           ExcNotMultiple(indices.size(), dofs_per_cell));
-  if (indices.size() <= 100)
-    {
-      Number dof_values[100];
-      for (unsigned int i=0; i<indices.size(); ++i)
-        dof_values[i] = get_vector_element (fe_function, indices[i]);
-      internal::do_function_derivatives(&dof_values[0], this->finite_element_output.shape_hessians,
-                                        *fe, this->finite_element_output.shape_function_to_row_table,
-                                        hessians, quadrature_points_fastest,
-                                        indices.size()/dofs_per_cell);
-    }
-  else
-    {
-      Vector<Number> dof_values(indices.size());
-      for (unsigned int i=0; i<indices.size(); ++i)
-        dof_values[i] = get_vector_element (fe_function, indices[i]);
-      internal::do_function_derivatives(dof_values.begin(),this->finite_element_output.shape_hessians,
-                                        *fe, this->finite_element_output.shape_function_to_row_table,
-                                        hessians, quadrature_points_fastest,
-                                        indices.size()/dofs_per_cell);
-    }
+
+  boost::container::small_vector<Number, 200> dof_values(indices.size());
+  for (unsigned int i=0; i<indices.size(); ++i)
+    dof_values[i] = get_vector_element (fe_function, indices[i]);
+  internal::do_function_derivatives(dof_values.data(), this->finite_element_output.shape_hessians,
+                                    *fe, this->finite_element_output.shape_function_to_row_table,
+                                    hessians, quadrature_points_fastest,
+                                    indices.size()/dofs_per_cell);
 }
 
 
@@ -3496,22 +3416,12 @@ void FEValuesBase<dim,spacedim>::get_function_laplacians (
           ExcAccessToUninitializedField("update_hessians"));
   AssertDimension (fe->n_components(), 1);
   AssertDimension (indices.size(), dofs_per_cell);
-  if (dofs_per_cell <= 100)
-    {
-      Number dof_values[100];
-      for (unsigned int i=0; i<dofs_per_cell; ++i)
-        dof_values[i] = get_vector_element (fe_function, indices[i]);
-      internal::do_function_laplacians(&dof_values[0], this->finite_element_output.shape_hessians,
-                                       laplacians);
-    }
-  else
-    {
-      Vector<Number> dof_values(dofs_per_cell);
-      for (unsigned int i=0; i<dofs_per_cell; ++i)
-        dof_values[i] = get_vector_element (fe_function, indices[i]);
-      internal::do_function_laplacians(dof_values.begin(), this->finite_element_output.shape_hessians,
-                                       laplacians);
-    }
+
+  boost::container::small_vector<Number, 200> dof_values(dofs_per_cell);
+  for (unsigned int i=0; i<dofs_per_cell; ++i)
+    dof_values[i] = get_vector_element (fe_function, indices[i]);
+  internal::do_function_laplacians(dof_values.data(), this->finite_element_output.shape_hessians,
+                                   laplacians);
 }
 
 
@@ -3553,26 +3463,14 @@ void FEValuesBase<dim,spacedim>::get_function_laplacians (
           ExcNotMultiple(indices.size(), dofs_per_cell));
   Assert (this->update_flags & update_hessians,
           ExcAccessToUninitializedField("update_hessians"));
-  if (indices.size() <= 100)
-    {
-      Number dof_values[100];
-      for (unsigned int i=0; i<indices.size(); ++i)
-        dof_values[i] = get_vector_element (fe_function, indices[i]);
-      internal::do_function_laplacians(&dof_values[0], this->finite_element_output.shape_hessians,
-                                       *fe, this->finite_element_output.shape_function_to_row_table,
-                                       laplacians, false,
-                                       indices.size()/dofs_per_cell);
-    }
-  else
-    {
-      Vector<Number> dof_values(indices.size());
-      for (unsigned int i=0; i<indices.size(); ++i)
-        dof_values[i] = get_vector_element (fe_function, indices[i]);
-      internal::do_function_laplacians(dof_values.begin(),this->finite_element_output.shape_hessians,
-                                       *fe, this->finite_element_output.shape_function_to_row_table,
-                                       laplacians, false,
-                                       indices.size()/dofs_per_cell);
-    }
+
+  boost::container::small_vector<Number, 200> dof_values(indices.size());
+  for (unsigned int i=0; i<indices.size(); ++i)
+    dof_values[i] = get_vector_element (fe_function, indices[i]);
+  internal::do_function_laplacians(dof_values.data(), this->finite_element_output.shape_hessians,
+                                   *fe, this->finite_element_output.shape_function_to_row_table,
+                                   laplacians, false,
+                                   indices.size()/dofs_per_cell);
 }
 
 
@@ -3590,26 +3488,14 @@ void FEValuesBase<dim,spacedim>::get_function_laplacians (
           ExcNotMultiple(indices.size(), dofs_per_cell));
   Assert (this->update_flags & update_hessians,
           ExcAccessToUninitializedField("update_hessians"));
-  if (indices.size() <= 100)
-    {
-      Number dof_values[100];
-      for (unsigned int i=0; i<indices.size(); ++i)
-        dof_values[i] = get_vector_element (fe_function, indices[i]);
-      internal::do_function_laplacians(&dof_values[0], this->finite_element_output.shape_hessians,
-                                       *fe, this->finite_element_output.shape_function_to_row_table,
-                                       laplacians, quadrature_points_fastest,
-                                       indices.size()/dofs_per_cell);
-    }
-  else
-    {
-      Vector<Number> dof_values(indices.size());
-      for (unsigned int i=0; i<indices.size(); ++i)
-        dof_values[i] = get_vector_element (fe_function, indices[i]);
-      internal::do_function_laplacians(dof_values.begin(),this->finite_element_output.shape_hessians,
-                                       *fe, this->finite_element_output.shape_function_to_row_table,
-                                       laplacians, quadrature_points_fastest,
-                                       indices.size()/dofs_per_cell);
-    }
+
+  boost::container::small_vector<Number, 200> dof_values(indices.size());
+  for (unsigned int i=0; i<indices.size(); ++i)
+    dof_values[i] = get_vector_element (fe_function, indices[i]);
+  internal::do_function_laplacians(dof_values.data(), this->finite_element_output.shape_hessians,
+                                   *fe, this->finite_element_output.shape_function_to_row_table,
+                                   laplacians, quadrature_points_fastest,
+                                   indices.size()/dofs_per_cell);
 }
 
 
@@ -3650,22 +3536,12 @@ void FEValuesBase<dim,spacedim>::get_function_third_derivatives (
           ExcAccessToUninitializedField("update_3rd_derivatives"));
   AssertDimension (fe_function.size(), present_cell->n_dofs_for_dof_handler());
   AssertDimension (indices.size(), dofs_per_cell);
-  if (dofs_per_cell <= 100)
-    {
-      Number dof_values[100];
-      for (unsigned int i=0; i<dofs_per_cell; ++i)
-        dof_values[i] = get_vector_element (fe_function, indices[i]);
-      internal::do_function_derivatives(&dof_values[0], this->finite_element_output.shape_3rd_derivatives,
-                                        third_derivatives);
-    }
-  else
-    {
-      Vector<Number> dof_values(dofs_per_cell);
-      for (unsigned int i=0; i<dofs_per_cell; ++i)
-        dof_values[i] = get_vector_element (fe_function, indices[i]);
-      internal::do_function_derivatives(dof_values.begin(), this->finite_element_output.shape_3rd_derivatives,
-                                        third_derivatives);
-    }
+
+  boost::container::small_vector<Number, 200> dof_values(dofs_per_cell);
+  for (unsigned int i=0; i<dofs_per_cell; ++i)
+    dof_values[i] = get_vector_element (fe_function, indices[i]);
+  internal::do_function_derivatives(dof_values.data(), this->finite_element_output.shape_3rd_derivatives,
+                                    third_derivatives);
 }
 
 
@@ -3710,26 +3586,14 @@ void FEValuesBase<dim, spacedim>::get_function_third_derivatives (
           ExcAccessToUninitializedField("update_3rd_derivatives"));
   Assert (indices.size() % dofs_per_cell == 0,
           ExcNotMultiple(indices.size(), dofs_per_cell));
-  if (indices.size() <= 100)
-    {
-      Number dof_values[100];
-      for (unsigned int i=0; i<indices.size(); ++i)
-        dof_values[i] = get_vector_element (fe_function, indices[i]);
-      internal::do_function_derivatives(&dof_values[0], this->finite_element_output.shape_3rd_derivatives,
-                                        *fe, this->finite_element_output.shape_function_to_row_table,
-                                        third_derivatives, quadrature_points_fastest,
-                                        indices.size()/dofs_per_cell);
-    }
-  else
-    {
-      Vector<Number> dof_values(indices.size());
-      for (unsigned int i=0; i<indices.size(); ++i)
-        dof_values[i] = get_vector_element (fe_function, indices[i]);
-      internal::do_function_derivatives(dof_values.begin(),this->finite_element_output.shape_3rd_derivatives,
-                                        *fe, this->finite_element_output.shape_function_to_row_table,
-                                        third_derivatives, quadrature_points_fastest,
-                                        indices.size()/dofs_per_cell);
-    }
+
+  boost::container::small_vector<Number, 200> dof_values(indices.size());
+  for (unsigned int i=0; i<indices.size(); ++i)
+    dof_values[i] = get_vector_element (fe_function, indices[i]);
+  internal::do_function_derivatives(dof_values.data(), this->finite_element_output.shape_3rd_derivatives,
+                                    *fe, this->finite_element_output.shape_function_to_row_table,
+                                    third_derivatives, quadrature_points_fastest,
+                                    indices.size()/dofs_per_cell);
 }
 
 
index da34a01a1cc71b134330a2639c574e654602ede8..6cc6f27786145d1382a9db400f06d7c75d753034 100644 (file)
 #include <deal.II/grid/tria_iterator.h>
 #include <deal.II/grid/tria_accessor.h>
 #include <deal.II/fe/fe_q.h>
+
 #include <cmath>
 
+#include <boost/container/small_vector.hpp>
+
 DEAL_II_NAMESPACE_OPEN
 
 using namespace Manifolds;
@@ -101,23 +104,9 @@ get_new_point (const std::vector<Point<spacedim> > &surrounding_points,
   // First sort points in the order of their weights. This is done to
   // produce unique points even if get_intermediate_points is not
   // associative (as for the SphericalManifold).
-  unsigned int permutation_short[30];
-  std::vector<unsigned int> permutation_long;
-  unsigned int *permutation;
-  if (n_points > 30)
-    {
-      permutation_long.resize(n_points);
-      permutation = &permutation_long[0];
-    }
-  else
-    permutation = &permutation_short[0];
-
-  for (unsigned int i=0; i<n_points; ++i)
-    permutation[i] = i;
-
-  std::sort(permutation,
-            permutation + n_points,
-            CompareWeights(weights));
+  boost::container::small_vector<unsigned int, 100> permutation(n_points);
+  std::iota(permutation.begin(), permutation.end(), 0);
+  std::sort(permutation.begin(), permutation.end(), CompareWeights(weights));
 
   // Now loop over points in the order of their associated weight
   Point<spacedim> p = surrounding_points[permutation[0]];
index c23306384330873ac8bec0f776347e55857bcde1..24ee2f794d5fb7a7257d2e26dff22370a0f6e26d 100644 (file)
@@ -35,6 +35,8 @@ DEAL_II_DISABLE_EXTRA_DIAGNOSTICS
 #  include <Teuchos_RCP.hpp>
 DEAL_II_ENABLE_EXTRA_DIAGNOSTICS
 
+#include <boost/container/small_vector.hpp>
+
 DEAL_II_NAMESPACE_OPEN
 
 namespace TrilinosWrappers
@@ -1395,11 +1397,8 @@ namespace TrilinosWrappers
     TrilinosScalar *col_value_ptr;
     TrilinosWrappers::types::int_type n_columns;
 
-    TrilinosScalar short_val_array[100];
-    TrilinosWrappers::types::int_type short_index_array[100];
-    std::vector<TrilinosScalar> long_val_array;
-    std::vector<TrilinosWrappers::types::int_type> long_index_array;
-
+    boost::container::small_vector<TrilinosScalar, 200> local_value_array(n_cols);
+    boost::container::small_vector<TrilinosWrappers::types::int_type, 200> local_index_array(n_cols);
 
     // If we don't elide zeros, the pointers are already available... need to
     // cast to non-const pointers as that is the format taken by Trilinos (but
@@ -1414,18 +1413,8 @@ namespace TrilinosWrappers
       {
         // Otherwise, extract nonzero values in each row and get the
         // respective indices.
-        if (n_cols > 100)
-          {
-            long_val_array.resize(n_cols);
-            long_index_array.resize(n_cols);
-            col_index_ptr = &long_index_array[0];
-            col_value_ptr = &long_val_array[0];
-          }
-        else
-          {
-            col_index_ptr = &short_index_array[0];
-            col_value_ptr = &short_val_array[0];
-          }
+        col_index_ptr = local_index_array.data();
+        col_value_ptr = local_value_array.data();
 
         n_columns = 0;
         for (size_type j=0; j<n_cols; ++j)
@@ -1586,10 +1575,8 @@ namespace TrilinosWrappers
     TrilinosScalar *col_value_ptr;
     TrilinosWrappers::types::int_type n_columns;
 
-    double short_val_array[100];
-    TrilinosWrappers::types::int_type short_index_array[100];
-    std::vector<TrilinosScalar> long_val_array;
-    std::vector<TrilinosWrappers::types::int_type> long_index_array;
+    boost::container::small_vector<TrilinosScalar, 100> local_value_array(n_cols);
+    boost::container::small_vector<TrilinosWrappers::types::int_type, 100> local_index_array(n_cols);
 
     // If we don't elide zeros, the pointers are already available... need to
     // cast to non-const pointers as that is the format taken by Trilinos (but
@@ -1608,18 +1595,8 @@ namespace TrilinosWrappers
       {
         // Otherwise, extract nonzero values in each row and the corresponding
         // index.
-        if (n_cols > 100)
-          {
-            long_val_array.resize(n_cols);
-            long_index_array.resize(n_cols);
-            col_index_ptr = &long_index_array[0];
-            col_value_ptr = &long_val_array[0];
-          }
-        else
-          {
-            col_index_ptr = &short_index_array[0];
-            col_value_ptr = &short_val_array[0];
-          }
+        col_index_ptr = local_index_array.data();
+        col_value_ptr = local_value_array.data();
 
         n_columns = 0;
         for (size_type j=0; j<n_cols; ++j)

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.