From 0bd7098169cbdf548811624fc11a68e7dfe93b40 Mon Sep 17 00:00:00 2001 From: David Wells Date: Fri, 23 Feb 2018 15:06:29 -0500 Subject: [PATCH] Deprecate fixed_int_power. We can use constexpr functions instead of recursive templates. --- include/deal.II/base/utilities.h | 33 +++++++++++++- .../deal.II/matrix_free/cuda_fe_evaluation.h | 6 +-- .../matrix_free/cuda_tensor_product_kernels.h | 6 +-- .../deal.II/matrix_free/evaluation_kernels.h | 24 +++++----- include/deal.II/matrix_free/fe_evaluation.h | 4 +- include/deal.II/matrix_free/operators.h | 4 +- .../matrix_free/tensor_product_kernels.h | 44 +++++++++---------- 7 files changed, 73 insertions(+), 48 deletions(-) diff --git a/include/deal.II/base/utilities.h b/include/deal.II/base/utilities.h index 717fd9662f..d6009340f5 100644 --- a/include/deal.II/base/utilities.h +++ b/include/deal.II/base/utilities.h @@ -302,9 +302,24 @@ namespace Utilities * * Use this class as in fixed_int_power@<5,2@>::%value to * compute 52. + * + * @deprecated This template has been deprecated in favor of C++11's support + * for constexpr calculations, e.g., use + * + * @code + * constexpr int value = Utilities::pow(2, dim); + * @endcode + * + * instead of + * + * @code + * const int value = Utilities::fixed_int_power<2, dim>::value; + * @endcode + * + * to obtain a constant expression for value. */ template - struct fixed_int_power + struct DEAL_II_DEPRECATED fixed_int_power { static const int value = a *fixed_int_power::value; }; @@ -312,13 +327,27 @@ namespace Utilities /** * Base case for the power operation with N=0, which gives the * result 1. + * + * @deprecated This template is deprecated: see the note in the general + * version of this template for more information. */ template - struct fixed_int_power + struct DEAL_II_DEPRECATED fixed_int_power { static const int value = 1; }; + /** + * A replacement for std::pow that allows compile-time + * calculations for constant expression arguments. + */ + constexpr + unsigned int + pow(const unsigned int base, const unsigned int iexp) + { + return iexp == 0 ? 1 : base*dealii::Utilities::pow(base, iexp - 1); + } + /** * Optimized replacement for std::lower_bound for searching within * the range of column indices. Slashes execution time by approximately one diff --git a/include/deal.II/matrix_free/cuda_fe_evaluation.h b/include/deal.II/matrix_free/cuda_fe_evaluation.h index 3465d7fefb..27bfb8d513 100644 --- a/include/deal.II/matrix_free/cuda_fe_evaluation.h +++ b/include/deal.II/matrix_free/cuda_fe_evaluation.h @@ -76,10 +76,8 @@ namespace CUDAWrappers typedef typename MatrixFree::Data data_type; static const unsigned int dimension = dim; static const unsigned int n_components = n_components_; - static const unsigned int n_q_points = - Utilities::fixed_int_power::value; - static const unsigned int tensor_dofs_per_cell = - Utilities::fixed_int_power::value; + static constexpr unsigned int n_q_points = Utilities::pow(n_q_points_1d, dim); + static constexpr unsigned int tensor_dofs_per_cell = Utilities::pow(fe_degree + 1, dim); /** * Constructor. diff --git a/include/deal.II/matrix_free/cuda_tensor_product_kernels.h b/include/deal.II/matrix_free/cuda_tensor_product_kernels.h index b2452d9434..d9a3558570 100644 --- a/include/deal.II/matrix_free/cuda_tensor_product_kernels.h +++ b/include/deal.II/matrix_free/cuda_tensor_product_kernels.h @@ -64,10 +64,8 @@ namespace CUDAWrappers struct EvaluatorTensorProduct { - static const unsigned int dofs_per_cell = - Utilities::fixed_int_power::value; - static const unsigned int n_q_points = - Utilities::fixed_int_power::value; + static constexpr unsigned int dofs_per_cell = Utilities::pow(fe_degree + 1, dim); + static constexpr unsigned int n_q_points = Utilities::pow(n_q_points_1d, dim); __device__ EvaluatorTensorProduct(); diff --git a/include/deal.II/matrix_free/evaluation_kernels.h b/include/deal.II/matrix_free/evaluation_kernels.h index 7f69067eb7..d299486286 100644 --- a/include/deal.II/matrix_free/evaluation_kernels.h +++ b/include/deal.II/matrix_free/evaluation_kernels.h @@ -783,12 +783,12 @@ namespace internal for (unsigned int q=basis_size_1; q!=0; --q) FEEvaluationImplBasisChange ::do_forward(transformation_matrix, - values_in + (q-1)*Utilities::fixed_int_power::value, - my_scratch + (q-1)*Utilities::fixed_int_power::value); + values_in + (q-1)*Utilities::pow(basis_size_1, dim-1), + my_scratch + (q-1)*Utilities::pow(basis_size_2, dim-1)); EvaluatorTensorProduct eval_val (transformation_matrix); const unsigned int n_inner_blocks = (dim > 1 && basis_size_2 < 10) ? basis_size_2 : 1; - const unsigned int n_blocks = Utilities::fixed_int_power::value; + const unsigned int n_blocks = Utilities::pow(basis_size_2, dim-1); for (unsigned int ii=0; ii ::do_backward(transformation_matrix, false, - my_scratch + q*Utilities::fixed_int_power::value, - values_out + q*Utilities::fixed_int_power::value); + my_scratch + q*Utilities::pow(basis_size_2, dim-1), + values_out + q*Utilities::pow(basis_size_1, dim-1)); } }; @@ -872,7 +872,7 @@ namespace internal eval(AlignedVector(), shape_info.shape_gradients_collocation_eo, shape_info.shape_hessians_collocation_eo); - constexpr unsigned int n_q_points = Utilities::fixed_int_power::value; + constexpr unsigned int n_q_points = Utilities::pow(fe_degree+1, dim); for (unsigned int c=0; c(), shape_info.shape_gradients_collocation_eo, shape_info.shape_hessians_collocation_eo); - constexpr unsigned int n_q_points = Utilities::fixed_int_power::value; + constexpr unsigned int n_q_points = Utilities::pow(fe_degree+1, dim); for (unsigned int c=0; c::value; + constexpr unsigned int n_q_points = Utilities::pow(n_q_points_1d, dim); for (unsigned int c=0; c::value; + constexpr unsigned int n_q_points = Utilities::pow(n_q_points_1d, dim); for (unsigned int c=0; c -1 ? - Utilities::fixed_int_power::value : + Utilities::pow(fe_degree+1, dim-1) : (dim > 1 ? Utilities::fixed_power(data.fe_degree+1) : 1); const unsigned int n_q_points = fe_degree > -1 ? - Utilities::fixed_int_power::value : data.n_q_points_face; + Utilities::pow(n_q_points_1d, dim-1) : data.n_q_points_face; if (evaluate_grad == false) for (unsigned int c=0; c -1 ? - Utilities::fixed_int_power::value : + Utilities::pow(fe_degree+1, dim-1) : (dim > 1 ? Utilities::fixed_power(data.fe_degree+1) : 1); const unsigned int n_q_points = fe_degree > -1 ? diff --git a/include/deal.II/matrix_free/fe_evaluation.h b/include/deal.II/matrix_free/fe_evaluation.h index adc1247592..d4e83961fb 100644 --- a/include/deal.II/matrix_free/fe_evaluation.h +++ b/include/deal.II/matrix_free/fe_evaluation.h @@ -1912,8 +1912,8 @@ public: typedef typename BaseClass::gradient_type gradient_type; static constexpr unsigned int dimension = dim; static constexpr unsigned int n_components = n_components_; - static constexpr unsigned int static_n_q_points = Utilities::fixed_int_power::value; - static constexpr unsigned int static_dofs_per_component = Utilities::fixed_int_power::value; + static constexpr unsigned int static_n_q_points = Utilities::pow(n_q_points_1d, dim); + static constexpr unsigned int static_dofs_per_component = Utilities::pow(fe_degree + 1, dim); static constexpr unsigned int tensor_dofs_per_cell = static_dofs_per_component *n_components; static constexpr unsigned int static_dofs_per_cell = static_dofs_per_component *n_components; diff --git a/include/deal.II/matrix_free/operators.h b/include/deal.II/matrix_free/operators.h index 88e0bbd0fb..ced0beb93e 100644 --- a/include/deal.II/matrix_free/operators.h +++ b/include/deal.II/matrix_free/operators.h @@ -851,7 +851,7 @@ namespace MatrixFreeOperators CellwiseInverseMassMatrix ::fill_inverse_JxW_values(AlignedVector > &inverse_jxw) const { - const unsigned int dofs_per_cell = Utilities::fixed_int_power::value; + constexpr unsigned int dofs_per_cell = Utilities::pow(fe_degree + 1, dim); Assert(inverse_jxw.size() > 0 && inverse_jxw.size() % dofs_per_cell == 0, ExcMessage("Expected diagonal to be a multiple of scalar dof per cells")); @@ -883,7 +883,7 @@ namespace MatrixFreeOperators const VectorizedArray *in_array, VectorizedArray *out_array) const { - const unsigned int dofs_per_cell = Utilities::fixed_int_power::value; + constexpr unsigned int dofs_per_cell = Utilities::pow(fe_degree + 1, dim); Assert(inverse_coefficients.size() > 0 && inverse_coefficients.size() % dofs_per_cell == 0, ExcMessage("Expected diagonal to be a multiple of scalar dof per cells")); diff --git a/include/deal.II/matrix_free/tensor_product_kernels.h b/include/deal.II/matrix_free/tensor_product_kernels.h index 91b3f21bc0..cae328832f 100644 --- a/include/deal.II/matrix_free/tensor_product_kernels.h +++ b/include/deal.II/matrix_free/tensor_product_kernels.h @@ -113,8 +113,8 @@ namespace internal template struct EvaluatorTensorProduct { - static const unsigned int n_rows_of_product = Utilities::fixed_int_power::value; - static const unsigned int n_columns_of_product = Utilities::fixed_int_power::value; + static constexpr unsigned int n_rows_of_product = Utilities::pow(n_rows, dim); + static constexpr unsigned int n_columns_of_product = Utilities::pow(n_columns, dim); /** * Empty constructor. Does nothing. Be careful when using 'values' and @@ -298,9 +298,9 @@ namespace internal constexpr int mm = contract_over_rows ? n_rows : n_columns, nn = contract_over_rows ? n_columns : n_rows; - constexpr int stride = Utilities::fixed_int_power::value; + constexpr int stride = Utilities::pow(n_columns, direction); constexpr int n_blocks1 = one_line ? 1 : stride; - constexpr int n_blocks2 = Utilities::fixed_int_power=dim)?0:(dim-direction-1)>::value; + constexpr int n_blocks2 = Utilities::pow(n_rows, (direction>=dim) ? 0 : (dim-direction-1)); for (int i2=0; i2 2 ? n_rows : 1; AssertIndexRange (face_direction, dim); - constexpr int stride = Utilities::fixed_int_power::value; - constexpr int out_stride = Utilities::fixed_int_power::value; + constexpr int stride = Utilities::pow(n_rows, face_direction); + constexpr int out_stride = Utilities::pow(n_rows, dim-1); const Number *DEAL_II_RESTRICT shape_values = this->shape_values; for (int i2=0; i2 struct EvaluatorTensorProduct { - static const unsigned int n_rows_of_product = Utilities::fixed_int_power::value; - static const unsigned int n_columns_of_product = Utilities::fixed_int_power::value; + static constexpr unsigned int n_rows_of_product = Utilities::pow(n_rows, dim); + static constexpr unsigned int n_columns_of_product = Utilities::pow(n_columns, dim); /** * Constructor, taking the data from ShapeInfo @@ -900,9 +900,9 @@ namespace internal constexpr int n_cols = nn / 2; constexpr int mid = mm / 2; - constexpr int stride = Utilities::fixed_int_power::value; + constexpr int stride = Utilities::pow(n_columns, direction); constexpr int n_blocks1 = stride; - constexpr int n_blocks2 = Utilities::fixed_int_power=dim)?0:(dim-direction-1)>::value; + constexpr int n_blocks2 = Utilities::pow(n_rows, (direction>=dim) ? 0 : (dim-direction-1)); for (int i2=0; i2::value; + constexpr int stride = Utilities::pow(n_columns, direction); constexpr int n_blocks1 = stride; - constexpr int n_blocks2 = Utilities::fixed_int_power=dim)?0:(dim-direction-1)>::value; + constexpr int n_blocks2 = Utilities::pow(n_rows, (direction>=dim) ? 0 : (dim-direction-1)); for (int i2=0; i2::value; + constexpr int stride = Utilities::pow(n_columns, direction); constexpr int n_blocks1 = stride; - constexpr int n_blocks2 = Utilities::fixed_int_power=dim)?0:(dim-direction-1)>::value; + constexpr int n_blocks2 = Utilities::pow(n_rows, (direction>=dim) ? 0 : (dim-direction-1)); for (int i2=0; i2 struct EvaluatorTensorProduct { - static const unsigned int n_rows_of_product = Utilities::fixed_int_power::value; - static const unsigned int n_columns_of_product = Utilities::fixed_int_power::value; + static constexpr unsigned int n_rows_of_product = Utilities::pow(n_rows, dim); + static constexpr unsigned int n_columns_of_product = Utilities::pow(n_columns, dim); /** * Empty constructor. Does nothing. Be careful when using 'values' and @@ -1531,9 +1531,9 @@ namespace internal constexpr int n_cols = nn / 2; constexpr int mid = mm / 2; - constexpr int stride = Utilities::fixed_int_power::value; + constexpr int stride = Utilities::pow(n_columns, direction); constexpr int n_blocks1 = one_line ? 1 : stride; - constexpr int n_blocks2 = Utilities::fixed_int_power=dim)?0:(dim-direction-1)>::value; + constexpr int n_blocks2 = Utilities::pow(n_rows, (direction>=dim) ? 0 : (dim-direction-1)); constexpr int offset = (n_columns+1)/2; @@ -1722,8 +1722,8 @@ namespace internal template struct EvaluatorTensorProduct { - static const unsigned int n_rows_of_product = Utilities::fixed_int_power::value; - static const unsigned int n_columns_of_product = Utilities::fixed_int_power::value; + static constexpr unsigned int n_rows_of_product = Utilities::pow(n_rows, dim); + static constexpr unsigned int n_columns_of_product = Utilities::pow(n_columns, dim); /** * Empty constructor. Does nothing. Be careful when using 'values' and @@ -1886,9 +1886,9 @@ namespace internal constexpr int n_cols = nn / 2; constexpr int mid = mm / 2; - constexpr int stride = Utilities::fixed_int_power::value; + constexpr int stride = Utilities::pow(n_columns, direction); constexpr int n_blocks1 = one_line ? 1 : stride; - constexpr int n_blocks2 = Utilities::fixed_int_power=dim)?0:(dim-direction-1)>::value; + constexpr int n_blocks2 = Utilities::pow(n_rows, (direction>=dim) ? 0 : (dim-direction-1)); // this code may look very inefficient at first sight due to the many // different cases with if's at the innermost loop part, but all of the -- 2.39.5