From: Wolfgang Bangerth Date: Mon, 27 Jul 2015 22:06:51 +0000 (-0500) Subject: Remove MappingQ1::compute_fill and simplify MappingQ1::compute_fill_face. X-Git-Tag: v8.4.0-rc2~701^2~3 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=256a8c0a7e4f54a4699435b0942ea3dbd4ebb9f6;p=dealii.git Remove MappingQ1::compute_fill and simplify MappingQ1::compute_fill_face. Also move around a few functions in the .cc file to make the flow easier to read. The old style used rather generic function names that don't really give away what concretely a function does. --- diff --git a/include/deal.II/fe/mapping_q1.h b/include/deal.II/fe/mapping_q1.h index f2b6f05242..81071b09a7 100644 --- a/include/deal.II/fe/mapping_q1.h +++ b/include/deal.II/fe/mapping_q1.h @@ -386,25 +386,14 @@ public: const unsigned int n_orig_q_points, InternalData &data) const; - /** - * Do the computation for the fill_* functions. - */ - void compute_fill (const typename Triangulation::cell_iterator &cell, - const unsigned int npts, - const DataSetDescriptor data_set, - const CellSimilarity::Similarity cell_similarity, - const InternalData &data, - std::vector > &quadrature_points) const; - /** * Do the computation for the fill_* functions. */ void compute_fill_face (const typename Triangulation::cell_iterator &cell, const unsigned int face_no, const unsigned int subface_no, - const unsigned int npts, const DataSetDescriptor data_set, - const std::vector &weights, + const Quadrature &quadrature, const InternalData &internal_data, FEValuesData &output_data) const; diff --git a/source/fe/mapping_q.cc b/source/fe/mapping_q.cc index 034c0540c8..edb1827ae4 100644 --- a/source/fe/mapping_q.cc +++ b/source/fe/mapping_q.cc @@ -334,14 +334,13 @@ fill_fe_face_values (const typename Triangulation::cell_iterator & const unsigned int n_q_points = quadrature.size(); this->compute_fill_face (cell, face_no, numbers::invalid_unsigned_int, - n_q_points, QProjector::DataSetDescriptor:: face (face_no, cell->face_orientation(face_no), cell->face_flip(face_no), cell->face_rotation(face_no), n_q_points), - quadrature.get_weights(), + quadrature, *p_data, output_data); } @@ -383,7 +382,6 @@ fill_fe_subface_values (const typename Triangulation::cell_iterato const unsigned int n_q_points = quadrature.size(); this->compute_fill_face (cell, face_no, subface_no, - n_q_points, QProjector::DataSetDescriptor:: subface (face_no, subface_no, cell->face_orientation(face_no), @@ -391,7 +389,7 @@ fill_fe_subface_values (const typename Triangulation::cell_iterato cell->face_rotation(face_no), n_q_points, cell->subface_case(face_no)), - quadrature.get_weights(), + quadrature, *p_data, output_data); } diff --git a/source/fe/mapping_q1.cc b/source/fe/mapping_q1.cc index b5bc80917d..7e090cfb08 100644 --- a/source/fe/mapping_q1.cc +++ b/source/fe/mapping_q1.cc @@ -603,121 +603,122 @@ MappingQ1::get_subface_data (const UpdateFlags update_flags, -template -void -MappingQ1::compute_fill (const typename Triangulation::cell_iterator &cell, - const unsigned int n_q_points, - const DataSetDescriptor data_set, - const CellSimilarity::Similarity cell_similarity, - const InternalData &data, - std::vector > &quadrature_points) const +namespace internal { - const UpdateFlags update_flags(data.current_update_flags()); - - // if necessary, recompute the - // support points of the - // transformation of this cell - // (note that we need to first - // check the triangulation pointer, - // since otherwise the second test - // might trigger an exception if - // the triangulations are not the - // same) - if ((data.mapping_support_points.size() == 0) - || - (&cell->get_triangulation() != - &data.cell_of_current_support_points->get_triangulation()) - || - (cell != data.cell_of_current_support_points)) - { - compute_mapping_support_points(cell, data.mapping_support_points); - data.cell_of_current_support_points = cell; - } - - // first compute quadrature points - if (update_flags & update_quadrature_points) - { - AssertDimension (quadrature_points.size(), n_q_points); - - for (unsigned int point=0; point result = (shape[0] * - data.mapping_support_points[0]); - for (unsigned int k=1; k + void + maybe_compute_q_points (const typename dealii::MappingQ1::DataSetDescriptor data_set, + const typename dealii::MappingQ1::InternalData &data, + std::vector > &quadrature_points) { - AssertDimension (data.contravariant.size(), n_q_points); + const UpdateFlags update_flags(data.current_update_flags()); - // if the current cell is just a - // translation of the previous one, no - // need to recompute jacobians... - if (cell_similarity != CellSimilarity::translation) + if (update_flags & update_quadrature_points) { - std::fill(data.contravariant.begin(), data.contravariant.end(), - DerivativeForm<1,dim,spacedim>()); - - Assert (data.n_shape_functions > 0, ExcInternalError()); - const Tensor<1,spacedim> *supp_pts = - &data.mapping_support_points[0]; - - for (unsigned int point=0; point *data_derv = - &data.derivative(point+data_set, 0); - - double result [spacedim][dim]; - - // peel away part of sum to avoid zeroing the - // entries and adding for the first time - for (unsigned int i=0; i result = (shape[0] * + data.mapping_support_points[0]); for (unsigned int k=1; k + void + maybe_update_Jacobians (const CellSimilarity::Similarity cell_similarity, + const typename dealii::MappingQ1::DataSetDescriptor data_set, + const typename dealii::MappingQ1::InternalData &data) + { + const UpdateFlags update_flags(data.current_update_flags()); + if (update_flags & update_contravariant_transformation) + // if the current cell is just a + // translation of the previous one, no + // need to recompute jacobians... + if (cell_similarity != CellSimilarity::translation) + { + const unsigned int n_q_points = data.contravariant.size(); + + std::fill(data.contravariant.begin(), data.contravariant.end(), + DerivativeForm<1,dim,spacedim>()); + + Assert (data.n_shape_functions > 0, ExcInternalError()); + const Tensor<1,spacedim> *supp_pts = + &data.mapping_support_points[0]; + + for (unsigned int point=0; point *data_derv = + &data.derivative(point+data_set, 0); + + double result [spacedim][dim]; + + // peel away part of sum to avoid zeroing the + // entries and adding for the first time + for (unsigned int i=0; i void MappingQ1::compute_mapping_support_points( @@ -748,10 +749,27 @@ fill_fe_values (const typename Triangulation::cell_iterator &cell, const unsigned int n_q_points=quadrature.size(); - compute_fill (cell, n_q_points, DataSetDescriptor::cell (), cell_similarity, - data, - output_data.quadrature_points); + // if necessary, recompute the support points of the transformation of this cell + // (note that we need to first check the triangulation pointer, since otherwise + // the second test might trigger an exception if the triangulations are not the + // same) + if ((data.mapping_support_points.size() == 0) + || + (&cell->get_triangulation() != + &data.cell_of_current_support_points->get_triangulation()) + || + (cell != data.cell_of_current_support_points)) + { + compute_mapping_support_points(cell, data.mapping_support_points); + data.cell_of_current_support_points = cell; + } + internal::maybe_compute_q_points (DataSetDescriptor::cell (), + data, + output_data.quadrature_points); + internal::maybe_update_Jacobians (cell_similarity, + DataSetDescriptor::cell (), + data); const UpdateFlags update_flags(data.current_update_flags()); const std::vector &weights=quadrature.get_weights(); @@ -915,16 +933,25 @@ namespace internal { namespace { + /** + * Depending on what information is called for in the update flags of the + * @p data object, compute the various pieces of information that is required + * by the fill_fe_face_values() and fill_fe_subface_values() functions. + * This function simply unifies the work that would be done by + * those two functions. + * + * The resulting data is put into the @p output_data argument. + */ template void - compute_fill_face (const dealii::MappingQ1 &mapping, - const typename dealii::Triangulation::cell_iterator &cell, - const unsigned int face_no, - const unsigned int subface_no, - const unsigned int n_q_points, - const std::vector &weights, - const typename dealii::MappingQ1::InternalData &data, - FEValuesData &output_data) + maybe_compute_face_data (const dealii::MappingQ1 &mapping, + const typename dealii::Triangulation::cell_iterator &cell, + const unsigned int face_no, + const unsigned int subface_no, + const unsigned int n_q_points, + const std::vector &weights, + const typename dealii::MappingQ1::InternalData &data, + FEValuesData &output_data) { const UpdateFlags update_flags(data.current_update_flags()); @@ -1060,21 +1087,38 @@ template 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 unsigned int face_no, + const unsigned int subface_no, + const DataSetDescriptor data_set, + const Quadrature &quadrature, const InternalData &internal_data, FEValuesData &output_data) const { - compute_fill (cell, n_q_points, data_set, CellSimilarity::none, - internal_data, - output_data.quadrature_points); - internal::compute_fill_face (*this, - cell, face_no, subface_no, n_q_points, - weights, internal_data, - output_data); + // if necessary, recompute the support points of the transformation of this cell + // (note that we need to first check the triangulation pointer, since otherwise + // the second test might trigger an exception if the triangulations are not the + // same) + if ((internal_data.mapping_support_points.size() == 0) + || + (&cell->get_triangulation() != + &internal_data.cell_of_current_support_points->get_triangulation()) + || + (cell != internal_data.cell_of_current_support_points)) + { + compute_mapping_support_points(cell, internal_data.mapping_support_points); + internal_data.cell_of_current_support_points = cell; + } + + internal::maybe_compute_q_points (data_set, + internal_data, + output_data.quadrature_points); + internal::maybe_update_Jacobians (CellSimilarity::none, + data_set, + internal_data); + internal::maybe_compute_face_data (*this, + cell, face_no, subface_no, quadrature.size(), + quadrature.get_weights(), internal_data, + output_data); } @@ -1096,13 +1140,12 @@ fill_fe_face_values (const typename Triangulation::cell_iterator & const unsigned int n_q_points = quadrature.size(); compute_fill_face (cell, face_no, numbers::invalid_unsigned_int, - n_q_points, DataSetDescriptor::face (face_no, cell->face_orientation(face_no), cell->face_flip(face_no), cell->face_rotation(face_no), n_q_points), - quadrature.get_weights(), + quadrature, data, output_data); } @@ -1128,14 +1171,13 @@ fill_fe_subface_values (const typename Triangulation::cell_iterato const unsigned int n_q_points = quadrature.size(); compute_fill_face (cell, face_no, subface_no, - n_q_points, DataSetDescriptor::subface (face_no, subface_no, cell->face_orientation(face_no), cell->face_flip(face_no), cell->face_rotation(face_no), n_q_points, cell->subface_case(face_no)), - quadrature.get_weights(), + quadrature, data, output_data); }