From: Martin Kronbichler Date: Tue, 14 Dec 2021 16:04:05 +0000 (+0100) Subject: Make many variable names more consistent X-Git-Tag: v9.4.0-rc1~740^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F13056%2Fhead;p=dealii.git Make many variable names more consistent --- diff --git a/include/deal.II/matrix_free/evaluation_kernels.h b/include/deal.II/matrix_free/evaluation_kernels.h index 72327cd070..c85113a867 100644 --- a/include/deal.II/matrix_free/evaluation_kernels.h +++ b/include/deal.II/matrix_free/evaluation_kernels.h @@ -124,13 +124,13 @@ namespace internal evaluate(const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, const Number * values_dofs_actual, - FEEvaluationData & eval); + FEEvaluationData & fe_eval); static void integrate(const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, Number * values_dofs_actual, - FEEvaluationData & eval, + FEEvaluationData & fe_eval, const bool add_into_values_array); static Eval @@ -170,13 +170,13 @@ namespace internal evaluate(const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, const Number * values_dofs_actual, - FEEvaluationData & eval); + FEEvaluationData & fe_eval); static void integrate(const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, Number * values_dofs_actual, - FEEvaluationData & eval, + FEEvaluationData & fe_eval, const bool add_into_values_array); }; @@ -192,7 +192,7 @@ namespace internal const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, const Number * values_dofs_actual, - FEEvaluationData & eval) + FEEvaluationData & fe_eval) { if (evaluation_flag == EvaluationFlags::nothing) return; @@ -200,9 +200,9 @@ namespace internal std::array *, 3> univariate_shape_data; - const auto &shape_data = eval.get_shape_info().data; + const auto &shape_data = fe_eval.get_shape_info().data; - univariate_shape_data.fill(&shape_data[0]); + univariate_shape_data.fill(&shape_data.front()); if (shape_data.size() == dim) for (int i = 1; i < dim; ++i) @@ -218,14 +218,14 @@ namespace internal (Eval::n_rows_of_product > Eval::n_columns_of_product ? Eval::n_rows_of_product : Eval::n_columns_of_product); - Number *temp1 = eval.get_scratch_data().begin(); + Number *temp1 = fe_eval.get_scratch_data().begin(); Number *temp2; if (temp_size == 0) { - temp2 = - temp1 + - std::max(Utilities::fixed_power(shape_data[0].fe_degree + 1), - Utilities::fixed_power(shape_data[0].n_q_points_1d)); + temp2 = temp1 + std::max(Utilities::fixed_power( + shape_data.front().fe_degree + 1), + Utilities::fixed_power( + shape_data.front().n_q_points_1d)); } else { @@ -233,21 +233,21 @@ namespace internal } const std::size_t n_q_points = temp_size == 0 ? - eval.get_shape_info().n_q_points : + fe_eval.get_shape_info().n_q_points : Eval::n_columns_of_product; const std::size_t dofs_per_comp = (type == MatrixFreeFunctions::truncated_tensor) ? - Utilities::pow(shape_data[0].fe_degree + 1, dim) : - eval.get_shape_info().dofs_per_component_on_cell; + Utilities::pow(shape_data.front().fe_degree + 1, dim) : + fe_eval.get_shape_info().dofs_per_component_on_cell; const Number *values_dofs = values_dofs_actual; if (type == MatrixFreeFunctions::truncated_tensor) { const std::size_t n_dofs_per_comp = - eval.get_shape_info().dofs_per_component_on_cell; + fe_eval.get_shape_info().dofs_per_component_on_cell; Number *values_dofs_tmp = temp1 + 2 * (std::max(n_dofs_per_comp, n_q_points)); const int degree = - fe_degree != -1 ? fe_degree : shape_data[0].fe_degree; + fe_degree != -1 ? fe_degree : shape_data.front().fe_degree; for (unsigned int c = 0; c < n_components; ++c) for (int i = 0, count_p = 0, count_q = 0; i < (dim > 2 ? degree + 1 : 1); @@ -270,9 +270,9 @@ namespace internal values_dofs = values_dofs_tmp; } - Number *values_quad = eval.begin_values(); - Number *gradients_quad = eval.begin_gradients(); - Number *hessians_quad = eval.begin_hessians(); + Number *values_quad = fe_eval.begin_values(); + Number *gradients_quad = fe_eval.begin_gradients(); + Number *hessians_quad = fe_eval.begin_hessians(); switch (dim) { @@ -465,14 +465,14 @@ namespace internal const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, Number * values_dofs_actual, - FEEvaluationData & eval, + FEEvaluationData & fe_eval, const bool add_into_values_array) { std::array *, 3> univariate_shape_data; - const auto &shape_data = eval.get_shape_info().data; - univariate_shape_data.fill(&shape_data[0]); + const auto &shape_data = fe_eval.get_shape_info().data; + univariate_shape_data.fill(&shape_data.front()); if (shape_data.size() == dim) for (int i = 1; i < dim; ++i) @@ -488,14 +488,14 @@ namespace internal (Eval::n_rows_of_product > Eval::n_columns_of_product ? Eval::n_rows_of_product : Eval::n_columns_of_product); - Number *temp1 = eval.get_scratch_data().begin(); + Number *temp1 = fe_eval.get_scratch_data().begin(); Number *temp2; if (temp_size == 0) { - temp2 = - temp1 + - std::max(Utilities::fixed_power(shape_data[0].fe_degree + 1), - Utilities::fixed_power(shape_data[0].n_q_points_1d)); + temp2 = temp1 + std::max(Utilities::fixed_power( + shape_data.front().fe_degree + 1), + Utilities::fixed_power( + shape_data.front().n_q_points_1d)); } else { @@ -503,23 +503,23 @@ namespace internal } const std::size_t n_q_points = temp_size == 0 ? - eval.get_shape_info().n_q_points : + fe_eval.get_shape_info().n_q_points : Eval::n_columns_of_product; const unsigned int dofs_per_comp = (type == MatrixFreeFunctions::truncated_tensor) ? - Utilities::fixed_power(shape_data[0].fe_degree + 1) : - eval.get_shape_info().dofs_per_component_on_cell; + Utilities::fixed_power(shape_data.front().fe_degree + 1) : + fe_eval.get_shape_info().dofs_per_component_on_cell; // expand dof_values to tensor product for truncated tensor products Number *values_dofs = (type == MatrixFreeFunctions::truncated_tensor) ? - temp1 + - 2 * (std::max( - eval.get_shape_info().dofs_per_component_on_cell, n_q_points)) : + temp1 + 2 * (std::max( + fe_eval.get_shape_info().dofs_per_component_on_cell, + n_q_points)) : values_dofs_actual; - Number *values_quad = eval.begin_values(); - Number *gradients_quad = eval.begin_gradients(); - Number *hessians_quad = eval.begin_hessians(); + Number *values_quad = fe_eval.begin_values(); + Number *gradients_quad = fe_eval.begin_gradients(); + Number *hessians_quad = fe_eval.begin_hessians(); switch (dim) { @@ -745,10 +745,10 @@ namespace internal if (type == MatrixFreeFunctions::truncated_tensor) { const std::size_t n_dofs_per_comp = - eval.get_shape_info().dofs_per_component_on_cell; + fe_eval.get_shape_info().dofs_per_component_on_cell; values_dofs -= dofs_per_comp * n_components; const int degree = - fe_degree != -1 ? fe_degree : shape_data[0].fe_degree; + fe_degree != -1 ? fe_degree : shape_data.front().fe_degree; for (unsigned int c = 0; c < n_components; ++c) for (int i = 0, count_p = 0, count_q = 0; i < (dim > 2 ? degree + 1 : 1); @@ -771,28 +771,29 @@ namespace internal template inline void - FEEvaluationImpl::evaluate(const unsigned int n_components, - const EvaluationFlags::EvaluationFlags - evaluation_flag, - const Number *values_dofs_actual, - FEEvaluationData &eval) + FEEvaluationImpl< + MatrixFreeFunctions::tensor_none, + dim, + fe_degree, + n_q_points_1d, + Number>::evaluate(const unsigned int n_components, + const EvaluationFlags::EvaluationFlags evaluation_flag, + const Number * values_dofs_actual, + FEEvaluationData & fe_eval) { - const std::size_t n_dofs = eval.get_shape_info().dofs_per_component_on_cell; - const std::size_t n_q_points = eval.get_shape_info().n_q_points; + const std::size_t n_dofs = + fe_eval.get_shape_info().dofs_per_component_on_cell; + const std::size_t n_q_points = fe_eval.get_shape_info().n_q_points; - const auto &shape_data = eval.get_shape_info().data; + const auto &shape_data = fe_eval.get_shape_info().data; using Eval = EvaluatorTensorProduct; if (evaluation_flag & EvaluationFlags::values) { - const auto shape_values = shape_data[0].shape_values.data(); - auto values_quad_ptr = eval.begin_values(); + const auto shape_values = shape_data.front().shape_values.data(); + auto values_quad_ptr = fe_eval.begin_values(); auto values_dofs_actual_ptr = values_dofs_actual; Eval eval(shape_values, nullptr, nullptr, n_dofs, n_q_points); @@ -808,8 +809,8 @@ namespace internal if (evaluation_flag & EvaluationFlags::gradients) { - const auto shape_gradients = shape_data[0].shape_gradients.data(); - auto gradients_quad_ptr = eval.begin_gradients(); + const auto shape_gradients = shape_data.front().shape_gradients.data(); + auto gradients_quad_ptr = fe_eval.begin_gradients(); auto values_dofs_actual_ptr = values_dofs_actual; for (unsigned int c = 0; c < n_components; ++c) @@ -847,25 +848,26 @@ namespace internal Number>::integrate(const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, Number * values_dofs_actual, - FEEvaluationData &eval, + FEEvaluationData &fe_eval, const bool add_into_values_array) { // TODO: implement hessians AssertThrow(!(integration_flag & EvaluationFlags::hessians), ExcNotImplemented()); - const std::size_t n_dofs = eval.get_shape_info().dofs_per_component_on_cell; - const std::size_t n_q_points = eval.get_shape_info().n_q_points; + const std::size_t n_dofs = + fe_eval.get_shape_info().dofs_per_component_on_cell; + const std::size_t n_q_points = fe_eval.get_shape_info().n_q_points; - const auto &shape_data = eval.get_shape_info().data; + const auto &shape_data = fe_eval.get_shape_info().data; using Eval = EvaluatorTensorProduct; if (integration_flag & EvaluationFlags::values) { - const auto shape_values = shape_data[0].shape_values.data(); - auto values_quad_ptr = eval.begin_values(); + const auto shape_values = shape_data.front().shape_values.data(); + auto values_quad_ptr = fe_eval.begin_values(); auto values_dofs_actual_ptr = values_dofs_actual; Eval eval(shape_values, nullptr, nullptr, n_dofs, n_q_points); @@ -885,8 +887,8 @@ namespace internal if (integration_flag & EvaluationFlags::gradients) { - const auto shape_gradients = shape_data[0].shape_gradients.data(); - auto gradients_quad_ptr = eval.begin_gradients(); + const auto shape_gradients = shape_data.front().shape_gradients.data(); + auto gradients_quad_ptr = fe_eval.begin_gradients(); auto values_dofs_actual_ptr = values_dofs_actual; for (unsigned int c = 0; c < n_components; ++c) @@ -1336,7 +1338,7 @@ namespace internal evaluate(const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, const Number * values_dofs, - FEEvaluationData & eval); + FEEvaluationData & fe_eval); static void do_evaluate(const MatrixFreeFunctions::UnivariateShapeData &shape, @@ -1349,7 +1351,7 @@ namespace internal integrate(const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, Number * values_dofs, - FEEvaluationData & eval, + FEEvaluationData & fe_eval, const bool add_into_values_array); static void @@ -1369,7 +1371,7 @@ namespace internal const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, const Number * values_dofs, - FEEvaluationData & eval) + FEEvaluationData & fe_eval) { constexpr std::size_t n_points = Utilities::pow(fe_degree + 1, dim); @@ -1377,14 +1379,15 @@ namespace internal { if (evaluation_flag & EvaluationFlags::values) for (unsigned int i = 0; i < n_points; ++i) - eval.begin_values()[n_points * c + i] = + fe_eval.begin_values()[n_points * c + i] = values_dofs[n_points * c + i]; - do_evaluate(eval.get_shape_info().data[0], + do_evaluate(fe_eval.get_shape_info().data.front(), evaluation_flag, values_dofs + c * n_points, - eval.begin_gradients() + c * dim * n_points, - eval.begin_hessians() + c * dim * (dim + 1) / 2 * n_points); + fe_eval.begin_gradients() + c * dim * n_points, + fe_eval.begin_hessians() + + c * dim * (dim + 1) / 2 * n_points); } } @@ -1457,7 +1460,7 @@ namespace internal const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, Number * values_dofs, - FEEvaluationData & eval, + FEEvaluationData & fe_eval, const bool add_into_values_array) { constexpr std::size_t n_points = Utilities::pow(fe_degree + 1, dim); @@ -1469,18 +1472,19 @@ namespace internal if (add_into_values_array) for (unsigned int i = 0; i < n_points; ++i) values_dofs[n_points * c + i] += - eval.begin_values()[n_points * c + i]; + fe_eval.begin_values()[n_points * c + i]; else for (unsigned int i = 0; i < n_points; ++i) values_dofs[n_points * c + i] = - eval.begin_values()[n_points * c + i]; + fe_eval.begin_values()[n_points * c + i]; } - do_integrate(eval.get_shape_info().data[0], + do_integrate(fe_eval.get_shape_info().data.front(), integration_flag, values_dofs + c * n_points, - eval.begin_gradients() + c * dim * n_points, - eval.begin_hessians() + c * dim * (dim + 1) / 2 * n_points, + fe_eval.begin_gradients() + c * dim * n_points, + fe_eval.begin_hessians() + + c * dim * (dim + 1) / 2 * n_points, add_into_values_array || (integration_flag & EvaluationFlags::values)); } @@ -1610,13 +1614,13 @@ namespace internal evaluate(const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, const Number * values_dofs, - FEEvaluationData & eval); + FEEvaluationData & fe_eval); static void integrate(const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, Number * values_dofs, - FEEvaluationData & eval, + FEEvaluationData & fe_eval, const bool add_into_values_array); }; @@ -1631,9 +1635,9 @@ namespace internal Number>::evaluate(const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, const Number * values_dofs, - FEEvaluationData & eval) + FEEvaluationData & fe_eval) { - const auto &shape_data = eval.get_shape_info().data[0]; + const auto &shape_data = fe_eval.get_shape_info().data.front(); Assert(n_q_points_1d > fe_degree, ExcMessage("You lose information when going to a collocation space " @@ -1655,7 +1659,7 @@ namespace internal Number>::do_forward(1, shape_data.shape_values_eo, values_dofs + c * n_dofs, - eval.begin_values() + c * n_q_points); + fe_eval.begin_values() + c * n_q_points); // apply derivatives in the collocation space if (evaluation_flag & @@ -1664,9 +1668,9 @@ namespace internal do_evaluate(shape_data, evaluation_flag & (EvaluationFlags::gradients | EvaluationFlags::hessians), - eval.begin_values() + c * n_q_points, - eval.begin_gradients() + c * dim * n_q_points, - eval.begin_hessians() + + fe_eval.begin_values() + c * n_q_points, + fe_eval.begin_gradients() + c * dim * n_q_points, + fe_eval.begin_hessians() + c * dim * (dim + 1) / 2 * n_q_points); } } @@ -1682,10 +1686,10 @@ namespace internal Number>::integrate(const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, Number * values_dofs, - FEEvaluationData & eval, + FEEvaluationData & fe_eval, const bool add_into_values_array) { - const auto &shape_data = eval.get_shape_info().data[0]; + const auto &shape_data = fe_eval.get_shape_info().data.front(); Assert(n_q_points_1d > fe_degree, ExcMessage("You lose information when going to a collocation space " @@ -1703,9 +1707,9 @@ namespace internal do_integrate(shape_data, integration_flag & (EvaluationFlags::gradients | EvaluationFlags::hessians), - eval.begin_values() + c * n_q_points, - eval.begin_gradients() + c * dim * n_q_points, - eval.begin_hessians() + + fe_eval.begin_values() + c * n_q_points, + fe_eval.begin_gradients() + c * dim * n_q_points, + fe_eval.begin_hessians() + c * dim * (dim + 1) / 2 * n_q_points, /*add_into_values_array=*/ integration_flag & EvaluationFlags::values); @@ -1721,7 +1725,7 @@ namespace internal Number>::do_backward(1, shape_data.shape_values_eo, add_into_values_array, - eval.begin_values() + c * n_q_points, + fe_eval.begin_values() + c * n_q_points, values_dofs + c * Utilities::pow(fe_degree + 1, dim)); } @@ -1768,13 +1772,13 @@ namespace internal run(const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, const Number * values_dofs, - FEEvaluationData & eval) + FEEvaluationData & fe_eval) { - const auto element_type = eval.get_shape_info().element_type; + const auto element_type = fe_eval.get_shape_info().element_type; using ElementType = MatrixFreeFunctions::ElementType; - Assert(eval.get_shape_info().data.size() == 1 || - (eval.get_shape_info().data.size() == dim && + Assert(fe_eval.get_shape_info().data.size() == 1 || + (fe_eval.get_shape_info().data.size() == dim && element_type == ElementType::tensor_general), ExcNotImplemented()); @@ -1782,7 +1786,7 @@ namespace internal element_type == ElementType::tensor_symmetric_collocation) { FEEvaluationImplCollocation::evaluate( - n_components, evaluation_flag, values_dofs, eval); + n_components, evaluation_flag, values_dofs, fe_eval); } // '<=' on type means tensor_symmetric or tensor_symmetric_hermite, see // shape_info.h for more details @@ -1794,7 +1798,10 @@ namespace internal dim, fe_degree, n_q_points_1d, - Number>::evaluate(n_components, evaluation_flag, values_dofs, eval); + Number>::evaluate(n_components, + evaluation_flag, + values_dofs, + fe_eval); } else if (fe_degree >= 0 && element_type <= ElementType::tensor_symmetric) { @@ -1805,7 +1812,7 @@ namespace internal Number>::evaluate(n_components, evaluation_flag, values_dofs, - eval); + fe_eval); } else if (element_type == ElementType::tensor_symmetric_plus_dg0) { @@ -1816,7 +1823,7 @@ namespace internal Number>::evaluate(n_components, evaluation_flag, values_dofs, - eval); + fe_eval); } else if (element_type == ElementType::truncated_tensor) { @@ -1827,7 +1834,7 @@ namespace internal Number>::evaluate(n_components, evaluation_flag, values_dofs, - eval); + fe_eval); } else if (element_type == ElementType::tensor_none) { @@ -1838,7 +1845,7 @@ namespace internal Number>::evaluate(n_components, evaluation_flag, values_dofs, - eval); + fe_eval); } else { @@ -1849,7 +1856,7 @@ namespace internal Number>::evaluate(n_components, evaluation_flag, values_dofs, - eval); + fe_eval); } return false; @@ -1881,14 +1888,14 @@ namespace internal run(const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, Number * values_dofs, - FEEvaluationData & eval, + FEEvaluationData & fe_eval, const bool sum_into_values_array) { - const auto element_type = eval.get_shape_info().element_type; + const auto element_type = fe_eval.get_shape_info().element_type; using ElementType = MatrixFreeFunctions::ElementType; - Assert(eval.get_shape_info().data.size() == 1 || - (eval.get_shape_info().data.size() == dim && + Assert(fe_eval.get_shape_info().data.size() == 1 || + (fe_eval.get_shape_info().data.size() == dim && element_type == ElementType::tensor_general), ExcNotImplemented()); @@ -1899,7 +1906,7 @@ namespace internal n_components, integration_flag, values_dofs, - eval, + fe_eval, sum_into_values_array); } // '<=' on type means tensor_symmetric or tensor_symmetric_hermite, see @@ -1915,7 +1922,7 @@ namespace internal Number>::integrate(n_components, integration_flag, values_dofs, - eval, + fe_eval, sum_into_values_array); } else if (fe_degree >= 0 && element_type <= ElementType::tensor_symmetric) @@ -1927,7 +1934,7 @@ namespace internal Number>::integrate(n_components, integration_flag, values_dofs, - eval, + fe_eval, sum_into_values_array); } else if (element_type == ElementType::tensor_symmetric_plus_dg0) @@ -1939,7 +1946,7 @@ namespace internal Number>::integrate(n_components, integration_flag, values_dofs, - eval, + fe_eval, sum_into_values_array); } else if (element_type == ElementType::truncated_tensor) @@ -1951,7 +1958,7 @@ namespace internal Number>::integrate(n_components, integration_flag, values_dofs, - eval, + fe_eval, sum_into_values_array); } else if (element_type == ElementType::tensor_none) @@ -1963,7 +1970,7 @@ namespace internal Number>::integrate(n_components, integration_flag, values_dofs, - eval, + fe_eval, sum_into_values_array); } else @@ -1975,7 +1982,7 @@ namespace internal Number>::integrate(n_components, integration_flag, values_dofs, - eval, + fe_eval, sum_into_values_array); } @@ -2479,12 +2486,13 @@ namespace internal static void interpolate(const unsigned int n_components, const EvaluationFlags::EvaluationFlags flags, - const MatrixFreeFunctions::ShapeInfo &data, + const MatrixFreeFunctions::ShapeInfo &shape_info, const Number * input, Number * output, const unsigned int face_no) { - Assert(static_cast(fe_degree) == data.data[0].fe_degree || + Assert(static_cast(fe_degree) == + shape_info.data.front().fe_degree || fe_degree == -1, ExcInternalError()); @@ -2494,10 +2502,10 @@ namespace internal output, flags, face_no, - data.data[0].fe_degree + 1, - data.data[0].shape_data_on_face, - data.dofs_per_component_on_cell, - 3 * data.dofs_per_component_on_face); + shape_info.data.front().fe_degree + 1, + shape_info.data.front().shape_data_on_face, + shape_info.dofs_per_component_on_cell, + 3 * shape_info.dofs_per_component_on_face); } /** @@ -2505,15 +2513,16 @@ namespace internal */ template static void - interpolate_quadrature(const unsigned int n_components, - const EvaluationFlags::EvaluationFlags flags, - const MatrixFreeFunctions::ShapeInfo &data, - const Number * input, - Number * output, - const unsigned int face_no) + interpolate_quadrature( + const unsigned int n_components, + const EvaluationFlags::EvaluationFlags flags, + const MatrixFreeFunctions::ShapeInfo &shape_info, + const Number * input, + Number * output, + const unsigned int face_no) { Assert(static_cast(fe_degree + 1) == - data.data[0].n_q_points_1d || + shape_info.data.front().n_q_points_1d || fe_degree == -1, ExcInternalError()); @@ -2523,10 +2532,10 @@ namespace internal output, flags, face_no, - data.data.front().quadrature.size(), - data.data.front().quadrature_data_on_face, - data.n_q_points, - data.n_q_points_face); + shape_info.data.front().quadrature.size(), + shape_info.data.front().quadrature_data_on_face, + shape_info.n_q_points, + shape_info.n_q_points_face); } private: @@ -2850,17 +2859,17 @@ namespace internal run(const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, const Number * values_dofs, - FEEvaluationData & eval) + FEEvaluationData & fe_eval) { - const auto & data = eval.get_shape_info(); - const auto & shape_data = data.data.front(); - const unsigned int face_no = eval.get_face_no(); - const unsigned int face_orientation = eval.get_face_orientation(); + const auto & shape_info = fe_eval.get_shape_info(); + const auto & shape_data = shape_info.data.front(); + const unsigned int face_no = fe_eval.get_face_no(); + const unsigned int face_orientation = fe_eval.get_face_orientation(); - if (data.element_type == MatrixFreeFunctions::tensor_none) + if (shape_info.element_type == MatrixFreeFunctions::tensor_none) { - const std::size_t n_dofs = data.dofs_per_component_on_cell; - const std::size_t n_q_points = data.n_q_points_faces[face_no]; + const std::size_t n_dofs = shape_info.dofs_per_component_on_cell; + const std::size_t n_q_points = shape_info.n_q_points_faces[face_no]; using Eval = EvaluatorTensorProduct; @@ -2870,7 +2879,7 @@ namespace internal const auto shape_values = &shape_data.shape_values_face(face_no, face_orientation, 0); - auto values_quad_ptr = eval.begin_values(); + auto values_quad_ptr = fe_eval.begin_values(); auto values_dofs_actual_ptr = values_dofs; Eval eval(shape_values, nullptr, nullptr, n_dofs, n_q_points); @@ -2886,7 +2895,7 @@ namespace internal if (evaluation_flag & EvaluationFlags::gradients) { - auto gradients_quad_ptr = eval.begin_gradients(); + auto gradients_quad_ptr = fe_eval.begin_gradients(); auto values_dofs_actual_ptr = values_dofs; std::array shape_gradients; @@ -2925,20 +2934,24 @@ namespace internal // Note: we always keep storage of values, 1st and 2nd derivatives in an // array, so reserve space for all three here - Number *temp = eval.get_scratch_data().begin(); + Number *temp = fe_eval.get_scratch_data().begin(); Number *scratch_data = temp + 3 * n_components * dofs_per_face; - FEFaceNormalEvaluationImpl::template interpolate< - true, - false>(n_components, evaluation_flag, data, values_dofs, temp, face_no); - - const unsigned int subface_index = eval.get_subface_index(); + FEFaceNormalEvaluationImpl:: + template interpolate(n_components, + evaluation_flag, + shape_info, + values_dofs, + temp, + face_no); + + const unsigned int subface_index = fe_eval.get_subface_index(); constexpr unsigned int n_q_points_1d_actual = fe_degree > -1 ? n_q_points_1d : 0; if (fe_degree > -1 && subface_index >= GeometryInfo::max_children_per_cell && - data.element_type <= MatrixFreeFunctions::tensor_symmetric) + shape_info.element_type <= MatrixFreeFunctions::tensor_symmetric) FEFaceEvaluationImpl & eval) + FEEvaluationData & fe_eval) { - const auto & data = eval.get_shape_info(); - const auto & shape_data = data.data.front(); - const unsigned int face_no = eval.get_face_no(); - const unsigned int face_orientation = eval.get_face_orientation(); + const auto & shape_info = fe_eval.get_shape_info(); + const auto & shape_data = shape_info.data.front(); + const unsigned int face_no = fe_eval.get_face_no(); + const unsigned int face_orientation = fe_eval.get_face_orientation(); - if (data.element_type == MatrixFreeFunctions::tensor_none) + if (shape_info.element_type == MatrixFreeFunctions::tensor_none) { - const std::size_t n_dofs = data.dofs_per_component_on_cell; - const std::size_t n_q_points = data.n_q_points_faces[face_no]; + const std::size_t n_dofs = shape_info.dofs_per_component_on_cell; + const std::size_t n_q_points = shape_info.n_q_points_faces[face_no]; using Eval = EvaluatorTensorProduct; @@ -3014,7 +3029,7 @@ namespace internal const auto shape_values = &shape_data.shape_values_face(face_no, face_orientation, 0); - auto values_quad_ptr = eval.begin_values(); + auto values_quad_ptr = fe_eval.begin_values(); auto values_dofs_actual_ptr = values_dofs; Eval eval(shape_values, nullptr, nullptr, n_dofs, n_q_points); @@ -3030,7 +3045,7 @@ namespace internal if (integration_flag & EvaluationFlags::gradients) { - auto gradients_quad_ptr = eval.begin_gradients(); + auto gradients_quad_ptr = fe_eval.begin_gradients(); auto values_dofs_actual_ptr = values_dofs; std::array shape_gradients; @@ -3072,7 +3087,7 @@ namespace internal fe_degree > -1 ? Utilities::pow(fe_degree + 1, dim - 1) : Utilities::pow(shape_data.fe_degree + 1, dim - 1); - Number *temp = eval.get_scratch_data().begin(); + Number *temp = fe_eval.get_scratch_data().begin(); Number *scratch_data = temp + 3 * n_components * dofs_per_face; if (face_orientation) @@ -3080,53 +3095,59 @@ namespace internal dim, n_components, integration_flag, - &eval.get_shape_info().face_orientations_quad(face_orientation, 0), + &fe_eval.get_shape_info().face_orientations_quad(face_orientation, 0), true, - data.n_q_points_face, + shape_info.n_q_points_face, temp, - eval.begin_values(), - eval.begin_gradients(), - eval.begin_hessians()); + fe_eval.begin_values(), + fe_eval.begin_gradients(), + fe_eval.begin_hessians()); const unsigned int n_q_points_1d_actual = fe_degree > -1 ? n_q_points_1d : 0; - const unsigned int subface_index = eval.get_subface_index(); + const unsigned int subface_index = fe_eval.get_subface_index(); if (fe_degree > -1 && - eval.get_subface_index() >= + fe_eval.get_subface_index() >= GeometryInfo::max_children_per_cell && - data.element_type <= MatrixFreeFunctions::tensor_symmetric) - FEFaceEvaluationImpl::integrate_in_face(n_components, - integration_flag, - shape_data, - temp, - eval.begin_values(), - eval.begin_gradients(), - eval.begin_hessians(), - scratch_data, - subface_index); + shape_info.element_type <= MatrixFreeFunctions::tensor_symmetric) + FEFaceEvaluationImpl< + true, + dim, + fe_degree, + n_q_points_1d_actual, + Number>::integrate_in_face(n_components, + integration_flag, + shape_data, + temp, + fe_eval.begin_values(), + fe_eval.begin_gradients(), + fe_eval.begin_hessians(), + scratch_data, + subface_index); else - FEFaceEvaluationImpl::integrate_in_face(n_components, - integration_flag, - shape_data, - temp, - eval.begin_values(), - eval.begin_gradients(), - eval.begin_hessians(), - scratch_data, - subface_index); + FEFaceEvaluationImpl< + false, + dim, + fe_degree, + n_q_points_1d_actual, + Number>::integrate_in_face(n_components, + integration_flag, + shape_data, + temp, + fe_eval.begin_values(), + fe_eval.begin_gradients(), + fe_eval.begin_hessians(), + scratch_data, + subface_index); FEFaceNormalEvaluationImpl:: - template interpolate( - n_components, integration_flag, data, temp, values_dofs, face_no); + template interpolate(n_components, + integration_flag, + shape_info, + temp, + values_dofs, + face_no); return false; } }; @@ -3143,7 +3164,7 @@ namespace internal const EvaluationFlags::EvaluationFlags evaluation_flag, typename Processor::Number2_ * global_vector_ptr, const std::vector> *sm_ptr, - const EvaluationData & eval, + const EvaluationData & fe_eval, typename Processor::VectorizedArrayType_ * temp1) { constexpr int dim = Processor::dim_; @@ -3154,13 +3175,13 @@ namespace internal using Number = typename Processor::Number_; using Number2_ = typename Processor::Number2_; - const auto & shape_data = eval.get_shape_info().data[0]; + const auto & shape_data = fe_eval.get_shape_info().data.front(); constexpr bool integrate = Processor::do_integrate; - const unsigned int face_no = eval.get_face_no(); - const auto & dof_info = eval.get_dof_info(); - const unsigned int cell = eval.get_cell_or_face_batch_id(); + const unsigned int face_no = fe_eval.get_face_no(); + const auto & dof_info = fe_eval.get_dof_info(); + const unsigned int cell = fe_eval.get_cell_or_face_batch_id(); const MatrixFreeFunctions::DoFInfo::DoFAccessIndex dof_access_index = - eval.get_dof_access_index(); + fe_eval.get_dof_access_index(); AssertIndexRange(cell, dof_info.index_storage_variants[dof_access_index].size()); constexpr unsigned int dofs_per_face = @@ -3171,9 +3192,10 @@ namespace internal bool all_faces_are_same = n_filled_lanes == n_lanes; if (n_face_orientations == n_lanes) for (unsigned int v = 1; v < n_lanes; ++v) - if (eval.get_all_face_numbers()[v] != eval.get_all_face_numbers()[0] || - eval.get_all_face_orientations()[v] != - eval.get_all_face_orientations()[0]) + if (fe_eval.get_all_face_numbers()[v] != + fe_eval.get_all_face_numbers()[0] || + fe_eval.get_all_face_orientations()[v] != + fe_eval.get_all_face_orientations()[0]) { all_faces_are_same = false; break; @@ -3192,23 +3214,23 @@ namespace internal if (n_face_orientations == n_lanes && dof_access_index == MatrixFreeFunctions::DoFInfo::dof_access_cell && - eval.get_is_interior_face() == 0) + fe_eval.get_is_interior_face() == 0) for (unsigned int v = 0; v < n_lanes; ++v) { // the loop breaks once an invalid_unsigned_int is hit for // all cases except the exterior faces in the ECL loop (where // some faces might be at the boundaries but others not) - if (eval.get_cell_ids()[v] == numbers::invalid_unsigned_int) + if (fe_eval.get_cell_ids()[v] == numbers::invalid_unsigned_int) continue; if (shape_data.nodal_at_cell_boundaries && - eval.get_all_face_orientations()[v] != 0) - orientation[v] = &eval.get_shape_info().face_orientations_dofs( - eval.get_all_face_orientations()[v], 0); + fe_eval.get_all_face_orientations()[v] != 0) + orientation[v] = &fe_eval.get_shape_info().face_orientations_dofs( + fe_eval.get_all_face_orientations()[v], 0); } - else if (eval.get_face_orientation() != 0) - orientation[0] = &eval.get_shape_info().face_orientations_dofs( - eval.get_face_orientation(), 0); + else if (fe_eval.get_face_orientation() != 0) + orientation[0] = &fe_eval.get_shape_info().face_orientations_dofs( + fe_eval.get_face_orientation(), 0); // face_to_cell_index_hermite std::array index_array_hermite = @@ -3217,13 +3239,13 @@ namespace internal { if (n_face_orientations == 1) index_array_hermite[0] = - &eval.get_shape_info().face_to_cell_index_hermite(face_no, 0); + &fe_eval.get_shape_info().face_to_cell_index_hermite(face_no, 0); else { - const auto &face_nos = eval.get_all_face_numbers(); + const auto &face_nos = fe_eval.get_all_face_numbers(); for (unsigned int v = 0; v < n_lanes; ++v) { - if (eval.get_cell_ids()[v] == numbers::invalid_unsigned_int) + if (fe_eval.get_cell_ids()[v] == numbers::invalid_unsigned_int) continue; grad_weight[v] = @@ -3232,8 +3254,8 @@ namespace internal (1 + (face_nos[v] % 2)))][0]; index_array_hermite[v] = - &eval.get_shape_info().face_to_cell_index_hermite(face_nos[v], - 0); + &fe_eval.get_shape_info().face_to_cell_index_hermite( + face_nos[v], 0); } } } @@ -3245,23 +3267,23 @@ namespace internal { if (n_face_orientations == 1) index_array_nodal[0] = - &eval.get_shape_info().face_to_cell_index_nodal(face_no, 0); + &fe_eval.get_shape_info().face_to_cell_index_nodal(face_no, 0); else { - const auto &face_nos = eval.get_all_face_numbers(); + const auto &face_nos = fe_eval.get_all_face_numbers(); for (unsigned int v = 0; v < n_lanes; ++v) { - if (eval.get_cell_ids()[v] == numbers::invalid_unsigned_int) + if (fe_eval.get_cell_ids()[v] == numbers::invalid_unsigned_int) continue; index_array_nodal[v] = - &eval.get_shape_info().face_to_cell_index_nodal(face_nos[v], - 0); + &fe_eval.get_shape_info().face_to_cell_index_nodal( + face_nos[v], 0); } } } - const unsigned int subface_index = eval.get_subface_index(); + const unsigned int subface_index = fe_eval.get_subface_index(); const auto reorientate = [&](const unsigned int v, const unsigned int i) { return (dim < 3 || orientation[v] == nullptr || subface_index < Utilities::pow(2U, dim)) ? @@ -3275,9 +3297,9 @@ namespace internal for (unsigned int comp = 0; comp < n_components; ++comp) { const std::size_t index_offset = - dof_info - .component_dof_indices_offset[eval.get_active_fe_index()] - [eval.get_first_selected_component()] + + dof_info.component_dof_indices_offset + [fe_eval.get_active_fe_index()] + [fe_eval.get_first_selected_component()] + comp * Utilities::pow(fe_degree + 1, dim); // case 1: contiguous and interleaved indices @@ -3514,7 +3536,7 @@ namespace internal } else if (n_face_orientations == n_lanes) { - const auto &cells = eval.get_cell_ids(); + const auto &cells = fe_eval.get_cell_ids(); for (unsigned int v = 0; v < n_lanes; ++v) if (cells[v] != numbers::invalid_unsigned_int) { @@ -3548,7 +3570,7 @@ namespace internal for (unsigned int v = 0; v < n_lanes; ++v) reordered_indices[v] = dof_info.dof_indices_contiguous[dof_access_index] - [eval.get_cell_ids()[v]]; + [fe_eval.get_cell_ids()[v]]; dof_indices = reordered_indices.data(); } @@ -3671,108 +3693,109 @@ namespace internal const EvaluationFlags::EvaluationFlags evaluation_flag, const Number2 * src_ptr, const std::vector> * sm_ptr, - FEEvaluationData &eval) + FEEvaluationData &fe_eval) { Assert(fe_degree > -1, ExcInternalError()); - Assert(eval.get_shape_info().element_type <= + Assert(fe_eval.get_shape_info().element_type <= MatrixFreeFunctions::tensor_symmetric, ExcInternalError()); const unsigned int dofs_per_face = Utilities::pow(fe_degree + 1, dim - 1); - VectorizedArrayType *temp = eval.get_scratch_data().begin(); + VectorizedArrayType *temp = fe_eval.get_scratch_data().begin(); VectorizedArrayType *scratch_data = temp + 3 * n_components * dofs_per_face; Processor p; - if (eval.get_dof_access_index() == + if (fe_eval.get_dof_access_index() == MatrixFreeFunctions::DoFInfo::dof_access_cell && - eval.get_is_interior_face() == false) + fe_eval.get_is_interior_face() == false) fe_face_evaluation_process_and_io( - p, n_components, evaluation_flag, src_ptr, sm_ptr, eval, temp); + p, n_components, evaluation_flag, src_ptr, sm_ptr, fe_eval, temp); else fe_face_evaluation_process_and_io<1>( - p, n_components, evaluation_flag, src_ptr, sm_ptr, eval, temp); + p, n_components, evaluation_flag, src_ptr, sm_ptr, fe_eval, temp); - const unsigned int subface_index = eval.get_subface_index(); + const unsigned int subface_index = fe_eval.get_subface_index(); if (subface_index >= GeometryInfo::max_children_per_cell) - FEFaceEvaluationImpl< - true, - dim, - fe_degree, - n_q_points_1d, - VectorizedArrayType>::evaluate_in_face(n_components, - evaluation_flag, - eval.get_shape_info().data[0], - temp, - eval.begin_values(), - eval.begin_gradients(), - eval.begin_hessians(), - scratch_data, - subface_index); + FEFaceEvaluationImpl:: + evaluate_in_face(n_components, + evaluation_flag, + fe_eval.get_shape_info().data.front(), + temp, + fe_eval.begin_values(), + fe_eval.begin_gradients(), + fe_eval.begin_hessians(), + scratch_data, + subface_index); else - FEFaceEvaluationImpl< - false, - dim, - fe_degree, - n_q_points_1d, - VectorizedArrayType>::evaluate_in_face(n_components, - evaluation_flag, - eval.get_shape_info().data[0], - temp, - eval.begin_values(), - eval.begin_gradients(), - eval.begin_hessians(), - scratch_data, - subface_index); + FEFaceEvaluationImpl:: + evaluate_in_face(n_components, + evaluation_flag, + fe_eval.get_shape_info().data.front(), + temp, + fe_eval.begin_values(), + fe_eval.begin_gradients(), + fe_eval.begin_hessians(), + scratch_data, + subface_index); // re-orientation for cases not possible with above algorithm if (subface_index < GeometryInfo::max_children_per_cell) { - if (eval.get_dof_access_index() == + if (fe_eval.get_dof_access_index() == MatrixFreeFunctions::DoFInfo::dof_access_cell && - eval.get_is_interior_face() == false) + fe_eval.get_is_interior_face() == false) { for (unsigned int v = 0; v < VectorizedArrayType::size(); ++v) { // the loop breaks once an invalid_unsigned_int is hit for // all cases except the exterior faces in the ECL loop (where // some faces might be at the boundaries but others not) - if (eval.get_cell_ids()[v] == numbers::invalid_unsigned_int) + if (fe_eval.get_cell_ids()[v] == + numbers::invalid_unsigned_int) continue; - if (eval.get_all_face_orientations()[v] != 0) + if (fe_eval.get_all_face_orientations()[v] != 0) adjust_for_face_orientation_per_lane( dim, n_components, v, evaluation_flag, - &eval.get_shape_info().face_orientations_quad( - eval.get_all_face_orientations()[v], 0), + &fe_eval.get_shape_info().face_orientations_quad( + fe_eval.get_all_face_orientations()[v], 0), false, Utilities::pow(n_q_points_1d, dim - 1), &temp[0][0], - eval.begin_values(), - eval.begin_gradients(), - eval.begin_hessians()); + fe_eval.begin_values(), + fe_eval.begin_gradients(), + fe_eval.begin_hessians()); } } - else if (eval.get_face_orientation() != 0) + else if (fe_eval.get_face_orientation() != 0) for (unsigned int c = 0; c < n_components; ++c) adjust_for_face_orientation( dim, n_components, evaluation_flag, - &eval.get_shape_info().face_orientations_quad( - eval.get_face_orientation(), 0), + &fe_eval.get_shape_info().face_orientations_quad( + fe_eval.get_face_orientation(), 0), false, Utilities::pow(n_q_points_1d, dim - 1), temp, - eval.begin_values(), - eval.begin_gradients(), - eval.begin_hessians()); + fe_eval.begin_values(), + fe_eval.begin_gradients(), + fe_eval.begin_hessians()); } return false; @@ -3785,15 +3808,15 @@ namespace internal const Number * vector_ptr, MatrixFreeFunctions::DoFInfo::IndexStorageVariants storage) { - const unsigned int fe_degree = shape_info.data[0].fe_degree; - if (fe_degree < 1 || !shape_info.data[0].nodal_at_cell_boundaries || + const unsigned int fe_degree = shape_info.data.front().fe_degree; + if (fe_degree < 1 || !shape_info.data.front().nodal_at_cell_boundaries || (evaluation_flag & EvaluationFlags::gradients && (fe_degree < 2 || - shape_info.data[0].element_type != + shape_info.data.front().element_type != MatrixFreeFunctions::tensor_symmetric_hermite)) || (evaluation_flag & EvaluationFlags::hessians) || vector_ptr == nullptr || - shape_info.data[0].element_type > + shape_info.data.front().element_type > MatrixFreeFunctions::tensor_symmetric || storage < MatrixFreeFunctions::DoFInfo::IndexStorageVariants::contiguous) @@ -3892,107 +3915,107 @@ namespace internal const EvaluationFlags::EvaluationFlags integration_flag, Number2 * dst_ptr, const std::vector> * sm_ptr, - FEEvaluationData &eval) + FEEvaluationData &fe_eval) { Assert(fe_degree > -1, ExcInternalError()); - Assert(eval.get_shape_info().element_type <= + Assert(fe_eval.get_shape_info().element_type <= MatrixFreeFunctions::tensor_symmetric, ExcInternalError()); const unsigned int dofs_per_face = Utilities::pow(fe_degree + 1, dim - 1); - VectorizedArrayType *temp = eval.get_scratch_data().begin(); + VectorizedArrayType *temp = fe_eval.get_scratch_data().begin(); VectorizedArrayType *scratch_data = temp + 3 * n_components * dofs_per_face; - const unsigned int subface_index = eval.get_subface_index(); + const unsigned int subface_index = fe_eval.get_subface_index(); // re-orientation for cases not possible with the io function below if (subface_index < GeometryInfo::max_children_per_cell) { - if (eval.get_dof_access_index() == + if (fe_eval.get_dof_access_index() == MatrixFreeFunctions::DoFInfo::dof_access_cell && - eval.get_is_interior_face() == false) + fe_eval.get_is_interior_face() == false) for (unsigned int v = 0; v < VectorizedArrayType::size(); ++v) { // the loop breaks once an invalid_unsigned_int is hit for // all cases except the exterior faces in the ECL loop (where // some faces might be at the boundaries but others not) - if (eval.get_cell_ids()[v] == numbers::invalid_unsigned_int) + if (fe_eval.get_cell_ids()[v] == numbers::invalid_unsigned_int) continue; - if (eval.get_all_face_orientations()[v] != 0) + if (fe_eval.get_all_face_orientations()[v] != 0) adjust_for_face_orientation_per_lane( dim, n_components, v, integration_flag, - &eval.get_shape_info().face_orientations_quad( - eval.get_all_face_orientations()[v], 0), + &fe_eval.get_shape_info().face_orientations_quad( + fe_eval.get_all_face_orientations()[v], 0), true, Utilities::pow(n_q_points_1d, dim - 1), &temp[0][0], - eval.begin_values(), - eval.begin_gradients(), - eval.begin_hessians()); + fe_eval.begin_values(), + fe_eval.begin_gradients(), + fe_eval.begin_hessians()); } - else if (eval.get_face_orientation() != 0) + else if (fe_eval.get_face_orientation() != 0) adjust_for_face_orientation( dim, n_components, integration_flag, - &eval.get_shape_info().face_orientations_quad( - eval.get_face_orientation(), 0), + &fe_eval.get_shape_info().face_orientations_quad( + fe_eval.get_face_orientation(), 0), true, Utilities::pow(n_q_points_1d, dim - 1), temp, - eval.begin_values(), - eval.begin_gradients(), - eval.begin_hessians()); + fe_eval.begin_values(), + fe_eval.begin_gradients(), + fe_eval.begin_hessians()); } - if (fe_degree > -1 && eval.get_subface_index() >= + if (fe_degree > -1 && fe_eval.get_subface_index() >= GeometryInfo::max_children_per_cell) - FEFaceEvaluationImpl< - true, - dim, - fe_degree, - n_q_points_1d, - VectorizedArrayType>::integrate_in_face(n_components, - integration_flag, - eval.get_shape_info().data[0], - temp, - eval.begin_values(), - eval.begin_gradients(), - eval.begin_hessians(), - scratch_data, - subface_index); + FEFaceEvaluationImpl:: + integrate_in_face(n_components, + integration_flag, + fe_eval.get_shape_info().data.front(), + temp, + fe_eval.begin_values(), + fe_eval.begin_gradients(), + fe_eval.begin_hessians(), + scratch_data, + subface_index); else - FEFaceEvaluationImpl< - false, - dim, - fe_degree, - n_q_points_1d, - VectorizedArrayType>::integrate_in_face(n_components, - integration_flag, - eval.get_shape_info().data[0], - temp, - eval.begin_values(), - eval.begin_gradients(), - eval.begin_hessians(), - scratch_data, - subface_index); + FEFaceEvaluationImpl:: + integrate_in_face(n_components, + integration_flag, + fe_eval.get_shape_info().data.front(), + temp, + fe_eval.begin_values(), + fe_eval.begin_gradients(), + fe_eval.begin_hessians(), + scratch_data, + subface_index); Processor p; - if (eval.get_dof_access_index() == + if (fe_eval.get_dof_access_index() == MatrixFreeFunctions::DoFInfo::dof_access_cell && - eval.get_is_interior_face() == false) + fe_eval.get_is_interior_face() == false) fe_face_evaluation_process_and_io( - p, n_components, integration_flag, dst_ptr, sm_ptr, eval, temp); + p, n_components, integration_flag, dst_ptr, sm_ptr, fe_eval, temp); else fe_face_evaluation_process_and_io<1>( - p, n_components, integration_flag, dst_ptr, sm_ptr, eval, temp); + p, n_components, integration_flag, dst_ptr, sm_ptr, fe_eval, temp); return false; } diff --git a/include/deal.II/matrix_free/evaluation_template_factory.h b/include/deal.II/matrix_free/evaluation_template_factory.h index d372e330f9..9eeb83700c 100644 --- a/include/deal.II/matrix_free/evaluation_template_factory.h +++ b/include/deal.II/matrix_free/evaluation_template_factory.h @@ -42,13 +42,13 @@ namespace internal evaluate(const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, const Number * values_dofs, - FEEvaluationData & eval); + FEEvaluationData & fe_eval); static void integrate(const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, Number * values_dofs, - FEEvaluationData & eval, + FEEvaluationData & fe_eval, const bool sum_into_values_array); static bool @@ -65,13 +65,13 @@ namespace internal evaluate(const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, const Number * values_dofs, - FEEvaluationData & eval); + FEEvaluationData & fe_eval); static void integrate(const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, Number * values_dofs, - FEEvaluationData & eval); + FEEvaluationData & fe_eval); }; @@ -84,14 +84,14 @@ namespace internal const EvaluationFlags::EvaluationFlags evaluation_flag, const Number * src_ptr, const std::vector> * sm_ptr, - FEEvaluationData &eval); + FEEvaluationData &fe_eval); static void integrate(const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, Number * dst_ptr, const std::vector> *sm_ptr, - FEEvaluationData &eval); + FEEvaluationData &fe_eval); }; @@ -101,7 +101,7 @@ namespace internal { static void apply(const unsigned int n_components, - const FEEvaluationData &fe_eval, + const FEEvaluationData &fe_fe_eval, const Number * in_array, Number * out_array); @@ -116,7 +116,7 @@ namespace internal static void transform_from_q_points_to_basis( const unsigned int n_components, - const FEEvaluationData &fe_eval, + const FEEvaluationData &fe_fe_eval, const Number * in_array, Number * out_array); }; diff --git a/include/deal.II/matrix_free/evaluation_template_factory.templates.h b/include/deal.II/matrix_free/evaluation_template_factory.templates.h index a644970b0a..cdaca44458 100644 --- a/include/deal.II/matrix_free/evaluation_template_factory.templates.h +++ b/include/deal.II/matrix_free/evaluation_template_factory.templates.h @@ -47,15 +47,15 @@ namespace internal const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, const Number * values_dofs, - FEEvaluationData & eval) + FEEvaluationData & fe_eval) { instantiation_helper_run<1, FEEvaluationImplEvaluateSelector>( - eval.get_shape_info().data[0].fe_degree, - eval.get_shape_info().data[0].n_q_points_1d, + fe_eval.get_shape_info().data[0].fe_degree, + fe_eval.get_shape_info().data[0].n_q_points_1d, n_components, evaluation_flag, values_dofs, - eval); + fe_eval); } @@ -66,16 +66,16 @@ namespace internal const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, Number * values_dofs, - FEEvaluationData & eval, + FEEvaluationData & fe_eval, const bool sum_into_values_array) { instantiation_helper_run<1, FEEvaluationImplIntegrateSelector>( - eval.get_shape_info().data[0].fe_degree, - eval.get_shape_info().data[0].n_q_points_1d, + fe_eval.get_shape_info().data[0].fe_degree, + fe_eval.get_shape_info().data[0].n_q_points_1d, n_components, integration_flag, values_dofs, - eval, + fe_eval, sum_into_values_array); } @@ -99,16 +99,16 @@ namespace internal const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, const Number * values_dofs, - FEEvaluationData & eval) + FEEvaluationData & fe_eval) { instantiation_helper_run<1, FEFaceEvaluationImplEvaluateSelector>( - eval.get_shape_info().data[0].fe_degree, - eval.get_shape_info().data[0].n_q_points_1d, + fe_eval.get_shape_info().data[0].fe_degree, + fe_eval.get_shape_info().data[0].n_q_points_1d, n_components, evaluation_flag, values_dofs, - eval); + fe_eval); } @@ -119,17 +119,17 @@ namespace internal const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, Number * values_dofs, - FEEvaluationData & eval) + FEEvaluationData & fe_eval) { instantiation_helper_run< 1, FEFaceEvaluationImplIntegrateSelector>( - eval.get_shape_info().data[0].fe_degree, - eval.get_shape_info().data[0].n_q_points_1d, + fe_eval.get_shape_info().data[0].fe_degree, + fe_eval.get_shape_info().data[0].n_q_points_1d, n_components, integration_flag, values_dofs, - eval); + fe_eval); } @@ -141,20 +141,20 @@ namespace internal const EvaluationFlags::EvaluationFlags evaluation_flag, const Number * src_ptr, const std::vector> * sm_ptr, - FEEvaluationData &eval) + FEEvaluationData &fe_eval) { instantiation_helper_run< 1, FEFaceEvaluationImplGatherEvaluateSelector>( - eval.get_shape_info().data[0].fe_degree, - eval.get_shape_info().data[0].n_q_points_1d, + fe_eval.get_shape_info().data[0].fe_degree, + fe_eval.get_shape_info().data[0].n_q_points_1d, n_components, evaluation_flag, src_ptr, sm_ptr, - eval); + fe_eval); } @@ -166,20 +166,20 @@ namespace internal const EvaluationFlags::EvaluationFlags integration_flag, Number * dst_ptr, const std::vector> * sm_ptr, - FEEvaluationData &eval) + FEEvaluationData &fe_eval) { instantiation_helper_run< 1, FEFaceEvaluationImplIntegrateScatterSelector>( - eval.get_shape_info().data[0].fe_degree, - eval.get_shape_info().data[0].n_q_points_1d, + fe_eval.get_shape_info().data[0].fe_degree, + fe_eval.get_shape_info().data[0].n_q_points_1d, n_components, integration_flag, dst_ptr, sm_ptr, - eval); + fe_eval); } diff --git a/include/deal.II/matrix_free/fe_evaluation.h b/include/deal.II/matrix_free/fe_evaluation.h index a21bf2b821..4840493292 100644 --- a/include/deal.II/matrix_free/fe_evaluation.h +++ b/include/deal.II/matrix_free/fe_evaluation.h @@ -771,7 +771,7 @@ protected: /** * A pointer to the underlying data. */ - const MatrixFree *matrix_info; + const MatrixFree *matrix_free; /** * A temporary data structure necessary to read degrees of freedom when no @@ -1386,25 +1386,25 @@ protected: * to quadrature points: * * @code - * FEEvaluation phi(matrix_free); + * FEEvaluation fe_eval(matrix_free); * for (unsigned int cell_index = cell_range.first; * cell_index < cell_range.second; ++cell_index) * { - * phi.reinit(cell_index); - * phi.read_dof_values(vector); - * phi.evaluate(EvaluationFlags::values); // interpolate values only - * for (unsigned int q=0; q val = phi.get_value(q); + * VectorizedArray val = fe_eval.get_value(q); * // do something with val * } * } * @endcode * - * Likewise, a gradient of the finite element solution represented by @p - * vector can be interpolated to the quadrature points by @p - * phi.get_gradient(q). The combination of read_dof_values(), evaluate() and - * get_value() is similar to what FEValues::get_function_values or + * Likewise, a gradient of the finite element solution represented by + * `vector` can be interpolated to the quadrature points by + * `fe_eval.get_gradient(q)`. The combination of read_dof_values(), evaluate() + * and get_value() is similar to what FEValues::get_function_values or * FEValues::get_function_gradients does, but it is in general much faster * because it makes use of the tensor product, see the description of the * evaluation routines below, and can do this operation for several cells at @@ -1421,16 +1421,16 @@ protected: * the following cell-wise integration: * * @code - * FEEvaluation phi(matrix_free); + * FEEvaluation fe_eval(matrix_free); * Function &function = ...; * for (unsigned int cell_index = cell_range.first; * cell_index < cell_range.second; ++cell_index) * { - * phi.reinit(cell_index); - * for (unsigned int q=0; q > p_vect = - * phi.quadrature_point(q); + * fe_eval.quadrature_point(q); * // Need to evaluate function for each component in VectorizedArray * VectorizedArray f_value = 0.0; * for (unsigned int v=0; v::size(); ++v) @@ -1440,14 +1440,14 @@ protected: * p[d] = p_vect[d][v]; * f_value[v] = function.value(p); * } - * phi.submit_value(f_value, q); + * fe_eval.submit_value(f_value, q); * } - * phi.integrate(EvaluationFlags::values); - * phi.distribute_local_to_global(dst); + * fe_eval.integrate(EvaluationFlags::values); + * fe_eval.distribute_local_to_global(dst); * } * @endcode * - * In this code, the call to @p phi.submit_value() prepares for the + * In this code, the call to @p fe_eval.submit_value() prepares for the * multiplication by the test function prior to the actual integration (inside * the submit call, the value to be tested is also multiplied by the * determinant of the Jacobian and the quadrature weight). In the @@ -1743,14 +1743,14 @@ protected: * template parameter on the number of components: * * @code - * FEEvaluation phi(matrix_free); + * FEEvaluation fe_eval(matrix_free); * @endcode * * If used this way, the components can be gathered from several components of * an @p std::vector through the call * * @code - * phi.read_dof_values(src, 0); + * fe_eval.read_dof_values(src, 0); * @endcode * * where the 0 means that the vectors starting from the zeroth vector in the @@ -1852,13 +1852,13 @@ protected: * points: * * @code - * phi1.evaluate(EvaluationFlags::values); - * phi2.evaluate(EvaluationFlags::gradients); - * for (unsigned int q=0; q val1 = phi1.get_value(q); - * Tensor<1,dim,VectorizedArray > grad2 = phi2.get_gradient(q); - * Point > point = phi1.quadrature_point(q); + * VectorizedArray val1 = fe_eval1.get_value(q); + * Tensor<1,dim,VectorizedArray > grad2 = fe_eval2.get_gradient(q); + * Point > point = fe_eval1.quadrature_point(q); * // ... some complicated formula combining those three... * } * @endcode @@ -2319,7 +2319,7 @@ public: */ template void - integrate_scatter(const EvaluationFlags::EvaluationFlags evaluation_flag, + integrate_scatter(const EvaluationFlags::EvaluationFlags integration_flag, VectorType & output_vector); /** @@ -2683,7 +2683,7 @@ public: * distribute_local_to_global() or set_dof_values() methods). */ void - integrate(const EvaluationFlags::EvaluationFlags evaluation_flag); + integrate(const EvaluationFlags::EvaluationFlags integration_flag); /** * @deprecated Please use the integrate() function with the EvaluationFlags argument. @@ -2700,7 +2700,7 @@ public: * call stores the result of the testing in the given array `values_array`. */ void - integrate(const EvaluationFlags::EvaluationFlags evaluation_flag, + integrate(const EvaluationFlags::EvaluationFlags integration_flag, VectorizedArrayType * values_array); /** @@ -2724,7 +2724,7 @@ public: */ template void - integrate_scatter(const EvaluationFlags::EvaluationFlags evaluation_flag, + integrate_scatter(const EvaluationFlags::EvaluationFlags integration_flag, VectorType & output_vector); /** @@ -2812,7 +2812,7 @@ namespace internal inline typename FEEvaluationData:: InitializationData extract_initialization_data( - const MatrixFree &data, + const MatrixFree &matrix_free, const unsigned int dof_no, const unsigned int first_selected_component, const unsigned int quad_no, @@ -2825,11 +2825,11 @@ namespace internal typename FEEvaluationData:: InitializationData init_data; - init_data.dof_info = &data.get_dof_info(dof_no); + init_data.dof_info = &matrix_free.get_dof_info(dof_no); init_data.mapping_data = &internal::MatrixFreeFunctions:: MappingInfoCellsOrFaces::get( - data.get_mapping_info(), quad_no); + matrix_free.get_mapping_info(), quad_no); init_data.active_fe_index = fe_degree != numbers::invalid_unsigned_int ? @@ -2847,7 +2847,7 @@ namespace internal 1)) : init_data.mapping_data->quad_index_from_n_q_points(n_q_points); - init_data.shape_info = &data.get_shape_info( + init_data.shape_info = &matrix_free.get_shape_info( dof_no, quad_no, init_data.dof_info->component_to_base_index[first_selected_component], @@ -2902,7 +2902,7 @@ inline FEEvaluationBaseset_data_pointers(scratch_data_array, n_components_); Assert( @@ -2965,7 +2965,7 @@ inline FEEvaluationBase()) - , matrix_info(nullptr) + , matrix_free(nullptr) { const unsigned int base_element_number = fe.component_to_base_index(first_selected_component).first; @@ -3005,12 +3005,12 @@ inline FEEvaluationBase &other) : FEEvaluationData(other) - , scratch_data_array(other.matrix_info == nullptr ? + , scratch_data_array(other.matrix_free == nullptr ? new AlignedVector() : - other.matrix_info->acquire_scratch_data()) - , matrix_info(other.matrix_info) + other.matrix_free->acquire_scratch_data()) + , matrix_free(other.matrix_free) { - if (other.matrix_info == nullptr) + if (other.matrix_free == nullptr) { Assert(other.mapped_geometry.get() != nullptr, ExcInternalError()); this->data = @@ -3056,21 +3056,21 @@ operator=(const FEEvaluationBase &other) { // release old memory - if (matrix_info == nullptr) + if (matrix_free == nullptr) { delete this->data; delete scratch_data_array; } else { - matrix_info->release_scratch_data(scratch_data_array); + matrix_free->release_scratch_data(scratch_data_array); } this->FEEvaluationData::operator=(other); - matrix_info = other.matrix_info; + matrix_free = other.matrix_free; - if (other.matrix_info == nullptr) + if (other.matrix_free == nullptr) { Assert(other.mapped_geometry.get() != nullptr, ExcInternalError()); this->data = @@ -3094,7 +3094,7 @@ operator=(const FEEvaluationBaseacquire_scratch_data(); + scratch_data_array = matrix_free->acquire_scratch_data(); } this->set_data_pointers(scratch_data_array, n_components_); @@ -3115,11 +3115,11 @@ inline FEEvaluationBase::~FEEvaluationBase() { - if (matrix_info != nullptr) + if (matrix_free != nullptr) { try { - matrix_info->release_scratch_data(scratch_data_array); + matrix_free->release_scratch_data(scratch_data_array); } catch (...) {} @@ -3142,10 +3142,10 @@ inline const MatrixFree & FEEvaluationBase:: get_matrix_free() const { - Assert(matrix_info != nullptr, + Assert(matrix_free != nullptr, ExcMessage( "FEEvaluation was not initialized with a MatrixFree object!")); - return *matrix_info; + return *matrix_free; } @@ -3165,19 +3165,19 @@ namespace internal n_components, Number, is_face, - VectorizedArrayType> & phi, - const MatrixFree *matrix_info, + VectorizedArrayType> & fe_eval, + const MatrixFree *matrix_free, GlobalVectorType & array, VectorizedArrayType2 & out, const FU & fu) { - (void)matrix_info; - Assert(matrix_info != nullptr, ExcNotImplemented()); + (void)matrix_free; + Assert(matrix_free != nullptr, ExcNotImplemented()); AssertDimension(array.size(), - matrix_info->get_task_info().cell_partition_data.back()); + matrix_free->get_task_info().cell_partition_data.back()); // 1) collect ids of cell - const auto cells = phi.get_cell_ids(); + const auto cells = fe_eval.get_cell_ids(); // 2) actually gather values for (unsigned int i = 0; i < VectorizedArrayType::size(); ++i) @@ -3201,7 +3201,7 @@ FEEvaluationBase:: { VectorizedArrayType out = Number(1.); internal::process_cell_data( - *this, this->matrix_info, array, out, [](auto &local, const auto &global) { + *this, this->matrix_free, array, out, [](auto &local, const auto &global) { local = global; }); return out; @@ -3220,7 +3220,7 @@ FEEvaluationBase:: const VectorizedArrayType & in) const { internal::process_cell_data( - *this, this->matrix_info, array, in, [](const auto &local, auto &global) { + *this, this->matrix_free, array, in, [](const auto &local, auto &global) { global = local; }); } @@ -3240,7 +3240,7 @@ FEEvaluationBase:: { std::array out; internal::process_cell_data( - *this, this->matrix_info, array, out, [](auto &local, const auto &global) { + *this, this->matrix_free, array, out, [](auto &local, const auto &global) { local = global; }); return out; @@ -3261,7 +3261,7 @@ FEEvaluationBase:: const std::array & in) const { internal::process_cell_data( - *this, this->matrix_info, array, in, [](const auto &local, auto &global) { + *this, this->matrix_free, array, in, [](const auto &local, auto &global) { global = local; }); } @@ -3366,14 +3366,14 @@ FEEvaluationBase:: // Case 1: No MatrixFree object given, simple case because we do not need to // process constraints and need not care about vectorization -> go to // separate function - if (this->matrix_info == nullptr) + if (this->matrix_free == nullptr) { read_write_operation_global(operation, src); return; } Assert(this->dof_info != nullptr, ExcNotInitialized()); - Assert(this->matrix_info->indices_initialized() == true, ExcNotInitialized()); + Assert(this->matrix_free->indices_initialized() == true, ExcNotInitialized()); if (this->n_fe_components == 1) for (unsigned int comp = 0; comp < n_components; ++comp) { @@ -3495,8 +3495,8 @@ FEEvaluationBase:: internal::MatrixFreeFunctions::DoFInfo::dof_access_cell ? &cells_copied[0] : (this->is_interior_face ? - &this->matrix_info->get_face_info(this->cell).cells_interior[0] : - &this->matrix_info->get_face_info(this->cell).cells_exterior[0]); + &this->matrix_free->get_face_info(this->cell).cells_interior[0] : + &this->matrix_free->get_face_info(this->cell).cells_exterior[0]); for (unsigned int v = 0; v < n_vectorization_actual; ++v) { Assert(cells[v] < this->dof_info->row_starts.size() - 1, @@ -3651,9 +3651,9 @@ FEEvaluationBase:: value[comp]); const Number *data_val = - this->matrix_info->constraint_pool_begin(indicator.second); + this->matrix_free->constraint_pool_begin(indicator.second); const Number *end_pool = - this->matrix_info->constraint_pool_end(indicator.second); + this->matrix_free->constraint_pool_end(indicator.second); for (; data_val != end_pool; ++data_val, ++dof_indices[v]) for (unsigned int comp = 0; comp < n_components; ++comp) operation.process_constraint(*dof_indices[v], @@ -3708,9 +3708,9 @@ FEEvaluationBase:: value); const Number *data_val = - this->matrix_info->constraint_pool_begin(indicator.second); + this->matrix_free->constraint_pool_begin(indicator.second); const Number *end_pool = - this->matrix_info->constraint_pool_end(indicator.second); + this->matrix_free->constraint_pool_end(indicator.second); for (; data_val != end_pool; ++data_val, ++dof_indices[v]) operation.process_constraint(*dof_indices[v], @@ -4274,8 +4274,8 @@ FEEvaluationBase:: internal::MatrixFreeFunctions::DoFInfo::dof_access_cell ? &cells_copied[0] : (this->is_interior_face ? - &this->matrix_info->get_face_info(this->cell).cells_interior[0] : - &this->matrix_info->get_face_info(this->cell).cells_exterior[0]); + &this->matrix_free->get_face_info(this->cell).cells_interior[0] : + &this->matrix_free->get_face_info(this->cell).cells_exterior[0]); } for (unsigned int v = 0; v < n_vectorization_actual; ++v) @@ -6715,12 +6715,12 @@ FEEvaluationmatrix_info->n_components(); + for (unsigned int no = 0; no < this->matrix_free->n_components(); ++no) for (unsigned int nf = 0; - nf < this->matrix_info->n_base_elements(no); + nf < this->matrix_free->n_base_elements(no); ++nf) - if (this->matrix_info + if (this->matrix_free ->get_shape_info(no, 0, nf, this->active_fe_index, 0) .data.front() .fe_degree == static_cast(fe_degree)) @@ -6735,9 +6735,9 @@ FEEvaluationquad_no; else for (unsigned int no = 0; - no < this->matrix_info->get_mapping_info().cell_data.size(); + no < this->matrix_free->get_mapping_info().cell_data.size(); ++no) - if (this->matrix_info->get_mapping_info() + if (this->matrix_free->get_mapping_info() .cell_data[no] .descriptor[this->active_quad_index] .n_q_points == n_q_points) @@ -6850,7 +6850,7 @@ FEEvaluationmapping_data != nullptr, ExcNotInitialized()); this->cell = cell_index; this->cell_type = - this->matrix_info->get_mapping_info().get_cell_type(cell_index); + this->matrix_free->get_mapping_info().get_cell_type(cell_index); const unsigned int offsets = this->mapping_data->data_index_offsets[cell_index]; @@ -6886,7 +6886,7 @@ FEEvaluation:: reinit(const TriaIterator> &cell) { - Assert(this->matrix_info == nullptr, + Assert(this->matrix_free == nullptr, ExcMessage("Cannot use initialization from cell iterator if " "initialized from MatrixFree object. Use variant for " "on the fly computation with arguments as for FEValues " @@ -6918,7 +6918,7 @@ FEEvaluation:: reinit(const typename Triangulation::cell_iterator &cell) { - Assert(this->matrix_info == 0, + Assert(this->matrix_free == 0, ExcMessage("Cannot use initialization from cell iterator if " "initialized from MatrixFree object. Use variant for " "on the fly computation with arguments as for FEValues " @@ -6943,7 +6943,7 @@ FEEvaluation::quadrature_point(const unsigned int q) const { - if (this->matrix_info == nullptr) + if (this->matrix_free == nullptr) { Assert((this->mapped_geometry->get_fe_values().get_update_flags() | update_quadrature_points), @@ -7094,14 +7094,18 @@ FEEvaluation -1) - SelectEvaluator:: - evaluate(n_components, evaluation_flag_actual, values_array, *this); + { + SelectEvaluator:: + evaluate(n_components, evaluation_flag_actual, values_array, *this); + } else - internal::FEEvaluationFactory::evaluate( - n_components, - evaluation_flag_actual, - const_cast(values_array), - *this); + { + internal::FEEvaluationFactory::evaluate( + n_components, + evaluation_flag_actual, + const_cast(values_array), + *this); + } # ifdef DEBUG if (evaluation_flag_actual & EvaluationFlags::values) @@ -7160,10 +7164,10 @@ namespace internal Number *>::value, VectorType>::type * = nullptr> VectorizedArrayType * - check_vector_access_inplace(const EvaluatorType &phi, VectorType &vector) + check_vector_access_inplace(const EvaluatorType &fe_eval, VectorType &vector) { - const unsigned int cell = phi.get_cell_or_face_batch_id(); - const auto & dof_info = phi.get_dof_info(); + const unsigned int cell = fe_eval.get_cell_or_face_batch_id(); + const auto & dof_info = fe_eval.get_dof_info(); // If the index storage is interleaved and contiguous and the vector // storage has the correct alignment, we can directly pass the pointer @@ -7189,7 +7193,8 @@ namespace internal [internal::MatrixFreeFunctions::DoFInfo::dof_access_cell] [cell * VectorizedArrayType::size()] + dof_info.component_dof_indices_offset - [phi.get_active_fe_index()][phi.get_first_selected_component()] * + [fe_eval.get_active_fe_index()] + [fe_eval.get_first_selected_component()] * VectorizedArrayType::size()); } else @@ -7349,7 +7354,7 @@ FEEvaluationhessians_quad_submitted == true, internal::ExcAccessToUninitializedField()); # endif - Assert(this->matrix_info != nullptr || + Assert(this->matrix_free != nullptr || this->mapped_geometry->is_initialized(), ExcNotInitialized()); @@ -7378,19 +7383,23 @@ FEEvaluation -1) - SelectEvaluator:: - integrate(n_components, - integration_flag_actual, - values_array, - *this, - sum_into_values_array); + { + SelectEvaluator:: + integrate(n_components, + integration_flag_actual, + values_array, + *this, + sum_into_values_array); + } else - internal::FEEvaluationFactory::integrate( - n_components, - integration_flag_actual, - values_array, - *this, - sum_into_values_array); + { + internal::FEEvaluationFactory::integrate( + n_components, + integration_flag_actual, + values_array, + *this, + sum_into_values_array); + } # ifdef DEBUG this->dof_values_initialized = true; @@ -7560,20 +7569,20 @@ FEFaceEvaluationmapping_data != nullptr, ExcNotInitialized()); if (face_index >= - this->matrix_info->get_task_info().face_partition_data.back() && + this->matrix_free->get_task_info().face_partition_data.back() && face_index < - this->matrix_info->get_task_info().boundary_partition_data.back()) + this->matrix_free->get_task_info().boundary_partition_data.back()) Assert(this->is_interior_face, ExcMessage( "Boundary faces do not have a neighbor. When looping over " "boundary faces use FEFaceEvaluation with the parameter " "is_interior_face set to true. ")); - this->reinit_face(this->matrix_info->get_face_info(face_index)); + this->reinit_face(this->matrix_free->get_face_info(face_index)); for (unsigned int i = 0; i < VectorizedArrayType::size(); ++i) this->cell_or_face_ids[i] = face_index * VectorizedArrayType::size() + i; - this->cell_type = this->matrix_info->get_mapping_info().face_type[face_index]; + this->cell_type = this->matrix_free->get_mapping_info().face_type[face_index]; const unsigned int offsets = this->mapping_data->data_index_offsets[face_index]; this->J_value = &this->mapping_data->JxW_values[offsets]; @@ -7611,20 +7620,20 @@ FEFaceEvaluationquad_no < - this->matrix_info->get_mapping_info().face_data_by_cells.size(), + this->matrix_free->get_mapping_info().face_data_by_cells.size(), ExcMessage( "You must set MatrixFree::AdditionalData::mapping_update_flags_faces_by_cells to use the present reinit method.")); AssertIndexRange(face_number, GeometryInfo::faces_per_cell); AssertIndexRange(cell_index, - this->matrix_info->get_mapping_info().cell_type.size()); + this->matrix_free->get_mapping_info().cell_type.size()); Assert(this->mapped_geometry == nullptr, ExcMessage("FEEvaluation was initialized without a matrix-free object." " Integer indexing is not possible")); if (this->mapped_geometry != nullptr) return; - Assert(this->matrix_info != nullptr, ExcNotInitialized()); + Assert(this->matrix_free != nullptr, ExcNotInitialized()); - this->cell_type = this->matrix_info->get_mapping_info().cell_type[cell_index]; + this->cell_type = this->matrix_free->get_mapping_info().cell_type[cell_index]; this->cell = cell_index; this->subface_index = GeometryInfo::max_children_per_cell; this->dof_access_index = @@ -7643,7 +7652,7 @@ FEFaceEvaluationmatrix_info->get_cell_and_face_to_plain_faces()(cell_index, + this->matrix_free->get_cell_and_face_to_plain_faces()(cell_index, face_number, i); @@ -7656,18 +7665,18 @@ FEFaceEvaluationmatrix_info->get_face_info(face_index / n_lanes); + this->matrix_free->get_face_info(face_index / n_lanes); // get cell ID on both sides of face auto cell_m = faces.cells_interior[face_index % n_lanes]; auto cell_p = faces.cells_exterior[face_index % n_lanes]; - const bool is_interior_face = cell_m != cell_this; + const bool face_identifies_as_interior = cell_m != cell_this; Assert(cell_m == cell_this || cell_p == cell_this, ExcInternalError()); // compare the IDs with the given cell ID - if (is_interior_face) + if (face_identifies_as_interior) { this->cell_ids[i] = cell_m; // neighbor has the other ID this->all_face_numbers[i] = faces.interior_face_no; @@ -7680,7 +7689,7 @@ FEFaceEvaluation= 8; unsigned int face_orientation = faces.face_orientation % 8; - if (is_interior_face != orientation_interior_face) + if (face_identifies_as_interior != orientation_interior_face) { constexpr std::array table{ {0, 1, 2, 3, 6, 5, 4, 7}}; @@ -7701,25 +7710,25 @@ FEFaceEvaluationmatrix_info->get_mapping_info() + this->matrix_free->get_mapping_info() .face_data_by_cells[this->quad_no] .data_index_offsets[cell_index * GeometryInfo::faces_per_cell + face_number]; AssertIndexRange(offsets, - this->matrix_info->get_mapping_info() + this->matrix_free->get_mapping_info() .face_data_by_cells[this->quad_no] .JxW_values.size()); - this->J_value = &this->matrix_info->get_mapping_info() + this->J_value = &this->matrix_free->get_mapping_info() .face_data_by_cells[this->quad_no] .JxW_values[offsets]; - this->normal_vectors = &this->matrix_info->get_mapping_info() + this->normal_vectors = &this->matrix_free->get_mapping_info() .face_data_by_cells[this->quad_no] .normal_vectors[offsets]; - this->jacobian = &this->matrix_info->get_mapping_info() + this->jacobian = &this->matrix_free->get_mapping_info() .face_data_by_cells[this->quad_no] .jacobians[!this->is_interior_face][offsets]; this->normal_x_jacobian = - &this->matrix_info->get_mapping_info() + &this->matrix_free->get_mapping_info() .face_data_by_cells[this->quad_no] .normals_times_jacobians[!this->is_interior_face][offsets]; @@ -7872,9 +7881,9 @@ FEFaceEvaluation:: - integrate(const EvaluationFlags::EvaluationFlags evaluation_flag) + integrate(const EvaluationFlags::EvaluationFlags integration_flag) { - integrate(evaluation_flag, this->values_dofs); + integrate(integration_flag, this->values_dofs); # ifdef DEBUG this->dof_values_initialized = true; @@ -8032,11 +8041,11 @@ FEFaceEvaluationdof_access_index][this->cell])) { if (fe_degree > -1) - internal::FEFaceEvaluationImplGatherEvaluateSelector< - dim, - Number, - VectorizedArrayType>::template run( - n_components, - evaluation_flag, - internal::get_beginning(input_vector), - shared_vector_data, - *this); + { + internal::FEFaceEvaluationImplGatherEvaluateSelector< + dim, + Number, + VectorizedArrayType>::template run( + n_components, + evaluation_flag, + internal::get_beginning(input_vector), + shared_vector_data, + *this); + } else - internal::FEFaceEvaluationGatherFactory< - dim, - Number, - VectorizedArrayType>::evaluate(n_components, - evaluation_flag, - internal::get_beginning( - input_vector), - shared_vector_data, - *this); + { + internal::FEFaceEvaluationGatherFactory< + dim, + Number, + VectorizedArrayType>::evaluate(n_components, + evaluation_flag, + internal::get_beginning( + input_vector), + shared_vector_data, + *this); + } } else { @@ -8169,26 +8182,30 @@ FEFaceEvaluationdof_access_index][this->cell])) { if (fe_degree > -1) - internal::FEFaceEvaluationImplIntegrateScatterSelector< - dim, - Number, - VectorizedArrayType>::template run( - n_components, - integration_flag, - internal::get_beginning(destination), - shared_vector_data, - *this); + { + internal::FEFaceEvaluationImplIntegrateScatterSelector< + dim, + Number, + VectorizedArrayType>::template run( + n_components, + integration_flag, + internal::get_beginning(destination), + shared_vector_data, + *this); + } else - internal::FEFaceEvaluationGatherFactory< - dim, - Number, - VectorizedArrayType>::integrate(n_components, - integration_flag, - internal::get_beginning( - destination), - shared_vector_data, - *this); + { + internal::FEFaceEvaluationGatherFactory< + dim, + Number, + VectorizedArrayType>::integrate(n_components, + integration_flag, + internal::get_beginning( + destination), + shared_vector_data, + *this); + } } else { @@ -8227,7 +8244,7 @@ FEFaceEvaluationmatrix_info->get_mapping_info() + Assert(this->matrix_free->get_mapping_info() .face_data_by_cells[this->quad_no] .quadrature_point_offsets.empty() == false, internal::ExcMatrixFreeAccessToUninitializedMappingField( @@ -8235,12 +8252,12 @@ FEFaceEvaluationcell * GeometryInfo::faces_per_cell + this->face_no; AssertIndexRange(index, - this->matrix_info->get_mapping_info() + this->matrix_free->get_mapping_info() .face_data_by_cells[this->quad_no] .quadrature_point_offsets.size()); - return this->matrix_info->get_mapping_info() + return this->matrix_free->get_mapping_info() .face_data_by_cells[this->quad_no] - .quadrature_points[this->matrix_info->get_mapping_info() + .quadrature_points[this->matrix_free->get_mapping_info() .face_data_by_cells[this->quad_no] .quadrature_point_offsets[index] + q]; diff --git a/include/deal.II/matrix_free/fe_evaluation_data.h b/include/deal.II/matrix_free/fe_evaluation_data.h index 32e84a9514..cd0b0c6226 100644 --- a/include/deal.II/matrix_free/fe_evaluation_data.h +++ b/include/deal.II/matrix_free/fe_evaluation_data.h @@ -133,7 +133,7 @@ public: * use, select one of the derived classes FEEvaluation or FEFaceEvaluation * with their respective arguments. */ - FEEvaluationData(const ShapeInfoType &info, + FEEvaluationData(const ShapeInfoType &shape_info, const bool is_interior_face = true); /** @@ -577,7 +577,7 @@ protected: * initialization. Used in derived classes when this class is initialized * from a MatrixFree object. */ - FEEvaluationData(const InitializationData &info, + FEEvaluationData(const InitializationData &initialization_data, const bool is_interior_face, const unsigned int quad_no, const unsigned int first_selected_component); diff --git a/include/deal.II/matrix_free/mapping_info_storage.h b/include/deal.II/matrix_free/mapping_info_storage.h index 98c8a3cdc1..e70c65033b 100644 --- a/include/deal.II/matrix_free/mapping_info_storage.h +++ b/include/deal.II/matrix_free/mapping_info_storage.h @@ -300,7 +300,7 @@ namespace internal static UpdateFlags compute_update_flags( const UpdateFlags update_flags, - const std::vector> &quad = + const std::vector> &quads = std::vector>()); /** diff --git a/include/deal.II/matrix_free/mapping_info_storage.templates.h b/include/deal.II/matrix_free/mapping_info_storage.templates.h index 9a7a035000..806c337104 100644 --- a/include/deal.II/matrix_free/mapping_info_storage.templates.h +++ b/include/deal.II/matrix_free/mapping_info_storage.templates.h @@ -124,7 +124,7 @@ namespace internal UpdateFlags MappingInfoStorage::compute_update_flags( const UpdateFlags update_flags, - const std::vector> &quad) + const std::vector> &quads) { // this class is build around the evaluation of jacobians, so compute // them in any case. The Jacobians will be inverted manually. Since we @@ -145,19 +145,19 @@ namespace internal // one quadrature point on the first component, but more points on later // components, we need to have Jacobian gradients anyway in order to // determine whether the Jacobian is constant throughout a cell - if (quad.empty() == false) + if (quads.empty() == false) { bool formula_with_one_point = false; - for (unsigned int i = 0; i < quad[0].size(); ++i) - if (quad[0][i].size() == 1) + for (unsigned int i = 0; i < quads[0].size(); ++i) + if (quads[0][i].size() == 1) { formula_with_one_point = true; break; } if (formula_with_one_point == true) - for (unsigned int comp = 1; comp < quad.size(); ++comp) - for (unsigned int i = 0; i < quad[comp].size(); ++i) - if (quad[comp][i].size() > 1) + for (unsigned int comp = 1; comp < quads.size(); ++comp) + for (unsigned int i = 0; i < quads[comp].size(); ++i) + if (quads[comp][i].size() > 1) { new_flags |= update_jacobian_grads; }