From: Julius Witte Date: Fri, 13 Mar 2020 08:02:34 +0000 (+0100) Subject: Apply suggested changes X-Git-Tag: v9.2.0-rc1~349^2~2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=35ea1e20687572c3a960e98283de058c529a3710;p=dealii.git Apply suggested changes --- diff --git a/include/deal.II/matrix_free/shape_info.h b/include/deal.II/matrix_free/shape_info.h index d143e40441..26d22a0f7c 100644 --- a/include/deal.II/matrix_free/shape_info.h +++ b/include/deal.II/matrix_free/shape_info.h @@ -97,6 +97,8 @@ namespace internal template struct UnivariateShapeData { + UnivariateShapeData(); + /** * Return the memory consumption of this class in bytes. */ @@ -108,7 +110,7 @@ namespace internal * will select the most efficient algorithm based on the given element * type. */ - ElementType element_type = ElementType::tensor_general; + ElementType element_type; /** * Stores the shape values of the 1D finite element evaluated on all 1D @@ -241,18 +243,18 @@ namespace internal /** * Stores the degree of the element. */ - unsigned int fe_degree = 0; + unsigned int fe_degree; /** * Stores the number of quadrature points per dimension. */ - unsigned int n_q_points_1d = 0; + unsigned int n_q_points_1d; /** * Indicates whether the basis functions are nodal in 0 and 1, i.e., the * end points of the unit cell. */ - bool nodal_at_cell_boundaries = false; + bool nodal_at_cell_boundaries; }; diff --git a/include/deal.II/matrix_free/shape_info.templates.h b/include/deal.II/matrix_free/shape_info.templates.h index da01d69018..37c15930bb 100644 --- a/include/deal.II/matrix_free/shape_info.templates.h +++ b/include/deal.II/matrix_free/shape_info.templates.h @@ -42,6 +42,18 @@ namespace internal { namespace MatrixFreeFunctions { + // ----------------- actual UnivariateShapeData functions + // -------------------- + + template + UnivariateShapeData::UnivariateShapeData() + : element_type(tensor_general) + , fe_degree(0) + , n_q_points_1d(0) + , nodal_at_cell_boundaries(false) + {} + + // ----------------- actual ShapeInfo functions -------------------- template @@ -87,31 +99,30 @@ namespace internal // assuming isotropy of dimensions and components data.resize(1); - UnivariateShapeData *univariate_shape_data = &data.front(); + UnivariateShapeData &univariate_shape_data = data.front(); data_access.reinit(n_dimensions, n_components); - data_access.fill(univariate_shape_data); - univariate_shape_data->quadrature = quad; - univariate_shape_data->fe_degree = fe->degree; - univariate_shape_data->n_q_points_1d = quad.size(); + data_access.fill(&univariate_shape_data); + univariate_shape_data.quadrature = quad; + univariate_shape_data.fe_degree = fe->degree; + univariate_shape_data.n_q_points_1d = quad.size(); // grant write access to common univariate shape data - auto &shape_values = univariate_shape_data->shape_values; - auto &shape_gradients = univariate_shape_data->shape_gradients; - auto &shape_hessians = univariate_shape_data->shape_hessians; + auto &shape_values = univariate_shape_data.shape_values; + auto &shape_gradients = univariate_shape_data.shape_gradients; + auto &shape_hessians = univariate_shape_data.shape_hessians; auto &shape_gradients_collocation = - univariate_shape_data->shape_gradients_collocation; + univariate_shape_data.shape_gradients_collocation; auto &shape_hessians_collocation = - univariate_shape_data->shape_hessians_collocation; - auto &inverse_shape_values = univariate_shape_data->inverse_shape_values; - auto &shape_data_on_face = univariate_shape_data->shape_data_on_face; - auto &values_within_subface = - univariate_shape_data->values_within_subface; + univariate_shape_data.shape_hessians_collocation; + auto &inverse_shape_values = univariate_shape_data.inverse_shape_values; + auto &shape_data_on_face = univariate_shape_data.shape_data_on_face; + auto &values_within_subface = univariate_shape_data.values_within_subface; auto &gradients_within_subface = - univariate_shape_data->gradients_within_subface; + univariate_shape_data.gradients_within_subface; auto &hessians_within_subface = - univariate_shape_data->hessians_within_subface; + univariate_shape_data.hessians_within_subface; auto &nodal_at_cell_boundaries = - univariate_shape_data->nodal_at_cell_boundaries; + univariate_shape_data.nodal_at_cell_boundaries; const unsigned int fe_degree = fe->degree; const unsigned int n_q_points_1d = quad.size(); @@ -418,9 +429,9 @@ namespace internal } if (element_type == tensor_general && - check_1d_shapes_symmetric(*univariate_shape_data)) + check_1d_shapes_symmetric(univariate_shape_data)) { - if (check_1d_shapes_collocation(*univariate_shape_data)) + if (check_1d_shapes_collocation(univariate_shape_data)) element_type = tensor_symmetric_collocation; else element_type = tensor_symmetric; @@ -439,7 +450,7 @@ namespace internal } } else if (element_type == tensor_symmetric_plus_dg0) - check_1d_shapes_symmetric(*univariate_shape_data); + check_1d_shapes_symmetric(univariate_shape_data); nodal_at_cell_boundaries = true; for (unsigned int i = 1; i < n_dofs_1d; ++i) @@ -522,8 +533,7 @@ namespace internal } } - // TODO !!! - univariate_shape_data->element_type = this->element_type; + univariate_shape_data.element_type = this->element_type; }