From: David Wells Date: Mon, 4 Sep 2017 21:26:35 +0000 (-0400) Subject: Use boost's small_vector in a few places. X-Git-Tag: v9.0.0-rc1~1122^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=7de0398a4277b6240b7d6f181f80763b74e436b1;p=dealii.git Use boost's small_vector in a few places. 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. --- diff --git a/source/fe/fe_values.cc b/source/fe/fe_values.cc index 9e4fc261c4..4892d82e31 100644 --- a/source/fe/fe_values.cc +++ b/source/fe/fe_values.cc @@ -42,6 +42,8 @@ #include +#include + DEAL_II_NAMESPACE_OPEN @@ -3096,22 +3098,10 @@ void FEValuesBase::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; ifinite_element_output.shape_values, values); - } - else - { - Vector dof_values(dofs_per_cell); - for (unsigned int i=0; ifinite_element_output.shape_values, - values); - } + boost::container::small_vector dof_values(dofs_per_cell); + for (unsigned int i=0; ifinite_element_output.shape_values, values); } @@ -3156,24 +3146,12 @@ void FEValuesBase::get_function_values ( ExcAccessToUninitializedField("update_values")); VectorSlice > > val(values); - if (indices.size() <= 100) - { - Number dof_values[100]; - for (unsigned int i=0; ifinite_element_output.shape_values, *fe, - this->finite_element_output.shape_function_to_row_table, val, - false, indices.size()/dofs_per_cell); - } - else - { - Vector dof_values(100); - for (unsigned int i=0; ifinite_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 dof_values(dofs_per_cell); + for (unsigned int i=0; ifinite_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::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; ifinite_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 dof_values(indices.size()); - for (unsigned int i=0; ifinite_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 dof_values(indices.size()); + for (unsigned int i=0; ifinite_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::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; ifinite_element_output.shape_gradients, - gradients); - } - else - { - Vector dof_values(dofs_per_cell); - for (unsigned int i=0; ifinite_element_output.shape_gradients, - gradients); - } + + boost::container::small_vector dof_values(dofs_per_cell); + for (unsigned int i=0; ifinite_element_output.shape_gradients, + gradients); } @@ -3317,26 +3272,13 @@ void FEValuesBase::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; ifinite_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 dof_values(indices.size()); - for (unsigned int i=0; ifinite_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 dof_values(indices.size()); + for (unsigned int i=0; ifinite_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::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; ifinite_element_output.shape_hessians, - hessians); - } - else - { - Vector dof_values(dofs_per_cell); - for (unsigned int i=0; ifinite_element_output.shape_hessians, - hessians); - } + + boost::container::small_vector dof_values(dofs_per_cell); + for (unsigned int i=0; ifinite_element_output.shape_hessians, + hessians); } @@ -3437,26 +3369,14 @@ void FEValuesBase::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; ifinite_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 dof_values(indices.size()); - for (unsigned int i=0; ifinite_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 dof_values(indices.size()); + for (unsigned int i=0; ifinite_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::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; ifinite_element_output.shape_hessians, - laplacians); - } - else - { - Vector dof_values(dofs_per_cell); - for (unsigned int i=0; ifinite_element_output.shape_hessians, - laplacians); - } + + boost::container::small_vector dof_values(dofs_per_cell); + for (unsigned int i=0; ifinite_element_output.shape_hessians, + laplacians); } @@ -3553,26 +3463,14 @@ void FEValuesBase::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; ifinite_element_output.shape_hessians, - *fe, this->finite_element_output.shape_function_to_row_table, - laplacians, false, - indices.size()/dofs_per_cell); - } - else - { - Vector dof_values(indices.size()); - for (unsigned int i=0; ifinite_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 dof_values(indices.size()); + for (unsigned int i=0; ifinite_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::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; ifinite_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 dof_values(indices.size()); - for (unsigned int i=0; ifinite_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 dof_values(indices.size()); + for (unsigned int i=0; ifinite_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::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; ifinite_element_output.shape_3rd_derivatives, - third_derivatives); - } - else - { - Vector dof_values(dofs_per_cell); - for (unsigned int i=0; ifinite_element_output.shape_3rd_derivatives, - third_derivatives); - } + + boost::container::small_vector dof_values(dofs_per_cell); + for (unsigned int i=0; ifinite_element_output.shape_3rd_derivatives, + third_derivatives); } @@ -3710,26 +3586,14 @@ void FEValuesBase::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; ifinite_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 dof_values(indices.size()); - for (unsigned int i=0; ifinite_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 dof_values(indices.size()); + for (unsigned int i=0; ifinite_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); } diff --git a/source/grid/manifold.cc b/source/grid/manifold.cc index da34a01a1c..6cc6f27786 100644 --- a/source/grid/manifold.cc +++ b/source/grid/manifold.cc @@ -20,8 +20,11 @@ #include #include #include + #include +#include + DEAL_II_NAMESPACE_OPEN using namespace Manifolds; @@ -101,23 +104,9 @@ get_new_point (const std::vector > &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 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 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 p = surrounding_points[permutation[0]]; diff --git a/source/lac/trilinos_sparse_matrix.cc b/source/lac/trilinos_sparse_matrix.cc index c233063843..24ee2f794d 100644 --- a/source/lac/trilinos_sparse_matrix.cc +++ b/source/lac/trilinos_sparse_matrix.cc @@ -35,6 +35,8 @@ DEAL_II_DISABLE_EXTRA_DIAGNOSTICS # include DEAL_II_ENABLE_EXTRA_DIAGNOSTICS +#include + 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 long_val_array; - std::vector long_index_array; - + boost::container::small_vector local_value_array(n_cols); + boost::container::small_vector 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 long_val_array; - std::vector long_index_array; + boost::container::small_vector local_value_array(n_cols); + boost::container::small_vector 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