From: Sebastian Kinnewig Date: Wed, 15 Nov 2023 12:08:20 +0000 (+0100) Subject: Introduce make_rcp(...) X-Git-Tag: relicensing~292^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F16270%2Fhead;p=dealii.git Introduce make_rcp(...) --- diff --git a/doc/news/changes/minor/20231115SebastianKinnewig b/doc/news/changes/minor/20231115SebastianKinnewig new file mode 100644 index 0000000000..881ab68b00 --- /dev/null +++ b/doc/news/changes/minor/20231115SebastianKinnewig @@ -0,0 +1,4 @@ +New: Introduce a new function, Utilities::Trilinos::internal::make_rcp(), +to create Teuchos::RCP objects, avoiding the usage of 'operator new'. +
+(Sebastian Kinnewig, 2023/11/15) diff --git a/include/deal.II/base/index_set.h b/include/deal.II/base/index_set.h index 28868eac53..8cc267ebba 100644 --- a/include/deal.II/base/index_set.h +++ b/include/deal.II/base/index_set.h @@ -21,6 +21,7 @@ #include #include #include +#include #include diff --git a/include/deal.II/base/trilinos_utilities.h b/include/deal.II/base/trilinos_utilities.h index 761a1683d0..bb3e1e5958 100644 --- a/include/deal.II/base/trilinos_utilities.h +++ b/include/deal.II/base/trilinos_utilities.h @@ -32,6 +32,10 @@ # endif #endif +#ifdef DEAL_II_TRILINOS_WITH_TPETRA +# include +#endif // DEAL_II_TRILINOS_WITH_TPETRA + DEAL_II_NAMESPACE_OPEN /** @@ -180,6 +184,35 @@ namespace Utilities duplicate_map(const Epetra_BlockMap &map, const Epetra_Comm &comm); } // namespace Trilinos #endif + + +#ifdef DEAL_II_TRILINOS_WITH_TPETRA + namespace Trilinos + { + namespace internal + { + /** + * Creates and returns a + * Teuchos::RCP + * object for type T. + * + * @note In Trilinos 14.0.0, the function + * Teuchos::make_rcp() + * was introduced, which should be preferred to this function. + */ +# if defined(DOXYGEN) || !DEAL_II_TRILINOS_VERSION_GTE(14, 0, 0) + template + Teuchos::RCP + make_rcp(Args &&...args); +# else + using Teuchos::make_rcp; +# endif // defined DOXYGEN || !DEAL_II_TRILINOS_VERSION_GTE(14, 0, 0) + } // namespace internal + } // namespace Trilinos +#endif // DEAL_II_TRILINOS_WITH_TPETRA + } // namespace Utilities DEAL_II_NAMESPACE_CLOSE diff --git a/include/deal.II/lac/trilinos_tpetra_vector.templates.h b/include/deal.II/lac/trilinos_tpetra_vector.templates.h index d5c4745fb2..06f32439f4 100644 --- a/include/deal.II/lac/trilinos_tpetra_vector.templates.h +++ b/include/deal.II/lac/trilinos_tpetra_vector.templates.h @@ -47,8 +47,11 @@ namespace LinearAlgebra template Vector::Vector() : Subscriptor() - , vector(new VectorType(Teuchos::rcp( - new MapType(0, 0, Utilities::Trilinos::tpetra_comm_self())))) + , vector(Utilities::Trilinos::internal::make_rcp( + Utilities::Trilinos::internal::make_rcp( + 0, + 0, + Utilities::Trilinos::tpetra_comm_self()))) {} @@ -56,7 +59,9 @@ namespace LinearAlgebra template Vector::Vector(const Vector &V) : Subscriptor() - , vector(new VectorType(V.trilinos_vector(), Teuchos::Copy)) + , vector(Utilities::Trilinos::internal::make_rcp( + V.trilinos_vector(), + Teuchos::Copy)) {} @@ -73,7 +78,7 @@ namespace LinearAlgebra Vector::Vector(const IndexSet ¶llel_partitioner, const MPI_Comm communicator) : Subscriptor() - , vector(new VectorType( + , vector(Utilities::Trilinos::internal::make_rcp( parallel_partitioner.make_tpetra_map_rcp(communicator, false))) {} @@ -89,7 +94,7 @@ namespace LinearAlgebra parallel_partitioner.make_tpetra_map_rcp(communicator, false); if (vector->getMap()->isSameAs(*input_map) == false) - vector = Teuchos::rcp(new VectorType(input_map)); + Utilities::Trilinos::internal::make_rcp(input_map); else if (omit_zeroing_entries == false) { vector->putScalar(0.); @@ -110,7 +115,7 @@ namespace LinearAlgebra Teuchos::RCP input_map = parallel_partitioner.make_tpetra_map_rcp(communicator, true); - vector = Teuchos::rcp(new VectorType(input_map)); + Utilities::Trilinos::internal::make_rcp(input_map); } @@ -188,7 +193,8 @@ namespace LinearAlgebra Tpetra::REPLACE); } else - vector = Teuchos::rcp(new VectorType(V.trilinos_vector())); + vector = Utilities::Trilinos::internal::make_rcp( + V.trilinos_vector()); } return *this; @@ -763,10 +769,10 @@ namespace LinearAlgebra const MPI_Comm mpi_comm) { source_stored_elements = source_index_set; - - tpetra_comm_pattern = - Teuchos::rcp(new TpetraWrappers::CommunicationPattern( - locally_owned_elements(), source_index_set, mpi_comm)); + tpetra_comm_pattern = Utilities::Trilinos::internal::make_rcp< + TpetraWrappers::CommunicationPattern>(locally_owned_elements(), + source_index_set, + mpi_comm); } } // namespace TpetraWrappers } // namespace LinearAlgebra diff --git a/source/base/index_set.cc b/source/base/index_set.cc index 5dd4b6d251..590a814b4a 100644 --- a/source/base/index_set.cc +++ b/source/base/index_set.cc @@ -995,17 +995,18 @@ IndexSet::make_tpetra_map_rcp(const MPI_Comm communicator, const bool linear = overlapping ? false : is_ascending_and_one_to_one(communicator); if (linear) - return Teuchos::RCP>( - new Tpetra::Map( - size(), - n_elements(), - 0, + return Utilities::Trilinos::internal::make_rcp< + Tpetra::Map>( + size(), + n_elements(), + 0, # ifdef DEAL_II_WITH_MPI - Teuchos::rcp(new Teuchos::MpiComm(communicator)) + Utilities::Trilinos::internal::make_rcp>( + communicator) # else - Teuchos::rcp(new Teuchos::Comm()) -# endif - )); + Utilities::Trilinos::internal::make_rcp>() +# endif // DEAL_WITH_MPI + ); else { const std::vector indices = get_index_vector(); @@ -1013,17 +1014,19 @@ IndexSet::make_tpetra_map_rcp(const MPI_Comm communicator, std::copy(indices.begin(), indices.end(), int_indices.begin()); const Teuchos::ArrayView arr_view( int_indices); - return Teuchos::RCP>( - new Tpetra::Map( - size(), - arr_view, - 0, + + return Utilities::Trilinos::internal::make_rcp< + Tpetra::Map>( + size(), + arr_view, + 0, # ifdef DEAL_II_WITH_MPI - Teuchos::rcp(new Teuchos::MpiComm(communicator)) + Utilities::Trilinos::internal::make_rcp>( + communicator) # else - Teuchos::rcp(new Teuchos::Comm()) -# endif - )); + Utilities::Trilinos::internal::make_rcp>() +# endif // DEAL_II_WITH_MPI + ); } } # endif diff --git a/source/base/trilinos_utilities.cc b/source/base/trilinos_utilities.cc index e9f2b105e5..9fb095e345 100644 --- a/source/base/trilinos_utilities.cc +++ b/source/base/trilinos_utilities.cc @@ -169,6 +169,23 @@ namespace Utilities } } // namespace Trilinos #endif + + +#ifdef DEAL_II_TRILINOS_WITH_TPETRA + namespace Trilinos + { +# if !DEAL_II_TRILINOS_VERSION_GTE(14, 0, 0) + template + Teuchos::RCP + make_rcp(Args &&...args) + { + return Teuchos::RCP(new T(std::forward(args)...)); + } +# endif // !DEAL_II_TRILINOS_VERSION_GTE(14, 0, 0) + + } // namespace Trilinos +#endif // DEAL_II_TRILINOS_WITH_TPETRA + } // namespace Utilities DEAL_II_NAMESPACE_CLOSE