From 8a20d60c1137b46c4da6306d7d449ece0892236e Mon Sep 17 00:00:00 2001 From: Daniel Arndt Date: Wed, 2 Jun 2021 11:52:53 -0400 Subject: [PATCH] Require MPI for Trilinos --- cmake/configure/configure_2_trilinos.cmake | 45 +++------ .../incompatibilities/20210602DanielArndt | 4 + include/deal.II/dofs/dof_accessor.templates.h | 2 +- include/deal.II/lac/read_write_vector.h | 10 +- .../deal.II/lac/read_write_vector.templates.h | 6 +- .../trilinos_epetra_communication_pattern.h | 10 +- include/deal.II/lac/trilinos_epetra_vector.h | 2 +- include/deal.II/lac/trilinos_precondition.h | 14 +-- include/deal.II/lac/trilinos_sparse_matrix.h | 12 +-- .../deal.II/lac/trilinos_sparsity_pattern.h | 8 +- .../trilinos_tpetra_communication_pattern.h | 2 +- include/deal.II/lac/trilinos_tpetra_vector.h | 2 +- .../lac/trilinos_tpetra_vector.templates.h | 18 ++-- include/deal.II/lac/trilinos_vector.h | 22 +---- include/deal.II/lac/vector_element_access.h | 2 +- .../numerics/data_out_dof_data.templates.h | 2 +- .../trilinos_epetra_communication_pattern.cc | 10 +- source/lac/trilinos_epetra_vector.cc | 42 ++++----- source/lac/trilinos_precondition.cc | 16 +--- source/lac/trilinos_sparse_matrix.cc | 93 +++---------------- source/lac/trilinos_sparsity_pattern.cc | 7 -- .../trilinos_tpetra_communication_pattern.cc | 10 +- source/lac/trilinos_tpetra_vector.cc | 6 +- source/lac/trilinos_vector.cc | 24 +---- 24 files changed, 93 insertions(+), 276 deletions(-) create mode 100644 doc/news/changes/incompatibilities/20210602DanielArndt diff --git a/cmake/configure/configure_2_trilinos.cmake b/cmake/configure/configure_2_trilinos.cmake index cc722180a5..fd415de1a1 100644 --- a/cmake/configure/configure_2_trilinos.cmake +++ b/cmake/configure/configure_2_trilinos.cmake @@ -17,6 +17,7 @@ # Configuration for the trilinos library: # +SET(FEATURE_TRILINOS_DEPENDS MPI) MACRO(FEATURE_TRILINOS_FIND_EXTERNAL var) FIND_PACKAGE(TRILINOS) @@ -87,41 +88,19 @@ MACRO(FEATURE_TRILINOS_FIND_EXTERNAL var) # Trilinos has to be configured with the same MPI configuration as # deal.II. # - IF( (TRILINOS_WITH_MPI AND NOT DEAL_II_WITH_MPI) - OR - (NOT TRILINOS_WITH_MPI AND DEAL_II_WITH_MPI)) + IF(NOT TRILINOS_WITH_MPI) MESSAGE(STATUS "Could not find a sufficient Trilinos installation: " - "Trilinos has to be configured with the same MPI configuration as deal.II." + "Trilinos has to have MPI support enabled." ) SET(TRILINOS_ADDITIONAL_ERROR_STRING ${TRILINOS_ADDITIONAL_ERROR_STRING} "The Trilinos installation (found at \"${TRILINOS_DIR}\")\n" - "has to be configured with the same MPI configuration as deal.II, but found:\n" - " DEAL_II_WITH_MPI = ${DEAL_II_WITH_MPI}\n" + "has to be configured with MPI support, but found:\n" " TRILINOS_WITH_MPI = ${TRILINOS_WITH_MPI}\n" ) SET(${var} FALSE) ENDIF() - # - # deal.II has to be configured with MPI if both Trilinos and PETSc are - # enabled. - # - IF(DEAL_II_WITH_TRILINOS AND DEAL_II_WITH_PETSC AND NOT DEAL_II_WITH_MPI) - MESSAGE(STATUS "Incompatible configuration settings: " - "MPI must be enabled to use both Trilinos and PETSc, as both libraries " - "provide mutually incompatible MPI stubs." - ) - SET(TRILINOS_ADDITIONAL_ERROR_STRING - ${TRILINOS_ADDITIONAL_ERROR_STRING} - "Incompatible Trilinos and PETSc libraries found. Both libraries were " - "configured without MPI support and cannot be used at the same time due " - "to incompatible MPI stub files. Either reconfigure deal.II, Trilinos, " - "and PETSc with MPI support, or disable one of the libraries.\n" - ) - SET(${var} FALSE) - ENDIF() - # # Trilinos has to be configured with 32bit indices if deal.II uses # unsigned int. @@ -391,15 +370,13 @@ MACRO(FEATURE_TRILINOS_CONFIGURE_EXTERNAL) "TrilinosWrappers::BlockSparseMatrix") SET(DEAL_II_EXPAND_TRILINOS_MPI_BLOCKVECTOR "TrilinosWrappers::MPI::BlockVector") SET(DEAL_II_EXPAND_TRILINOS_MPI_VECTOR "TrilinosWrappers::MPI::Vector") - IF (TRILINOS_WITH_MPI) - SET(DEAL_II_EXPAND_EPETRA_VECTOR "LinearAlgebra::EpetraWrappers::Vector") - IF (${DEAL_II_TRILINOS_WITH_TPETRA}) - SET(DEAL_II_EXPAND_TPETRA_VECTOR_DOUBLE "LinearAlgebra::TpetraWrappers::Vector") - SET(DEAL_II_EXPAND_TPETRA_VECTOR_FLOAT "LinearAlgebra::TpetraWrappers::Vector") - IF (${DEAL_II_WITH_COMPLEX_NUMBERS}) - SET(DEAL_II_EXPAND_TPETRA_VECTOR_COMPLEX_DOUBLE "LinearAlgebra::TpetraWrappers::Vector>") - SET(DEAL_II_EXPAND_TPETRA_VECTOR_COMPLEX_FLOAT "LinearAlgebra::TpetraWrappers::Vector>") - ENDIF() + SET(DEAL_II_EXPAND_EPETRA_VECTOR "LinearAlgebra::EpetraWrappers::Vector") + IF (${DEAL_II_TRILINOS_WITH_TPETRA}) + SET(DEAL_II_EXPAND_TPETRA_VECTOR_DOUBLE "LinearAlgebra::TpetraWrappers::Vector") + SET(DEAL_II_EXPAND_TPETRA_VECTOR_FLOAT "LinearAlgebra::TpetraWrappers::Vector") + IF (${DEAL_II_WITH_COMPLEX_NUMBERS}) + SET(DEAL_II_EXPAND_TPETRA_VECTOR_COMPLEX_DOUBLE "LinearAlgebra::TpetraWrappers::Vector>") + SET(DEAL_II_EXPAND_TPETRA_VECTOR_COMPLEX_FLOAT "LinearAlgebra::TpetraWrappers::Vector>") ENDIF() ENDIF() IF(${DEAL_II_TRILINOS_WITH_SACADO}) diff --git a/doc/news/changes/incompatibilities/20210602DanielArndt b/doc/news/changes/incompatibilities/20210602DanielArndt new file mode 100644 index 0000000000..1cccd041d7 --- /dev/null +++ b/doc/news/changes/incompatibilities/20210602DanielArndt @@ -0,0 +1,4 @@ +Changed: Trilinos support in deal.II requires both deal.II and Trilinos to be +configured with MPI support. +
+(Daniel Arndt, 2021/06/02) diff --git a/include/deal.II/dofs/dof_accessor.templates.h b/include/deal.II/dofs/dof_accessor.templates.h index b6f11ba067..7bf7c0333d 100644 --- a/include/deal.II/dofs/dof_accessor.templates.h +++ b/include/deal.II/dofs/dof_accessor.templates.h @@ -705,7 +705,7 @@ namespace internal -#if defined(DEAL_II_WITH_TRILINOS) && defined(DEAL_II_WITH_MPI) +#if defined(DEAL_II_WITH_TRILINOS) static std::vector sort_indices(const types::global_dof_index *v_begin, const types::global_dof_index *v_end) diff --git a/include/deal.II/lac/read_write_vector.h b/include/deal.II/lac/read_write_vector.h index 48d62e8e83..2e2fecad4e 100644 --- a/include/deal.II/lac/read_write_vector.h +++ b/include/deal.II/lac/read_write_vector.h @@ -218,7 +218,6 @@ namespace LinearAlgebra #ifdef DEAL_II_WITH_TRILINOS -# ifdef DEAL_II_WITH_MPI /** * Initialize this ReadWriteVector by supplying access to all locally * available entries in the given ghosted or non-ghosted vector. @@ -232,7 +231,6 @@ namespace LinearAlgebra */ void reinit(const TrilinosWrappers::MPI::Vector &trilinos_vec); -# endif #endif /** @@ -365,8 +363,7 @@ namespace LinearAlgebra const std::shared_ptr &communication_pattern = {}); -# ifdef DEAL_II_WITH_MPI -# ifdef DEAL_II_TRILINOS_WITH_TPETRA +# ifdef DEAL_II_TRILINOS_WITH_TPETRA /** * Imports all the elements present in the vector's IndexSet from the input * vector @p tpetra_vec. VectorOperation::values @p operation is used to @@ -380,7 +377,7 @@ namespace LinearAlgebra VectorOperation::values operation, const std::shared_ptr &communication_pattern = {}); -# endif +# endif /** * Imports all the elements present in the vector's IndexSet from the input @@ -395,7 +392,6 @@ namespace LinearAlgebra VectorOperation::values operation, const std::shared_ptr &communication_pattern = {}); -# endif #endif #ifdef DEAL_II_WITH_CUDA @@ -697,7 +693,7 @@ namespace LinearAlgebra void resize_val(const size_type new_allocated_size); -#if defined(DEAL_II_WITH_TRILINOS) && defined(DEAL_II_WITH_MPI) +#if defined(DEAL_II_WITH_TRILINOS) # ifdef DEAL_II_TRILINOS_WITH_TPETRA /** * Return a TpetraWrappers::CommunicationPattern and store it for future diff --git a/include/deal.II/lac/read_write_vector.templates.h b/include/deal.II/lac/read_write_vector.templates.h index 7334d813c8..59f92e4479 100644 --- a/include/deal.II/lac/read_write_vector.templates.h +++ b/include/deal.II/lac/read_write_vector.templates.h @@ -299,7 +299,7 @@ namespace LinearAlgebra -#if defined(DEAL_II_WITH_TRILINOS) && defined(DEAL_II_WITH_MPI) +#if defined(DEAL_II_WITH_TRILINOS) template void ReadWriteVector::reinit( @@ -570,7 +570,7 @@ namespace LinearAlgebra -#if defined(DEAL_II_WITH_TRILINOS) && defined(DEAL_II_WITH_MPI) +#if defined(DEAL_II_WITH_TRILINOS) # ifdef DEAL_II_TRILINOS_WITH_TPETRA template void @@ -1033,7 +1033,7 @@ namespace LinearAlgebra -#if defined(DEAL_II_WITH_TRILINOS) && defined(DEAL_II_WITH_MPI) +#if defined(DEAL_II_WITH_TRILINOS) # ifdef DEAL_II_TRILINOS_WITH_TPETRA template TpetraWrappers::CommunicationPattern diff --git a/include/deal.II/lac/trilinos_epetra_communication_pattern.h b/include/deal.II/lac/trilinos_epetra_communication_pattern.h index 9a55623a73..2168c0af61 100644 --- a/include/deal.II/lac/trilinos_epetra_communication_pattern.h +++ b/include/deal.II/lac/trilinos_epetra_communication_pattern.h @@ -21,13 +21,11 @@ #ifdef DEAL_II_WITH_TRILINOS -# ifdef DEAL_II_WITH_MPI +# include -# include +# include -# include - -# include +# include DEAL_II_NAMESPACE_OPEN @@ -90,8 +88,6 @@ namespace LinearAlgebra DEAL_II_NAMESPACE_CLOSE -# endif - #endif #endif diff --git a/include/deal.II/lac/trilinos_epetra_vector.h b/include/deal.II/lac/trilinos_epetra_vector.h index 8251cbc882..400dc2ec5d 100644 --- a/include/deal.II/lac/trilinos_epetra_vector.h +++ b/include/deal.II/lac/trilinos_epetra_vector.h @@ -19,7 +19,7 @@ #include -#if defined(DEAL_II_WITH_TRILINOS) && defined(DEAL_II_WITH_MPI) +#if defined(DEAL_II_WITH_TRILINOS) # include # include diff --git a/include/deal.II/lac/trilinos_precondition.h b/include/deal.II/lac/trilinos_precondition.h index b9bf3e112b..e02048e45f 100644 --- a/include/deal.II/lac/trilinos_precondition.h +++ b/include/deal.II/lac/trilinos_precondition.h @@ -26,19 +26,15 @@ # include # include -# include - -# ifdef DEAL_II_WITH_MPI -# include -# else -# include -# endif # include +# include # include # include # include # include +# include + // forward declarations # ifndef DOXYGEN class Ifpack_Preconditioner; @@ -239,11 +235,7 @@ namespace TrilinosWrappers * Internal communication pattern in case the matrix needs to be copied * from deal.II format. */ -# ifdef DEAL_II_WITH_MPI Epetra_MpiComm communicator; -# else - Epetra_SerialComm communicator; -# endif /** * Internal Trilinos map in case the matrix needs to be copied from diff --git a/include/deal.II/lac/trilinos_sparse_matrix.h b/include/deal.II/lac/trilinos_sparse_matrix.h index 204641022a..c73d4e4a36 100644 --- a/include/deal.II/lac/trilinos_sparse_matrix.h +++ b/include/deal.II/lac/trilinos_sparse_matrix.h @@ -38,19 +38,15 @@ # include # include # include +# include # include # include +# include # include # include # include # include -# ifdef DEAL_II_WITH_MPI -# include -# include -# else -# include -# endif DEAL_II_NAMESPACE_OPEN @@ -2348,11 +2344,7 @@ namespace TrilinosWrappers * Internal communication pattern in case the matrix needs to be copied * from deal.II format. */ -# ifdef DEAL_II_WITH_MPI Epetra_MpiComm communicator; -# else - Epetra_SerialComm communicator; -# endif /** * Epetra_Map that sets the partitioning of the domain space of diff --git a/include/deal.II/lac/trilinos_sparsity_pattern.h b/include/deal.II/lac/trilinos_sparsity_pattern.h index badaa307fd..4241dfa9d0 100644 --- a/include/deal.II/lac/trilinos_sparsity_pattern.h +++ b/include/deal.II/lac/trilinos_sparsity_pattern.h @@ -28,16 +28,12 @@ # include # include +# include +# include # include # include # include -# ifdef DEAL_II_WITH_MPI -# include -# include -# else -# include -# endif DEAL_II_NAMESPACE_OPEN diff --git a/include/deal.II/lac/trilinos_tpetra_communication_pattern.h b/include/deal.II/lac/trilinos_tpetra_communication_pattern.h index 4d6ee23800..fe14d69af0 100644 --- a/include/deal.II/lac/trilinos_tpetra_communication_pattern.h +++ b/include/deal.II/lac/trilinos_tpetra_communication_pattern.h @@ -19,7 +19,7 @@ #include -#if defined(DEAL_II_TRILINOS_WITH_TPETRA) && defined(DEAL_II_WITH_MPI) +#if defined(DEAL_II_TRILINOS_WITH_TPETRA) # include diff --git a/include/deal.II/lac/trilinos_tpetra_vector.h b/include/deal.II/lac/trilinos_tpetra_vector.h index 08857a4726..82fd1aa606 100644 --- a/include/deal.II/lac/trilinos_tpetra_vector.h +++ b/include/deal.II/lac/trilinos_tpetra_vector.h @@ -19,7 +19,7 @@ #include -#if defined(DEAL_II_TRILINOS_WITH_TPETRA) && defined(DEAL_II_WITH_MPI) +#if defined(DEAL_II_TRILINOS_WITH_TPETRA) # include # include diff --git a/include/deal.II/lac/trilinos_tpetra_vector.templates.h b/include/deal.II/lac/trilinos_tpetra_vector.templates.h index 25aec238d7..20ff03151c 100644 --- a/include/deal.II/lac/trilinos_tpetra_vector.templates.h +++ b/include/deal.II/lac/trilinos_tpetra_vector.templates.h @@ -22,19 +22,17 @@ #ifdef DEAL_II_TRILINOS_WITH_TPETRA -# ifdef DEAL_II_WITH_MPI +# include -# include +# include -# include +# include -# include +# include +# include +# include -# include -# include -# include - -# include +# include DEAL_II_NAMESPACE_OPEN @@ -684,8 +682,6 @@ namespace LinearAlgebra DEAL_II_NAMESPACE_CLOSE -# endif - #endif #endif diff --git a/include/deal.II/lac/trilinos_vector.h b/include/deal.II/lac/trilinos_vector.h index 4a73bb9a8f..c7a7239f8e 100644 --- a/include/deal.II/lac/trilinos_vector.h +++ b/include/deal.II/lac/trilinos_vector.h @@ -31,19 +31,15 @@ # include # include +# include +# include +# include +# include +# include # include # include # include -# ifdef DEAL_II_WITH_MPI // only if MPI is installed -# include -# include -# else -# include -# endif -# include -# include -# include DEAL_II_NAMESPACE_OPEN @@ -2167,18 +2163,10 @@ namespace TrilinosWrappers { static MPI_Comm comm; -# ifdef DEAL_II_WITH_MPI - const Epetra_MpiComm *mpi_comm = dynamic_cast(&vector->Map().Comm()); comm = mpi_comm->Comm(); -# else - - comm = MPI_COMM_SELF; - -# endif - return comm; } diff --git a/include/deal.II/lac/vector_element_access.h b/include/deal.II/lac/vector_element_access.h index 46087fe43d..ff0c2aa9af 100644 --- a/include/deal.II/lac/vector_element_access.h +++ b/include/deal.II/lac/vector_element_access.h @@ -78,7 +78,7 @@ namespace internal -#if defined(DEAL_II_WITH_TRILINOS) && defined(DEAL_II_WITH_MPI) +#if defined(DEAL_II_WITH_TRILINOS) template <> inline void ElementAccess::add( diff --git a/include/deal.II/numerics/data_out_dof_data.templates.h b/include/deal.II/numerics/data_out_dof_data.templates.h index ba4077e2e0..547603ce43 100644 --- a/include/deal.II/numerics/data_out_dof_data.templates.h +++ b/include/deal.II/numerics/data_out_dof_data.templates.h @@ -665,7 +665,7 @@ namespace internal -#if defined(DEAL_II_TRILINOS_WITH_TPETRA) && defined(DEAL_II_WITH_MPI) +#if defined(DEAL_II_TRILINOS_WITH_TPETRA) template <> inline void VectorHelper>::extract( diff --git a/source/lac/trilinos_epetra_communication_pattern.cc b/source/lac/trilinos_epetra_communication_pattern.cc index 868527d663..7c34e3ed0e 100644 --- a/source/lac/trilinos_epetra_communication_pattern.cc +++ b/source/lac/trilinos_epetra_communication_pattern.cc @@ -17,13 +17,11 @@ #ifdef DEAL_II_WITH_TRILINOS -# ifdef DEAL_II_WITH_MPI +# include -# include +# include -# include - -# include +# include DEAL_II_NAMESPACE_OPEN @@ -85,6 +83,4 @@ namespace LinearAlgebra DEAL_II_NAMESPACE_CLOSE -# endif - #endif diff --git a/source/lac/trilinos_epetra_vector.cc b/source/lac/trilinos_epetra_vector.cc index 57980bb36d..5671daf544 100644 --- a/source/lac/trilinos_epetra_vector.cc +++ b/source/lac/trilinos_epetra_vector.cc @@ -17,19 +17,17 @@ #ifdef DEAL_II_WITH_TRILINOS -# ifdef DEAL_II_WITH_MPI +# include -# include +# include -# include +# include -# include +# include +# include +# include -# include -# include -# include - -# include +# include DEAL_II_NAMESPACE_OPEN @@ -243,7 +241,7 @@ namespace LinearAlgebra Assert(this->size() == down_V.size(), ExcDimensionMismatch(this->size(), down_V.size())); -# if DEAL_II_TRILINOS_VERSION_GTE(11, 11, 0) +# if DEAL_II_TRILINOS_VERSION_GTE(11, 11, 0) Epetra_Import data_exchange(vector->Map(), down_V.trilinos_vector().Map()); const int ierr = vector->Import(down_V.trilinos_vector(), @@ -251,7 +249,7 @@ namespace LinearAlgebra Epetra_AddLocalAlso); Assert(ierr == 0, ExcTrilinosError(ierr)); (void)ierr; -# else +# else // In versions older than 11.11 the Import function is broken for // adding Hence, we provide a workaround in this case @@ -266,7 +264,7 @@ namespace LinearAlgebra ierr = vector->Update(1.0, dummy, 1.0); Assert(ierr == 0, ExcTrilinosError(ierr)); (void)ierr; -# endif +# endif } return *this; @@ -532,11 +530,11 @@ namespace LinearAlgebra Vector::size_type Vector::size() const { -# ifndef DEAL_II_WITH_64BIT_INDICES +# ifndef DEAL_II_WITH_64BIT_INDICES return vector->GlobalLength(); -# else +# else return vector->GlobalLength64(); -# endif +# endif } @@ -568,23 +566,23 @@ namespace LinearAlgebra // easy case: local range is contiguous if (vector->Map().LinearMap()) { -# ifndef DEAL_II_WITH_64BIT_INDICES +# ifndef DEAL_II_WITH_64BIT_INDICES is.add_range(vector->Map().MinMyGID(), vector->Map().MaxMyGID() + 1); -# else +# else is.add_range(vector->Map().MinMyGID64(), vector->Map().MaxMyGID64() + 1); -# endif +# endif } else if (vector->Map().NumMyElements() > 0) { const size_type n_indices = vector->Map().NumMyElements(); -# ifndef DEAL_II_WITH_64BIT_INDICES +# ifndef DEAL_II_WITH_64BIT_INDICES unsigned int *vector_indices = reinterpret_cast(vector->Map().MyGlobalElements()); -# else +# else size_type *vector_indices = reinterpret_cast(vector->Map().MyGlobalElements64()); -# endif +# endif is.add_indices(vector_indices, vector_indices + n_indices); } is.compress(); @@ -673,6 +671,4 @@ namespace LinearAlgebra DEAL_II_NAMESPACE_CLOSE -# endif - #endif diff --git a/source/lac/trilinos_precondition.cc b/source/lac/trilinos_precondition.cc index af0451586d..f52a1b0f5a 100644 --- a/source/lac/trilinos_precondition.cc +++ b/source/lac/trilinos_precondition.cc @@ -32,9 +32,7 @@ DEAL_II_NAMESPACE_OPEN namespace TrilinosWrappers { PreconditionBase::PreconditionBase() -# ifdef DEAL_II_WITH_MPI : communicator(MPI_COMM_SELF) -# endif {} @@ -42,12 +40,8 @@ namespace TrilinosWrappers PreconditionBase::PreconditionBase(const PreconditionBase &base) : Subscriptor() , preconditioner(base.preconditioner) - , -# ifdef DEAL_II_WITH_MPI - communicator(base.communicator) - , -# endif - vector_distributor(new Epetra_Map(*base.vector_distributor)) + , communicator(base.communicator) + , vector_distributor(new Epetra_Map(*base.vector_distributor)) {} @@ -56,9 +50,7 @@ namespace TrilinosWrappers PreconditionBase::clear() { preconditioner.reset(); -# ifdef DEAL_II_WITH_MPI communicator = MPI_COMM_SELF; -# endif vector_distributor.reset(); } @@ -66,11 +58,7 @@ namespace TrilinosWrappers MPI_Comm PreconditionBase::get_mpi_communicator() const { -# ifdef DEAL_II_WITH_MPI return communicator.Comm(); -# else - return MPI_COMM_SELF; -# endif } diff --git a/source/lac/trilinos_sparse_matrix.cc b/source/lac/trilinos_sparse_matrix.cc index 56520cb0bf..5fe22dadf9 100644 --- a/source/lac/trilinos_sparse_matrix.cc +++ b/source/lac/trilinos_sparse_matrix.cc @@ -76,7 +76,6 @@ namespace TrilinosWrappers return V.end(); } -# ifdef DEAL_II_WITH_MPI template <> double * begin(LinearAlgebra::EpetraWrappers::Vector &V) @@ -105,7 +104,7 @@ namespace TrilinosWrappers return V.trilinos_vector()[0] + V.trilinos_vector().MyLength(); } -# ifdef DEAL_II_TRILINOS_WITH_TPETRA +# ifdef DEAL_II_TRILINOS_WITH_TPETRA template Number * begin(LinearAlgebra::TpetraWrappers::Vector &V) @@ -135,7 +134,6 @@ namespace TrilinosWrappers return V.trilinos_vector().getData().get() + V.trilinos_vector().getLocalLength(); } -# endif # endif } // namespace internal @@ -2306,17 +2304,10 @@ namespace TrilinosWrappers MPI_Comm SparseMatrix::get_mpi_communicator() const { -# ifdef DEAL_II_WITH_MPI - const Epetra_MpiComm *mpi_comm = dynamic_cast(&matrix->RangeMap().Comm()); Assert(mpi_comm != nullptr, ExcInternalError()); return mpi_comm->Comm(); -# else - - return MPI_COMM_SELF; - -# endif } } // namespace TrilinosWrappers @@ -2325,35 +2316,13 @@ namespace TrilinosWrappers { namespace internal { - namespace - { -# ifndef DEAL_II_WITH_MPI - Epetra_Map - make_serial_Epetra_map(const IndexSet &serial_partitioning) - { - // See IndexSet::make_trilinos_map - return Epetra_Map( - TrilinosWrappers::types::int_type(serial_partitioning.size()), - TrilinosWrappers::types::int_type(serial_partitioning.n_elements()), - 0, - Epetra_SerialComm()); - } -# endif - } // namespace - namespace LinearOperatorImplementation { TrilinosPayload::TrilinosPayload() : use_transpose(false) - , -# ifdef DEAL_II_WITH_MPI - communicator(MPI_COMM_SELF) + , communicator(MPI_COMM_SELF) , domain_map(IndexSet().make_trilinos_map(communicator.Comm())) , range_map(IndexSet().make_trilinos_map(communicator.Comm())) -# else - domain_map(internal::make_serial_Epetra_map(IndexSet())) - , range_map(internal::make_serial_Epetra_map(IndexSet())) -# endif { vmult = [](Range &, const Domain &) { Assert(false, @@ -2386,21 +2355,13 @@ namespace TrilinosWrappers const TrilinosWrappers::SparseMatrix &matrix_exemplar, const TrilinosWrappers::SparseMatrix &matrix) : use_transpose(matrix_exemplar.trilinos_matrix().UseTranspose()) - , -# ifdef DEAL_II_WITH_MPI - communicator(matrix_exemplar.get_mpi_communicator()) + , communicator(matrix_exemplar.get_mpi_communicator()) , domain_map( matrix_exemplar.locally_owned_domain_indices().make_trilinos_map( communicator.Comm())) , range_map( matrix_exemplar.locally_owned_range_indices().make_trilinos_map( communicator.Comm())) -# else - domain_map(internal::make_serial_Epetra_map( - matrix_exemplar.locally_owned_domain_indices())) - , range_map(internal::make_serial_Epetra_map( - matrix_exemplar.locally_owned_range_indices())) -# endif { vmult = [&matrix_exemplar, &matrix](Range & tril_dst, const Domain &tril_src) { @@ -2471,21 +2432,13 @@ namespace TrilinosWrappers const TrilinosWrappers::SparseMatrix & matrix_exemplar, const TrilinosWrappers::PreconditionBase &preconditioner) : use_transpose(matrix_exemplar.trilinos_matrix().UseTranspose()) - , -# ifdef DEAL_II_WITH_MPI - communicator(matrix_exemplar.get_mpi_communicator()) + , communicator(matrix_exemplar.get_mpi_communicator()) , domain_map( matrix_exemplar.locally_owned_domain_indices().make_trilinos_map( communicator.Comm())) , range_map( matrix_exemplar.locally_owned_range_indices().make_trilinos_map( communicator.Comm())) -# else - domain_map(internal::make_serial_Epetra_map( - matrix_exemplar.locally_owned_domain_indices())) - , range_map(internal::make_serial_Epetra_map( - matrix_exemplar.locally_owned_range_indices())) -# endif { vmult = [&matrix_exemplar, &preconditioner](Range & tril_dst, const Domain &tril_src) { @@ -2592,19 +2545,11 @@ namespace TrilinosWrappers const TrilinosWrappers::PreconditionBase &preconditioner) : use_transpose( preconditioner_exemplar.trilinos_operator().UseTranspose()) - , -# ifdef DEAL_II_WITH_MPI - communicator(preconditioner_exemplar.get_mpi_communicator()) + , communicator(preconditioner_exemplar.get_mpi_communicator()) , domain_map(preconditioner_exemplar.locally_owned_domain_indices() .make_trilinos_map(communicator.Comm())) , range_map(preconditioner_exemplar.locally_owned_range_indices() .make_trilinos_map(communicator.Comm())) -# else - domain_map(internal::make_serial_Epetra_map( - preconditioner_exemplar.locally_owned_domain_indices())) - , range_map(internal::make_serial_Epetra_map( - preconditioner_exemplar.locally_owned_range_indices())) -# endif { vmult = [&preconditioner_exemplar, &preconditioner](Range &tril_dst, const Domain &tril_src) { @@ -2828,11 +2773,7 @@ namespace TrilinosWrappers MPI_Comm TrilinosPayload::get_mpi_communicator() const { -# ifdef DEAL_II_WITH_MPI return communicator.Comm(); -# else - return MPI_COMM_SELF; -# endif } @@ -3330,8 +3271,7 @@ namespace TrilinosWrappers dealii::LinearAlgebra::distributed::Vector &, const dealii::LinearAlgebra::distributed::Vector &) const; -# ifdef DEAL_II_WITH_MPI -# ifdef DEAL_II_TRILINOS_WITH_TPETRA +# ifdef DEAL_II_TRILINOS_WITH_TPETRA template void SparseMatrix::vmult( dealii::LinearAlgebra::TpetraWrappers::Vector &, @@ -3341,13 +3281,12 @@ namespace TrilinosWrappers SparseMatrix::vmult( dealii::LinearAlgebra::TpetraWrappers::Vector &, const dealii::LinearAlgebra::TpetraWrappers::Vector &) const; -# endif +# endif template void SparseMatrix::vmult( dealii::LinearAlgebra::EpetraWrappers::Vector &, const dealii::LinearAlgebra::EpetraWrappers::Vector &) const; -# endif template void SparseMatrix::Tvmult(MPI::Vector &, const MPI::Vector &) const; @@ -3361,8 +3300,7 @@ namespace TrilinosWrappers dealii::LinearAlgebra::distributed::Vector &, const dealii::LinearAlgebra::distributed::Vector &) const; -# ifdef DEAL_II_WITH_MPI -# ifdef DEAL_II_TRILINOS_WITH_TPETRA +# ifdef DEAL_II_TRILINOS_WITH_TPETRA template void SparseMatrix::Tvmult( dealii::LinearAlgebra::TpetraWrappers::Vector &, @@ -3372,13 +3310,12 @@ namespace TrilinosWrappers SparseMatrix::Tvmult( dealii::LinearAlgebra::TpetraWrappers::Vector &, const dealii::LinearAlgebra::TpetraWrappers::Vector &) const; -# endif +# endif template void SparseMatrix::Tvmult( dealii::LinearAlgebra::EpetraWrappers::Vector &, const dealii::LinearAlgebra::EpetraWrappers::Vector &) const; -# endif template void SparseMatrix::vmult_add(MPI::Vector &, const MPI::Vector &) const; @@ -3392,8 +3329,7 @@ namespace TrilinosWrappers dealii::LinearAlgebra::distributed::Vector &, const dealii::LinearAlgebra::distributed::Vector &) const; -# ifdef DEAL_II_WITH_MPI -# ifdef DEAL_II_TRILINOS_WITH_TPETRA +# ifdef DEAL_II_TRILINOS_WITH_TPETRA template void SparseMatrix::vmult_add( dealii::LinearAlgebra::TpetraWrappers::Vector &, @@ -3403,13 +3339,12 @@ namespace TrilinosWrappers SparseMatrix::vmult_add( dealii::LinearAlgebra::TpetraWrappers::Vector &, const dealii::LinearAlgebra::TpetraWrappers::Vector &) const; -# endif +# endif template void SparseMatrix::vmult_add( dealii::LinearAlgebra::EpetraWrappers::Vector &, const dealii::LinearAlgebra::EpetraWrappers::Vector &) const; -# endif template void SparseMatrix::Tvmult_add(MPI::Vector &, const MPI::Vector &) const; @@ -3423,8 +3358,7 @@ namespace TrilinosWrappers dealii::LinearAlgebra::distributed::Vector &, const dealii::LinearAlgebra::distributed::Vector &) const; -# ifdef DEAL_II_WITH_MPI -# ifdef DEAL_II_TRILINOS_WITH_TPETRA +# ifdef DEAL_II_TRILINOS_WITH_TPETRA template void SparseMatrix::Tvmult_add( dealii::LinearAlgebra::TpetraWrappers::Vector &, @@ -3434,13 +3368,12 @@ namespace TrilinosWrappers SparseMatrix::Tvmult_add( dealii::LinearAlgebra::TpetraWrappers::Vector &, const dealii::LinearAlgebra::TpetraWrappers::Vector &) const; -# endif +# endif template void SparseMatrix::Tvmult_add( dealii::LinearAlgebra::EpetraWrappers::Vector &, const dealii::LinearAlgebra::EpetraWrappers::Vector &) const; -# endif } // namespace TrilinosWrappers # endif // DOXYGEN diff --git a/source/lac/trilinos_sparsity_pattern.cc b/source/lac/trilinos_sparsity_pattern.cc index 719978f904..dfbcc6fe0d 100644 --- a/source/lac/trilinos_sparsity_pattern.cc +++ b/source/lac/trilinos_sparsity_pattern.cc @@ -942,17 +942,10 @@ namespace TrilinosWrappers MPI_Comm SparsityPattern::get_mpi_communicator() const { -# ifdef DEAL_II_WITH_MPI - const Epetra_MpiComm *mpi_comm = dynamic_cast(&graph->RangeMap().Comm()); Assert(mpi_comm != nullptr, ExcInternalError()); return mpi_comm->Comm(); -# else - - return MPI_COMM_SELF; - -# endif } diff --git a/source/lac/trilinos_tpetra_communication_pattern.cc b/source/lac/trilinos_tpetra_communication_pattern.cc index fb2c270904..cb6973cc2d 100644 --- a/source/lac/trilinos_tpetra_communication_pattern.cc +++ b/source/lac/trilinos_tpetra_communication_pattern.cc @@ -17,13 +17,11 @@ #ifdef DEAL_II_TRILINOS_WITH_TPETRA -# ifdef DEAL_II_WITH_MPI +# include -# include +# include -# include - -# include +# include DEAL_II_NAMESPACE_OPEN @@ -99,6 +97,4 @@ namespace LinearAlgebra DEAL_II_NAMESPACE_CLOSE -# endif - #endif diff --git a/source/lac/trilinos_tpetra_vector.cc b/source/lac/trilinos_tpetra_vector.cc index b269bdabcf..6308274e8c 100644 --- a/source/lac/trilinos_tpetra_vector.cc +++ b/source/lac/trilinos_tpetra_vector.cc @@ -16,7 +16,6 @@ #include #ifdef DEAL_II_TRILINOS_WITH_TPETRA -# ifdef DEAL_II_WITH_MPI DEAL_II_NAMESPACE_OPEN @@ -26,14 +25,13 @@ namespace LinearAlgebra { template class Vector; template class Vector; -# ifdef DEAL_II_WITH_COMPLEX_VALUES +# ifdef DEAL_II_WITH_COMPLEX_VALUES template class Vector>; template class Vector>; -# endif +# endif } // namespace TpetraWrappers } // namespace LinearAlgebra DEAL_II_NAMESPACE_CLOSE -# endif #endif diff --git a/source/lac/trilinos_vector.cc b/source/lac/trilinos_vector.cc index 53b4a28fc1..021da0fd94 100644 --- a/source/lac/trilinos_vector.cc +++ b/source/lac/trilinos_vector.cc @@ -140,11 +140,7 @@ namespace TrilinosWrappers { // When we clear the vector, reset the pointer and generate an empty // vector. -# ifdef DEAL_II_WITH_MPI Epetra_Map map(0, 0, Epetra_MpiComm(MPI_COMM_SELF)); -# else - Epetra_Map map(0, 0, Epetra_SerialComm()); -# endif has_ghosts = false; vector = std::make_unique(map); @@ -210,7 +206,6 @@ namespace TrilinosWrappers // version in case the underlying Epetra_MpiComm object is the same, // otherwise we might access an MPI_Comm object that has been // deleted -# ifdef DEAL_II_WITH_MPI const Epetra_MpiComm *my_comm = dynamic_cast(&vector->Comm()); const Epetra_MpiComm *v_comm = @@ -218,9 +213,6 @@ namespace TrilinosWrappers const bool same_communicators = my_comm != nullptr && v_comm != nullptr && my_comm->DataPtr() == v_comm->DataPtr(); -# else - const bool same_communicators = true; -# endif if (!same_communicators || vector->Map().SameAs(v.vector->Map()) == false) { @@ -267,7 +259,7 @@ namespace TrilinosWrappers last_action = Insert; } -# if defined(DEBUG) && defined(DEAL_II_WITH_MPI) +# if defined(DEBUG) const Epetra_MpiComm *comm_ptr = dynamic_cast(&(v.vector->Comm())); Assert(comm_ptr != nullptr, ExcInternalError()); @@ -343,7 +335,7 @@ namespace TrilinosWrappers } else vector = std::move(actual_vec); -# if defined(DEBUG) && defined(DEAL_II_WITH_MPI) +# if defined(DEBUG) const Epetra_MpiComm *comm_ptr = dynamic_cast(&(vector->Comm())); Assert(comm_ptr != nullptr, ExcInternalError()); @@ -422,7 +414,6 @@ namespace TrilinosWrappers // check equality for MPI communicators to avoid accessing a possibly // invalid MPI_Comm object -# ifdef DEAL_II_WITH_MPI const Epetra_MpiComm *my_comm = dynamic_cast(&vector->Comm()); const Epetra_MpiComm *v_comm = @@ -450,9 +441,6 @@ namespace TrilinosWrappers // else // same_communicators = true; // } -# else - const bool same_communicators = true; -# endif // distinguish three cases. First case: both vectors have the same // layout (just need to copy the local data, not reset the memory and @@ -612,7 +600,6 @@ namespace TrilinosWrappers # ifdef DEBUG -# ifdef DEAL_II_WITH_MPI // check that every process has decided to use the same mode. This will // otherwise result in undefined behavior in the call to // GlobalAssemble(). @@ -628,7 +615,6 @@ namespace TrilinosWrappers "this vector was an addition or a set operation. This will " "prevent the compress() operation from succeeding.")); -# endif # endif // Now pass over the information about what we did last to the vector. @@ -756,7 +742,6 @@ namespace TrilinosWrappers ++ptr; } -# ifdef DEAL_II_WITH_MPI // in parallel, check that the vector // is zero on _all_ processors. const Epetra_MpiComm *mpi_comm = @@ -764,9 +749,6 @@ namespace TrilinosWrappers Assert(mpi_comm != nullptr, ExcInternalError()); unsigned int num_nonzero = Utilities::MPI::sum(flag, mpi_comm->Comm()); return num_nonzero == 0; -# else - return flag == 0; -# endif } @@ -774,14 +756,12 @@ namespace TrilinosWrappers bool Vector::is_non_negative() const { -# ifdef DEAL_II_WITH_MPI // if this vector is a parallel one, then // we need to communicate to determine // the answer to the current // function. this still has to be // implemented AssertThrow(local_size() == size(), ExcNotImplemented()); -# endif // get a representation of the vector and // loop over all the elements TrilinosScalar *start_ptr; -- 2.39.5