From 5fd1ea95eaeadb1eaf4567d27351ed3ac398bf06 Mon Sep 17 00:00:00 2001 From: Peter Munch Date: Mon, 17 Jan 2022 21:27:18 +0100 Subject: [PATCH] Move jacobian_gradients into FEEvalData --- include/deal.II/matrix_free/fe_evaluation.h | 23 +++++++++---------- .../deal.II/matrix_free/fe_evaluation_data.h | 7 ++++++ 2 files changed, 18 insertions(+), 12 deletions(-) diff --git a/include/deal.II/matrix_free/fe_evaluation.h b/include/deal.II/matrix_free/fe_evaluation.h index 86e54949e2..7b1d2a1482 100644 --- a/include/deal.II/matrix_free/fe_evaluation.h +++ b/include/deal.II/matrix_free/fe_evaluation.h @@ -4811,10 +4811,7 @@ FEEvaluationBase:: // cell with general Jacobian else { - const auto &jac_grad = - this->mapping_data->jacobian_gradients - [1 - this->is_interior_face] - [this->mapping_data->data_index_offsets[this->cell] + q_point]; + const auto &jac_grad = this->jacobian_gradients[q_point]; for (unsigned int comp = 0; comp < n_components; ++comp) { VectorizedArrayType tmp[dim][dim]; @@ -4915,10 +4912,7 @@ FEEvaluationBase:: // cell with general Jacobian else { - const Tensor<1, dim *(dim + 1) / 2, Tensor<1, dim, VectorizedArrayType>> - &jac_grad = - this->mapping_data->jacobian_gradients - [0][this->mapping_data->data_index_offsets[this->cell] + q_point]; + const auto &jac_grad = this->jacobian_gradients[q_point]; for (unsigned int comp = 0; comp < n_components; ++comp) { // compute laplacian before the gradient because it needs to access @@ -5244,10 +5238,7 @@ FEEvaluationBase:: { const Tensor<2, dim, VectorizedArrayType> jac = this->jacobian[q_point]; const VectorizedArrayType JxW = this->J_value[q_point]; - const auto & jac_grad = - this->mapping_data->jacobian_gradients - [1 - this->is_interior_face] - [this->mapping_data->data_index_offsets[this->cell] + q_point]; + const auto &jac_grad = this->jacobian_gradients[q_point]; for (unsigned int comp = 0; comp < n_components; ++comp) { // 1. tmp = hessian_in(u) * J @@ -6902,6 +6893,8 @@ FEEvaluationmapping_data->data_index_offsets[cell_index]; this->jacobian = &this->mapping_data->jacobians[0][offsets]; this->J_value = &this->mapping_data->JxW_values[offsets]; + this->jacobian_gradients = + this->mapping_data->jacobian_gradients[0].data() + offsets; for (unsigned int i = 0; i < VectorizedArrayType::size(); ++i) this->cell_ids[i] = cell_index * VectorizedArrayType::size() + i; @@ -7638,6 +7631,9 @@ FEFaceEvaluationnormal_x_jacobian = &this->mapping_data ->normals_times_jacobians[!this->is_interior_face][offsets]; + this->jacobian_gradients = + this->mapping_data->jacobian_gradients[!this->is_interior_face].data() + + offsets; # ifdef DEBUG this->dof_values_initialized = false; @@ -7774,6 +7770,9 @@ FEFaceEvaluationmatrix_free->get_mapping_info() .face_data_by_cells[this->quad_no] .normals_times_jacobians[!this->is_interior_face][offsets]; + this->jacobian_gradients = + this->mapping_data->jacobian_gradients[!this->is_interior_face].data() + + offsets; # ifdef DEBUG this->dof_values_initialized = false; diff --git a/include/deal.II/matrix_free/fe_evaluation_data.h b/include/deal.II/matrix_free/fe_evaluation_data.h index de3cb37535..6d6f90cc55 100644 --- a/include/deal.II/matrix_free/fe_evaluation_data.h +++ b/include/deal.II/matrix_free/fe_evaluation_data.h @@ -614,6 +614,13 @@ protected: */ const Tensor<2, dim, Number> *jacobian; + /** + * A pointer to the gradients of the inverse Jacobian transformation of the + * present cell. + */ + const Tensor<1, dim *(dim + 1) / 2, Tensor<1, dim, Number>> + *jacobian_gradients; + /** * A pointer to the Jacobian determinant of the present cell. If on a * Cartesian cell or on a cell with constant Jacobian, this is just the -- 2.39.5