From: Maximilian Bergbauer Date: Wed, 7 Jun 2023 09:25:27 +0000 (+0200) Subject: Reuse MappingInfo X-Git-Tag: v9.5.0-rc1~116^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F15358%2Fhead;p=dealii.git Reuse MappingInfo --- diff --git a/include/deal.II/non_matching/mapping_info.h b/include/deal.II/non_matching/mapping_info.h index 85c7c512bd..aed1d2a12d 100644 --- a/include/deal.II/non_matching/mapping_info.h +++ b/include/deal.II/non_matching/mapping_info.h @@ -434,6 +434,12 @@ namespace NonMatching unsigned int compute_n_q_points(const unsigned int n_q_points_unvectorized); + /** + * Clears fields to make the object reusable. + */ + void + clear(); + /** * Resize the unit_point data field. */ @@ -673,6 +679,17 @@ namespace NonMatching + template + void + MappingInfo::clear() + { + n_q_points_unvectorized.clear(); + unit_points_index.clear(); + data_index_offsets.clear(); + } + + + template void MappingInfo::reinit( @@ -778,6 +795,8 @@ namespace NonMatching const std::vector> &quadrature_vector, const unsigned int n_unfiltered_cells) { + clear(); + do_cell_index_compression = n_unfiltered_cells != numbers::invalid_unsigned_int; @@ -873,6 +892,8 @@ namespace NonMatching const std::vector> &quadrature_vector, const unsigned int n_unfiltered_cells) { + clear(); + Assert( additional_data.use_global_weights == false, ExcMessage( @@ -976,6 +997,8 @@ namespace NonMatching const std::vector>> &quadrature_vector, const unsigned int n_unfiltered_cells) { + clear(); + do_cell_index_compression = n_unfiltered_cells != numbers::invalid_unsigned_int; diff --git a/tests/non_matching/mapping_info_03.cc b/tests/non_matching/mapping_info_03.cc index cda12083aa..795bad9a17 100644 --- a/tests/non_matching/mapping_info_03.cc +++ b/tests/non_matching/mapping_info_03.cc @@ -43,15 +43,14 @@ test() FE_Q fe(degree); MappingQ mapping(degree); + NonMatching::MappingInfo::AdditionalData additional_data; + additional_data.use_global_weights = true; + NonMatching::MappingInfo mapping_info(mapping, + update_JxW_values, + additional_data); deallog << "Check JxW faces..." << std::endl; { - NonMatching::MappingInfo::AdditionalData additional_data; - additional_data.use_global_weights = true; - NonMatching::MappingInfo mapping_info(mapping, - update_JxW_values, - additional_data); - // 1) build vector of quadratures std::vector>> quad_vec; // prescribe JxW (there is no meaning in the actual values, they just have @@ -96,13 +95,6 @@ test() deallog << "\n\nCheck JxW cells..." << std::endl; { - NonMatching::MappingInfo::AdditionalData additional_data; - additional_data.use_global_weights = true; - NonMatching::MappingInfo mapping_info(mapping, - update_JxW_values, - additional_data); - - // 1) build vector of quadratures std::vector> quad_vec; // prescribe JxW (there is no meaning in the actual values, they just have