From e712e34272260b31bfd33b78ee0fd7f874e2c072 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Mon, 27 Jul 2015 13:51:30 -0500 Subject: [PATCH] Revise internal functions of MappingQ1. The compute_fill* functions now also take FEValuesData as arguments, just like the fill_fe_values() functions. --- include/deal.II/fe/mapping.h | 15 ++++-- include/deal.II/fe/mapping_q1.h | 9 +--- source/fe/mapping_q.cc | 14 +----- source/fe/mapping_q1.cc | 86 +++++++++++++-------------------- 4 files changed, 50 insertions(+), 74 deletions(-) diff --git a/include/deal.II/fe/mapping.h b/include/deal.II/fe/mapping.h index d4b01051b3..6fc969415f 100644 --- a/include/deal.II/fe/mapping.h +++ b/include/deal.II/fe/mapping.h @@ -740,7 +740,10 @@ private: * return value of this function). * @param[in] quadrature A reference to the quadrature formula in use * for the current evaluation. This quadrature object is the same - * as the one used when creating the @p internal_data object + * as the one used when creating the @p internal_data object. The + * object is used both to map the location of quadrature points, + * as well as to compute the JxW values for each quadrature + * point (which involves the quadrature weights). * @param[in] internal_data A reference to an object previously * created by get_data() and that may be used to store information * the mapping can compute once on the reference cell. See the @@ -789,7 +792,10 @@ private: * information is requested. * @param[in] quadrature A reference to the quadrature formula in use * for the current evaluation. This quadrature object is the same - * as the one used when creating the @p internal_data object + * as the one used when creating the @p internal_data object. The + * object is used both to map the location of quadrature points, + * as well as to compute the JxW values for each quadrature + * point (which involves the quadrature weights). * @param[in] internal_data A reference to an object previously * created by get_data() and that may be used to store information * the mapping can compute once on the reference cell. See the @@ -822,7 +828,10 @@ private: * given cell for which information is requested. * @param[in] quadrature A reference to the quadrature formula in use * for the current evaluation. This quadrature object is the same - * as the one used when creating the @p internal_data object + * as the one used when creating the @p internal_data object. The + * object is used both to map the location of quadrature points, + * as well as to compute the JxW values for each quadrature + * point (which involves the quadrature weights). * @param[in] internal_data A reference to an object previously * created by get_data() and that may be used to store information * the mapping can compute once on the reference cell. See the diff --git a/include/deal.II/fe/mapping_q1.h b/include/deal.II/fe/mapping_q1.h index d3eb9006a9..f2b6f05242 100644 --- a/include/deal.II/fe/mapping_q1.h +++ b/include/deal.II/fe/mapping_q1.h @@ -405,13 +405,8 @@ public: const unsigned int npts, const DataSetDescriptor data_set, const std::vector &weights, - const InternalData &mapping_data, - std::vector > &quadrature_points, - std::vector &JxW_values, - std::vector > &boundary_form, - std::vector > &normal_vectors, - std::vector > &jacobians, - std::vector > &inverse_jacobians) const; + const InternalData &internal_data, + FEValuesData &output_data) const; /** * Compute shape values and/or derivatives. diff --git a/source/fe/mapping_q.cc b/source/fe/mapping_q.cc index f3bcb5394f..034c0540c8 100644 --- a/source/fe/mapping_q.cc +++ b/source/fe/mapping_q.cc @@ -343,12 +343,7 @@ fill_fe_face_values (const typename Triangulation::cell_iterator & n_q_points), quadrature.get_weights(), *p_data, - output_data.quadrature_points, - output_data.JxW_values, - output_data.boundary_forms, - output_data.normal_vectors, - output_data.jacobians, - output_data.inverse_jacobians); + output_data); } @@ -398,12 +393,7 @@ fill_fe_subface_values (const typename Triangulation::cell_iterato cell->subface_case(face_no)), quadrature.get_weights(), *p_data, - output_data.quadrature_points, - output_data.JxW_values, - output_data.boundary_forms, - output_data.normal_vectors, - output_data.jacobians, - output_data.inverse_jacobians); + output_data); } diff --git a/source/fe/mapping_q1.cc b/source/fe/mapping_q1.cc index c4cb6be7fb..b5bc80917d 100644 --- a/source/fe/mapping_q1.cc +++ b/source/fe/mapping_q1.cc @@ -924,21 +924,17 @@ namespace internal const unsigned int n_q_points, const std::vector &weights, const typename dealii::MappingQ1::InternalData &data, - std::vector &JxW_values, - std::vector > &boundary_forms, - std::vector > &normal_vectors, - std::vector > &jacobians, - std::vector > &inverse_jacobians) + FEValuesData &output_data) { const UpdateFlags update_flags(data.current_update_flags()); if (update_flags & update_boundary_forms) { - AssertDimension (boundary_forms.size(), n_q_points); + AssertDimension (output_data.boundary_forms.size(), n_q_points); if (update_flags & update_normal_vectors) - AssertDimension (normal_vectors.size(), n_q_points); + AssertDimension (output_data.normal_vectors.size(), n_q_points); if (update_flags & update_JxW_values) - AssertDimension (JxW_values.size(), n_q_points); + AssertDimension (output_data.JxW_values.size(), n_q_points); // map the unit tangentials to the real cell. checking for d!=dim-1 // eliminates compiler warnings regarding unsigned int expressions < @@ -970,14 +966,14 @@ namespace internal // fields (because it has only dim-1 components), but we // can still compute the boundary form by simply // looking at the number of the face - boundary_forms[i][0] = (face_no == 0 ? - -1 : +1); + output_data.boundary_forms[i][0] = (face_no == 0 ? + -1 : +1); break; case 2: - cross_product (boundary_forms[i], data.aux[0][i]); + cross_product (output_data.boundary_forms[i], data.aux[0][i]); break; case 3: - cross_product (boundary_forms[i], data.aux[0][i], data.aux[1][i]); + cross_product (output_data.boundary_forms[i], data.aux[0][i], data.aux[1][i]); break; default: Assert(false, ExcNotImplemented()); @@ -998,9 +994,9 @@ namespace internal if (dim==1) { // J is a tangent vector - boundary_forms[point] = data.contravariant[point].transpose()[0]; - boundary_forms[point] /= - (face_no == 0 ? -1. : +1.) * boundary_forms[point].norm(); + output_data.boundary_forms[point] = data.contravariant[point].transpose()[0]; + output_data.boundary_forms[point] /= + (face_no == 0 ? -1. : +1.) * output_data.boundary_forms[point].norm(); } @@ -1014,7 +1010,7 @@ namespace internal // then compute the face normal from the face tangent // and the cell normal: - cross_product (boundary_forms[point], + cross_product (output_data.boundary_forms[point], data.aux[0][point], cell_normal); } @@ -1026,31 +1022,32 @@ namespace internal if (update_flags & (update_normal_vectors | update_JxW_values)) - for (unsigned int i=0; i::subface_ratio( cell->subface_case(face_no), subface_no); - JxW_values[i] *= area_ratio; + output_data.JxW_values[i] *= area_ratio; } } if (update_flags & update_normal_vectors) - normal_vectors[i] = Point(boundary_forms[i] / boundary_forms[i].norm()); + output_data.normal_vectors[i] = Point(output_data.boundary_forms[i] / + output_data.boundary_forms[i].norm()); } if (update_flags & update_jacobians) for (unsigned int point=0; point void -MappingQ1::compute_fill_face ( - const typename Triangulation::cell_iterator &cell, - const unsigned int face_no, - const unsigned int subface_no, - const unsigned int n_q_points, - const DataSetDescriptor data_set, - const std::vector &weights, - const InternalData &data, - std::vector > &quadrature_points, - std::vector &JxW_values, - std::vector > &boundary_forms, - std::vector > &normal_vectors, - std::vector > &jacobians, - std::vector > &inverse_jacobians) const +MappingQ1:: +compute_fill_face (const typename Triangulation::cell_iterator &cell, + const unsigned int face_no, + const unsigned int subface_no, + const unsigned int n_q_points, + const DataSetDescriptor data_set, + const std::vector &weights, + const InternalData &internal_data, + FEValuesData &output_data) const { compute_fill (cell, n_q_points, data_set, CellSimilarity::none, - data, quadrature_points); + internal_data, + output_data.quadrature_points); internal::compute_fill_face (*this, cell, face_no, subface_no, n_q_points, - weights, data, - JxW_values, boundary_forms, normal_vectors, - jacobians, inverse_jacobians); + weights, internal_data, + output_data); } @@ -1112,12 +1104,7 @@ fill_fe_face_values (const typename Triangulation::cell_iterator & n_q_points), quadrature.get_weights(), data, - output_data.quadrature_points, - output_data.JxW_values, - output_data.boundary_forms, - output_data.normal_vectors, - output_data.jacobians, - output_data.inverse_jacobians); + output_data); } @@ -1150,12 +1137,7 @@ fill_fe_subface_values (const typename Triangulation::cell_iterato cell->subface_case(face_no)), quadrature.get_weights(), data, - output_data.quadrature_points, - output_data.JxW_values, - output_data.boundary_forms, - output_data.normal_vectors, - output_data.jacobians, - output_data.inverse_jacobians); + output_data); } -- 2.39.5