From ee0147853e589e5a3a8b5086d7ec0f57f7a974a3 Mon Sep 17 00:00:00 2001 From: Simon Sticko Date: Fri, 13 Sep 2019 10:54:02 +0200 Subject: [PATCH] Add function MappingCartesian::fill_jacobian_derivatives(..) Presently all of the fill_fe_*values functions in MappingCartesian contain loops that fill the derivatives of the Jacobian with zeros. In order to remove duplication, add a function which does this instead. --- include/deal.II/fe/mapping_cartesian.h | 12 ++ source/fe/mapping_cartesian.cc | 169 ++++++++----------------- 2 files changed, 67 insertions(+), 114 deletions(-) diff --git a/include/deal.II/fe/mapping_cartesian.h b/include/deal.II/fe/mapping_cartesian.h index 7b18a9e518..89f115096c 100644 --- a/include/deal.II/fe/mapping_cartesian.h +++ b/include/deal.II/fe/mapping_cartesian.h @@ -328,6 +328,18 @@ private: const unsigned int face_no, const InternalData & data, std::vector> &normal_vectors) const; + + /** + * Since the Jacobian is constant for this mapping all derivatives of the + * Jacobian are identically zero. Fill these quantities with zeros if the + * corresponding update flags say that they should be updated. + */ + void + maybe_update_jacobian_derivatives( + const InternalData & data, + const CellSimilarity::Similarity cell_similarity, + internal::FEValuesImplementation::MappingRelatedData + &output_data) const; }; /*@}*/ diff --git a/source/fe/mapping_cartesian.cc b/source/fe/mapping_cartesian.cc index 96aa3d3574..f26f1b1aff 100644 --- a/source/fe/mapping_cartesian.cc +++ b/source/fe/mapping_cartesian.cc @@ -315,6 +315,58 @@ MappingCartesian::maybe_update_normal_vectors( +template +void +MappingCartesian::maybe_update_jacobian_derivatives( + const InternalData & data, + const CellSimilarity::Similarity cell_similarity, + internal::FEValuesImplementation::MappingRelatedData + &output_data) const +{ + if (cell_similarity != CellSimilarity::translation) + { + if (data.update_each & update_jacobian_grads) + for (unsigned int i = 0; i < output_data.jacobian_grads.size(); ++i) + output_data.jacobian_grads[i] = DerivativeForm<2, dim, spacedim>(); + + if (data.update_each & update_jacobian_pushed_forward_grads) + for (unsigned int i = 0; + i < output_data.jacobian_pushed_forward_grads.size(); + ++i) + output_data.jacobian_pushed_forward_grads[i] = Tensor<3, spacedim>(); + + if (data.update_each & update_jacobian_2nd_derivatives) + for (unsigned int i = 0; + i < output_data.jacobian_2nd_derivatives.size(); + ++i) + output_data.jacobian_2nd_derivatives[i] = + DerivativeForm<3, dim, spacedim>(); + + if (data.update_each & update_jacobian_pushed_forward_2nd_derivatives) + for (unsigned int i = 0; + i < output_data.jacobian_pushed_forward_2nd_derivatives.size(); + ++i) + output_data.jacobian_pushed_forward_2nd_derivatives[i] = + Tensor<4, spacedim>(); + + if (data.update_each & update_jacobian_3rd_derivatives) + for (unsigned int i = 0; + i < output_data.jacobian_3rd_derivatives.size(); + ++i) + output_data.jacobian_3rd_derivatives[i] = + DerivativeForm<4, dim, spacedim>(); + + if (data.update_each & update_jacobian_pushed_forward_3rd_derivatives) + for (unsigned int i = 0; + i < output_data.jacobian_pushed_forward_3rd_derivatives.size(); + ++i) + output_data.jacobian_pushed_forward_3rd_derivatives[i] = + Tensor<5, spacedim>(); + } +} + + + template CellSimilarity::Similarity MappingCartesian::fill_fe_values( @@ -361,51 +413,8 @@ MappingCartesian::fill_fe_values( for (unsigned int j = 0; j < dim; ++j) output_data.jacobians[i][j][j] = data.cell_extents[j]; } - // "compute" the derivative of the Jacobian at the quadrature - // points, which are all zero of course - if (data.update_each & update_jacobian_grads) - if (cell_similarity != CellSimilarity::translation) - for (unsigned int i = 0; i < output_data.jacobian_grads.size(); ++i) - output_data.jacobian_grads[i] = DerivativeForm<2, dim, spacedim>(); - - if (data.update_each & update_jacobian_pushed_forward_grads) - if (cell_similarity != CellSimilarity::translation) - for (unsigned int i = 0; - i < output_data.jacobian_pushed_forward_grads.size(); - ++i) - output_data.jacobian_pushed_forward_grads[i] = Tensor<3, spacedim>(); - - // "compute" the hessian of the Jacobian at the quadrature points, - // which are all also zero of course - if (data.update_each & update_jacobian_2nd_derivatives) - if (cell_similarity != CellSimilarity::translation) - for (unsigned int i = 0; i < output_data.jacobian_2nd_derivatives.size(); - ++i) - output_data.jacobian_2nd_derivatives[i] = - DerivativeForm<3, dim, spacedim>(); - if (data.update_each & update_jacobian_pushed_forward_2nd_derivatives) - if (cell_similarity != CellSimilarity::translation) - for (unsigned int i = 0; - i < output_data.jacobian_pushed_forward_2nd_derivatives.size(); - ++i) - output_data.jacobian_pushed_forward_2nd_derivatives[i] = - Tensor<4, spacedim>(); - - if (data.update_each & update_jacobian_3rd_derivatives) - if (cell_similarity != CellSimilarity::translation) - for (unsigned int i = 0; i < output_data.jacobian_3rd_derivatives.size(); - ++i) - output_data.jacobian_3rd_derivatives[i] = - DerivativeForm<4, dim, spacedim>(); - - if (data.update_each & update_jacobian_pushed_forward_3rd_derivatives) - if (cell_similarity != CellSimilarity::translation) - for (unsigned int i = 0; - i < output_data.jacobian_pushed_forward_3rd_derivatives.size(); - ++i) - output_data.jacobian_pushed_forward_3rd_derivatives[i] = - Tensor<5, spacedim>(); + maybe_update_jacobian_derivatives(data, cell_similarity, output_data); // "compute" inverse Jacobian at the quadrature points, which are // all the same @@ -481,41 +490,7 @@ MappingCartesian::fill_fe_face_values( output_data.jacobians[i][d][d] = data.cell_extents[d]; } - if (data.update_each & update_jacobian_grads) - for (unsigned int i = 0; i < output_data.jacobian_grads.size(); ++i) - output_data.jacobian_grads[i] = DerivativeForm<2, dim, spacedim>(); - - if (data.update_each & update_jacobian_pushed_forward_grads) - for (unsigned int i = 0; - i < output_data.jacobian_pushed_forward_grads.size(); - ++i) - output_data.jacobian_pushed_forward_grads[i] = Tensor<3, spacedim>(); - - if (data.update_each & update_jacobian_2nd_derivatives) - for (unsigned int i = 0; i < output_data.jacobian_2nd_derivatives.size(); - ++i) - output_data.jacobian_2nd_derivatives[i] = - DerivativeForm<3, dim, spacedim>(); - - if (data.update_each & update_jacobian_pushed_forward_2nd_derivatives) - for (unsigned int i = 0; - i < output_data.jacobian_pushed_forward_2nd_derivatives.size(); - ++i) - output_data.jacobian_pushed_forward_2nd_derivatives[i] = - Tensor<4, spacedim>(); - - if (data.update_each & update_jacobian_3rd_derivatives) - for (unsigned int i = 0; i < output_data.jacobian_3rd_derivatives.size(); - ++i) - output_data.jacobian_3rd_derivatives[i] = - DerivativeForm<4, dim, spacedim>(); - - if (data.update_each & update_jacobian_pushed_forward_3rd_derivatives) - for (unsigned int i = 0; - i < output_data.jacobian_pushed_forward_3rd_derivatives.size(); - ++i) - output_data.jacobian_pushed_forward_3rd_derivatives[i] = - Tensor<5, spacedim>(); + maybe_update_jacobian_derivatives(data, CellSimilarity::none, output_data); if (data.update_each & update_inverse_jacobians) for (unsigned int i = 0; i < output_data.inverse_jacobians.size(); ++i) @@ -593,41 +568,7 @@ MappingCartesian::fill_fe_subface_values( output_data.jacobians[i][d][d] = data.cell_extents[d]; } - if (data.update_each & update_jacobian_grads) - for (unsigned int i = 0; i < output_data.jacobian_grads.size(); ++i) - output_data.jacobian_grads[i] = DerivativeForm<2, dim, spacedim>(); - - if (data.update_each & update_jacobian_pushed_forward_grads) - for (unsigned int i = 0; - i < output_data.jacobian_pushed_forward_grads.size(); - ++i) - output_data.jacobian_pushed_forward_grads[i] = Tensor<3, spacedim>(); - - if (data.update_each & update_jacobian_2nd_derivatives) - for (unsigned int i = 0; i < output_data.jacobian_2nd_derivatives.size(); - ++i) - output_data.jacobian_2nd_derivatives[i] = - DerivativeForm<3, dim, spacedim>(); - - if (data.update_each & update_jacobian_pushed_forward_2nd_derivatives) - for (unsigned int i = 0; - i < output_data.jacobian_pushed_forward_2nd_derivatives.size(); - ++i) - output_data.jacobian_pushed_forward_2nd_derivatives[i] = - Tensor<4, spacedim>(); - - if (data.update_each & update_jacobian_3rd_derivatives) - for (unsigned int i = 0; i < output_data.jacobian_3rd_derivatives.size(); - ++i) - output_data.jacobian_3rd_derivatives[i] = - DerivativeForm<4, dim, spacedim>(); - - if (data.update_each & update_jacobian_pushed_forward_3rd_derivatives) - for (unsigned int i = 0; - i < output_data.jacobian_pushed_forward_3rd_derivatives.size(); - ++i) - output_data.jacobian_pushed_forward_3rd_derivatives[i] = - Tensor<5, spacedim>(); + maybe_update_jacobian_derivatives(data, CellSimilarity::none, output_data); if (data.update_each & update_inverse_jacobians) for (unsigned int i = 0; i < output_data.inverse_jacobians.size(); ++i) -- 2.39.5