From: Wolfgang Bangerth Date: Thu, 7 Mar 2024 17:01:55 +0000 (-0700) Subject: Also remove the specialization structure of MGTwoLevelTransferNonNested. X-Git-Tag: v9.6.0-rc1~493^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=1c466388a02e06dc6482017e991e459777ec1d6e;p=dealii.git Also remove the specialization structure of MGTwoLevelTransferNonNested. --- diff --git a/include/deal.II/multigrid/mg_transfer_global_coarsening.h b/include/deal.II/multigrid/mg_transfer_global_coarsening.h index f84b2ea795..6f6e637272 100644 --- a/include/deal.II/multigrid/mg_transfer_global_coarsening.h +++ b/include/deal.II/multigrid/mg_transfer_global_coarsening.h @@ -206,9 +206,10 @@ namespace MGTransferGlobalCoarseningTools } // namespace MGTransferGlobalCoarseningTools + /** * An abstract base class for transfer operators between two multigrid levels. - * Specialization for LinearAlgebra::distributed::Vector. The implementation of + * The implementation of * restriction and prolongation between levels is delegated to derived classes, * which implement prolongate_and_add_internal() and restrict_and_add_internal() * accordingly. @@ -651,54 +652,19 @@ private: /** * Class for transfer between two non-nested multigrid levels. - * */ template class MGTwoLevelTransferNonNested : public MGTwoLevelTransferBase { -public: - /** - * Perform prolongation. - */ - void - prolongate_and_add(VectorType &dst, const VectorType &src) const override; - - /** - * Perform restriction. - */ - void - restrict_and_add(VectorType &dst, const VectorType &src) const override; - - /** - * Perform interpolation of a solution vector from the fine level to the - * coarse level. This function is different from restriction, where a - * weighted residual is transferred to a coarser level (transposition of - * prolongation matrix). - */ - void - interpolate(VectorType &dst, const VectorType &src) const override; - - /** - * Return the memory consumption of the allocated memory in this class. - */ - std::size_t - memory_consumption() const override; -}; - - - -/** - * Class for transfer between two non-nested multigrid levels. - * - * Specialization for LinearAlgebra::distributed::Vector. - * - */ -template -class MGTwoLevelTransferNonNested> - : public MGTwoLevelTransferBase> -{ private: + static_assert( + std::is_same_v< + VectorType, + LinearAlgebra::distributed::Vector>, + "This class is currently only implemented for vectors of " + "type LinearAlgebra::distributed::Vector."); + + using Number = typename VectorType::value_type; using VectorizedArrayType = VectorizedArray; mg::SignalsNonNested signals_non_nested; @@ -752,6 +718,9 @@ public: bool enforce_all_points_found; }; + /** + * Constructor. + */ MGTwoLevelTransferNonNested(const AdditionalData &data = AdditionalData()); /** @@ -775,9 +744,7 @@ public: * prolongation matrix). */ void - interpolate( - LinearAlgebra::distributed::Vector &dst, - const LinearAlgebra::distributed::Vector &src) const override; + interpolate(VectorType &dst, const VectorType &src) const override; /** * Enable inplace vector operations if external and internal vectors @@ -826,17 +793,15 @@ protected: * Perform prolongation. */ void - prolongate_and_add_internal( - LinearAlgebra::distributed::Vector &dst, - const LinearAlgebra::distributed::Vector &src) const override; + prolongate_and_add_internal(VectorType &dst, + const VectorType &src) const override; /** * Perform restriction. */ void - restrict_and_add_internal( - LinearAlgebra::distributed::Vector &dst, - const LinearAlgebra::distributed::Vector &src) const override; + restrict_and_add_internal(VectorType &dst, + const VectorType &src) const override; private: /** @@ -844,18 +809,15 @@ private: */ template void - prolongate_and_add_internal_comp( - LinearAlgebra::distributed::Vector &dst, - const LinearAlgebra::distributed::Vector &src) const; + prolongate_and_add_internal_comp(VectorType &dst, + const VectorType &src) const; /** * Perform restriction for correct number of components. */ template void - restrict_and_add_internal_comp( - LinearAlgebra::distributed::Vector &dst, - const LinearAlgebra::distributed::Vector &src) const; + restrict_and_add_internal_comp(VectorType &dst, const VectorType &src) const; /** * Pointer to the DoFHandler object used during initialization. diff --git a/include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h b/include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h index 3025bd18a9..69cabe4a3b 100644 --- a/include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h +++ b/include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h @@ -4716,9 +4716,9 @@ namespace internal -template -MGTwoLevelTransferNonNested>:: - MGTwoLevelTransferNonNested(const AdditionalData &data) +template +MGTwoLevelTransferNonNested::MGTwoLevelTransferNonNested( + const AdditionalData &data) : additional_data(data) , rpe(typename Utilities::MPI::RemotePointEvaluation::AdditionalData( data.tolerance, @@ -4727,15 +4727,15 @@ MGTwoLevelTransferNonNested>:: {})) {} -template +template void -MGTwoLevelTransferNonNested>:: - reinit(const DoFHandler &dof_handler_fine, - const DoFHandler &dof_handler_coarse, - const Mapping &mapping_fine, - const Mapping &mapping_coarse, - const AffineConstraints &constraint_fine, - const AffineConstraints &constraint_coarse) +MGTwoLevelTransferNonNested::reinit( + const DoFHandler &dof_handler_fine, + const DoFHandler &dof_handler_coarse, + const Mapping &mapping_fine, + const Mapping &mapping_coarse, + const AffineConstraints &constraint_fine, + const AffineConstraints &constraint_coarse) { AssertThrow(dof_handler_coarse.get_fe().has_support_points(), ExcNotImplemented()); @@ -4889,13 +4889,12 @@ namespace internal -template +template template void -MGTwoLevelTransferNonNested>:: - prolongate_and_add_internal_comp( - LinearAlgebra::distributed::Vector &dst, - const LinearAlgebra::distributed::Vector &src) const +MGTwoLevelTransferNonNested::prolongate_and_add_internal_comp( + VectorType &dst, + const VectorType &src) const { using Traits = internal::FEPointEvaluation::EvaluatorTypeTraits; @@ -5001,12 +5000,11 @@ MGTwoLevelTransferNonNested>:: -template +template void -MGTwoLevelTransferNonNested>:: - prolongate_and_add_internal( - LinearAlgebra::distributed::Vector &dst, - const LinearAlgebra::distributed::Vector &src) const +MGTwoLevelTransferNonNested::prolongate_and_add_internal( + VectorType &dst, + const VectorType &src) const { if (this->fe_coarse->n_components() == 1) prolongate_and_add_internal_comp<1>(dst, src); @@ -5018,13 +5016,12 @@ MGTwoLevelTransferNonNested>:: -template +template template void -MGTwoLevelTransferNonNested>:: - restrict_and_add_internal_comp( - LinearAlgebra::distributed::Vector &dst, - const LinearAlgebra::distributed::Vector &src) const +MGTwoLevelTransferNonNested::restrict_and_add_internal_comp( + VectorType &dst, + const VectorType &src) const { using Traits = internal::FEPointEvaluation::EvaluatorTypeTraits; @@ -5126,12 +5123,11 @@ MGTwoLevelTransferNonNested>:: -template +template void -MGTwoLevelTransferNonNested>:: - restrict_and_add_internal( - LinearAlgebra::distributed::Vector &dst, - const LinearAlgebra::distributed::Vector &src) const +MGTwoLevelTransferNonNested::restrict_and_add_internal( + VectorType &dst, + const VectorType &src) const { if (this->fe_coarse->n_components() == 1) restrict_and_add_internal_comp<1>(dst, src); @@ -5143,11 +5139,11 @@ MGTwoLevelTransferNonNested>:: -template +template void -MGTwoLevelTransferNonNested>:: - interpolate(LinearAlgebra::distributed::Vector &dst, - const LinearAlgebra::distributed::Vector &src) const +MGTwoLevelTransferNonNested::interpolate( + VectorType &dst, + const VectorType &src) const { AssertThrow(false, ExcNotImplemented()); (void)dst; @@ -5156,9 +5152,9 @@ MGTwoLevelTransferNonNested>:: -template +template void -MGTwoLevelTransferNonNested>:: +MGTwoLevelTransferNonNested:: enable_inplace_operations_if_possible( const std::shared_ptr &external_partitioner_coarse, @@ -5175,10 +5171,9 @@ MGTwoLevelTransferNonNested>:: -template +template std::size_t -MGTwoLevelTransferNonNested>:: - memory_consumption() const +MGTwoLevelTransferNonNested::memory_consumption() const { std::size_t size = 0; @@ -5192,40 +5187,40 @@ MGTwoLevelTransferNonNested>:: -template +template boost::signals2::connection -MGTwoLevelTransferNonNested>:: - connect_prolongation_cell_loop(const std::function &slot) +MGTwoLevelTransferNonNested::connect_prolongation_cell_loop( + const std::function &slot) { return signals_non_nested.prolongation_cell_loop.connect(slot); } -template +template boost::signals2::connection -MGTwoLevelTransferNonNested>:: - connect_restriction_cell_loop(const std::function &slot) +MGTwoLevelTransferNonNested::connect_restriction_cell_loop( + const std::function &slot) { return signals_non_nested.restriction_cell_loop.connect(slot); } -template +template boost::signals2::connection -MGTwoLevelTransferNonNested>:: - connect_prolongation(const std::function &slot) +MGTwoLevelTransferNonNested::connect_prolongation( + const std::function &slot) { return signals_non_nested.prolongation.connect(slot); } -template +template boost::signals2::connection -MGTwoLevelTransferNonNested>:: - connect_restriction(const std::function &slot) +MGTwoLevelTransferNonNested::connect_restriction( + const std::function &slot) { return signals_non_nested.restriction.connect(slot); }