]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fix compatibility issues of TpetraWrappers with older Trilinos versions. 16448/head
authorSebastian Kinnewig <kinnewig@ifam.uni-hannover.de>
Wed, 10 Jan 2024 12:06:29 +0000 (13:06 +0100)
committerSebastian Kinnewig <kinnewig@ifam.uni-hannover.de>
Wed, 10 Jan 2024 16:02:13 +0000 (17:02 +0100)
include/deal.II/lac/trilinos_tpetra_sparse_matrix.h
include/deal.II/lac/trilinos_tpetra_sparse_matrix.templates.h
include/deal.II/lac/trilinos_tpetra_sparsity_pattern.h
source/lac/trilinos_tpetra_sparsity_pattern.cc

index 6af4632bbfa71d91189b4f746748efc2625fd2d3..858115ca85fac56a181f89a3bccd2bb411555131 100644 (file)
@@ -107,9 +107,15 @@ namespace LinearAlgebra
       /**
        * Typedef for the NodeType
        */
+#  if DEAL_II_TRILINOS_VERSION_GTE(14, 2, 0)
       using NodeType = Tpetra::KokkosCompat::KokkosDeviceWrapperNode<
         typename MemorySpace::kokkos_space::execution_space,
         typename MemorySpace::kokkos_space>;
+#  else
+      using NodeType = Kokkos::Compat::KokkosDeviceWrapperNode<
+        typename MemorySpace::kokkos_space::execution_space,
+        typename MemorySpace::kokkos_space>;
+#  endif
 
       /**
        * Typedef for Tpetra::CrsMatrix
index e3f709d5b1bacd5a0973e1d41339aa685ad81796..451878c8d2425136505054f1385c4365b8d68991 100644 (file)
@@ -112,7 +112,11 @@ namespace LinearAlgebra
 
         std::vector<TrilinosWrappers::types::int_type> ghost_rows;
         Teuchos::Array<size_t>                         n_entries_per_row(
+#  if DEAL_II_TRILINOS_VERSION_GTE(14, 0, 0)
           row_space_map->getLocalNumElements());
+#  else
+          row_space_map->getNodeNumElements());
+#  endif
         {
           size_type own = 0;
           for (const auto global_row : relevant_rows)
@@ -132,8 +136,13 @@ namespace LinearAlgebra
         // doubles aside.
         Teuchos::RCP<GraphType<NodeType>> graph;
 
+#  if DEAL_II_TRILINOS_VERSION_GTE(12, 16, 0)
         graph = Utilities::Trilinos::internal::make_rcp<GraphType<NodeType>>(
-          row_space_map, n_entries_per_row());
+          row_space_map, n_entries_per_row);
+#  else
+        graph = Utilities::Trilinos::internal::make_rcp<GraphType<NodeType>>(
+          row_space_map, Teuchos::arcpFromArray(n_entries_per_row));
+#  endif
 
         // This functions assumes that the sparsity pattern sits on all
         // processors (completely). The parallel version uses a Tpetra graph
@@ -315,8 +324,13 @@ namespace LinearAlgebra
     {
       Teuchos::Array<size_t> n_entries_per_row_array(n_entries_per_row.begin(),
                                                      n_entries_per_row.end());
+#  if DEAL_II_TRILINOS_VERSION_GTE(12, 16, 0)
+      matrix = Utilities::Trilinos::internal::make_rcp<MatrixType>(
+        column_space_map, n_entries_per_row_array);
+#  else
       matrix = Utilities::Trilinos::internal::make_rcp<MatrixType>(
-        column_space_map, n_entries_per_row_array());
+        column_space_map, Teuchos::arcpFromArray(n_entries_per_row_array));
+#  endif
     }
 
 
@@ -349,9 +363,15 @@ namespace LinearAlgebra
     {
       Teuchos::Array<size_t> n_entries_per_row_array(n_entries_per_row.begin(),
                                                      n_entries_per_row.end());
+#  if DEAL_II_TRILINOS_VERSION_GTE(12, 16, 0)
+      matrix = Utilities::Trilinos::internal::make_rcp<MatrixType>(
+        row_parallel_partitioning.make_tpetra_map_rcp(communicator, false),
+        n_entries_per_row_array);
+#  else
       matrix = Utilities::Trilinos::internal::make_rcp<MatrixType>(
         row_parallel_partitioning.make_tpetra_map_rcp(communicator, false),
-        n_entries_per_row_array());
+        Teuchos::arcpFromArray(n_entries_per_row_array));
+#  endif
     }
 
 
@@ -406,7 +426,11 @@ namespace LinearAlgebra
     inline unsigned int
     SparseMatrix<Number, MemorySpace>::local_size() const
     {
+#  if DEAL_II_TRILINOS_VERSION_GTE(14, 0, 0)
       return matrix->getLocalNumRows();
+#  else
+      return matrix->getNodeNumRows();
+#  endif
     }
 
 
@@ -645,10 +669,19 @@ namespace LinearAlgebra
         }
       else
         {
+#  if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
           typename MatrixType::values_host_view_type     values;
           typename MatrixType::local_inds_host_view_type indices;
+#  else
+          Teuchos::ArrayView<const Number> values;
+          Teuchos::ArrayView<const int>    indices;
+#  endif
 
+#  if DEAL_II_TRILINOS_VERSION_GTE(14, 0, 0)
           for (size_t i = 0; i < matrix->getLocalNumRows(); ++i)
+#  else
+          for (size_t i = 0; i < matrix->getNodeNumRows(); ++i)
+#  endif
             {
               matrix->getLocalRowView(i, indices, values);
 
index 9271f41ae27c491ddbd22d88205b3e6589a41872..613890c63ef59ee5c6ec90582c743205b518fcc5 100644 (file)
@@ -89,9 +89,15 @@ namespace LinearAlgebra
         /**
          * Typedef for the NodeType
          */
+#  if DEAL_II_TRILINOS_VERSION_GTE(14, 2, 0)
         using NodeType = Tpetra::KokkosCompat::KokkosDeviceWrapperNode<
           typename MemorySpace::kokkos_space::execution_space,
           typename MemorySpace::kokkos_space>;
+#  else
+        using NodeType = Kokkos::Compat::KokkosDeviceWrapperNode<
+          typename MemorySpace::kokkos_space::execution_space,
+          typename MemorySpace::kokkos_space>;
+#  endif
 
         /**
          * Constructor.
@@ -195,9 +201,15 @@ namespace LinearAlgebra
         /**
          * Typedef for the NodeType
          */
+#  if DEAL_II_TRILINOS_VERSION_GTE(14, 2, 0)
         using NodeType = Tpetra::KokkosCompat::KokkosDeviceWrapperNode<
           typename MemorySpace::kokkos_space::execution_space,
           typename MemorySpace::kokkos_space>;
+#  else
+        using NodeType = Kokkos::Compat::KokkosDeviceWrapperNode<
+          typename MemorySpace::kokkos_space::execution_space,
+          typename MemorySpace::kokkos_space>;
+#  endif
 
         /**
          * Constructor. Create an iterator into the matrix @p matrix for the
@@ -307,9 +319,15 @@ namespace LinearAlgebra
       /**
        * Typedef for the NodeType
        */
+#  if DEAL_II_TRILINOS_VERSION_GTE(14, 2, 0)
       using NodeType = Tpetra::KokkosCompat::KokkosDeviceWrapperNode<
         typename MemorySpace::kokkos_space::execution_space,
         typename MemorySpace::kokkos_space>;
+#  else
+      using NodeType = Kokkos::Compat::KokkosDeviceWrapperNode<
+        typename MemorySpace::kokkos_space::execution_space,
+        typename MemorySpace::kokkos_space>;
+#  endif
 
       /**
        * Declare an alias for the iterator class.
index b831cff3b6611963feec5a8cc47e2d1fe4d351f0..84732f9fe923f5b1a08faa42a4695311e65d3561 100644 (file)
@@ -60,11 +60,16 @@ namespace LinearAlgebra
           {
             // get a representation of the present row
             std::size_t ncols;
+#  if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
             typename Tpetra::CrsGraph<
               int,
               dealii::types::signed_global_dof_index,
               NodeType>::nonconst_global_inds_host_view_type
               column_indices_view(colnum_cache->data(), colnum_cache->size());
+#  else
+            Teuchos::ArrayView<long long> column_indices_view(
+              colnum_cache->data(), colnum_cache->size());
+#  endif
             sparsity_pattern->graph->getGlobalRowCopy(this->a_row,
                                                       column_indices_view,
                                                       ncols);
@@ -407,14 +412,26 @@ namespace LinearAlgebra
             "that. Either use more MPI processes or recompile Trilinos "
             "with 'local ordinate = long long' "));
 
+#  if DEAL_II_TRILINOS_VERSION_GTE(12, 16, 0)
+        if (row_map->getComm()->getSize() > 1)
+          graph =
+            Utilities::Trilinos::internal::make_rcp<GraphType<MemorySpace>>(
+              row_map, n_entries_per_row);
+        else
+          graph =
+            Utilities::Trilinos::internal::make_rcp<GraphType<MemorySpace>>(
+              row_map, col_map, n_entries_per_row);
+#  else
         if (row_map->getComm()->getSize() > 1)
           graph =
             Utilities::Trilinos::internal::make_rcp<GraphType<MemorySpace>>(
-              row_map, n_entries_per_row());
+              row_map, Teuchos::arcpFromArray(n_entries_per_row));
         else
           graph =
             Utilities::Trilinos::internal::make_rcp<GraphType<MemorySpace>>(
-              row_map, col_map, n_entries_per_row());
+              row_map, col_map, Teuchos::arcpFromArray(n_entries_per_row));
+
+#  endif
 
         AssertDimension(sp.n_rows(), graph->getGlobalNumRows());
         AssertDimension(sp.n_cols(), graph->getGlobalNumEntries());
@@ -717,8 +734,13 @@ namespace LinearAlgebra
       Assert(column_space_map.get(), ExcInternalError());
       if (nonlocal_graph.get() != nullptr)
         {
+#  if DEAL_II_TRILINOS_VERSION_GTE(14, 0, 0)
           if (nonlocal_graph->getRowMap()->getLocalNumElements() > 0 &&
               column_space_map->getGlobalNumElements() > 0)
+#  else
+          if (nonlocal_graph->getRowMap()->getNodeNumElements() > 0 &&
+              column_space_map->getGlobalNumElements() > 0)
+#  endif
             {
               // Insert dummy element at (row, column) that corresponds to row 0
               // in local index counting.
@@ -731,16 +753,27 @@ namespace LinearAlgebra
               if (column_space_map->getGlobalNumElements() ==
                   graph->getRangeMap()->getGlobalNumElements())
                 column = row;
-              // if not, take a column index that we have ourselves since we
-              // know for sure it is there (and it will not create spurious
-              // messages to many ranks like putting index 0 on many processors)
+                // if not, take a column index that we have ourselves since we
+                // know for sure it is there (and it will not create spurious
+                // messages to many ranks like putting index 0 on many
+                // processors)
+#  if DEAL_II_TRILINOS_VERSION_GTE(14, 0, 0)
               else if (column_space_map->getLocalNumElements() > 0)
+#  else
+              else if (column_space_map->getNodeNumElements() > 0)
+#  endif
                 column = column_space_map->getGlobalElement(0);
               nonlocal_graph->insertGlobalIndices(row, 1, &column);
             }
+#  if DEAL_II_TRILINOS_VERSION_GTE(14, 0, 0)
           Assert(nonlocal_graph->getRowMap()->getLocalNumElements() == 0 ||
                    column_space_map->getGlobalNumElements() == 0,
                  ExcInternalError());
+#  else
+          Assert(nonlocal_graph->getRowMap()->getNodeNumElements() == 0 ||
+                   column_space_map->getGlobalNumElements() == 0,
+                 ExcInternalError());
+#  endif
 
           nonlocal_graph->fillComplete(column_space_map, graph->getRangeMap());
           graph->fillComplete(column_space_map, graph->getRangeMap());
@@ -776,6 +809,7 @@ namespace LinearAlgebra
     SparsityPattern<MemorySpace>::exists(const size_type i,
                                          const size_type j) const
     {
+#  if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
       if (!row_is_stored_locally(i))
         return false;
 
@@ -796,6 +830,10 @@ namespace LinearAlgebra
         col_indices.data();
 
       return static_cast<size_t>(local_col_index) != col_indices.size();
+#  else
+      Assert(false, ExcNotImplemented());
+      return false;
+#  endif
     }
 
 
@@ -807,7 +845,12 @@ namespace LinearAlgebra
       size_type local_b = 0;
       for (int i = 0; i < static_cast<int>(local_size()); ++i)
         {
+#  if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
           typename GraphType::local_inds_host_view_type indices;
+#  else
+          Teuchos::ArrayView<const int> indices;
+#  endif
+
           graph->getLocalRowView(i, indices);
           const auto num_entries = indices.size();
           for (unsigned int j = 0; j < static_cast<unsigned int>(num_entries);
@@ -831,7 +874,11 @@ namespace LinearAlgebra
     unsigned int
     SparsityPattern<MemorySpace>::local_size() const
     {
+#  if DEAL_II_TRILINOS_VERSION_GTE(14, 0, 0)
       return graph->getLocalNumRows();
+#  else
+      return graph->getNodeNumRows();
+#  endif
     }
 
 
@@ -862,7 +909,11 @@ namespace LinearAlgebra
     unsigned int
     SparsityPattern<MemorySpace>::max_entries_per_row() const
     {
+#  if DEAL_II_TRILINOS_VERSION_GTE(14, 0, 0)
       return graph->getLocalMaxNumRowEntries();
+#  else
+      return graph->getNodeMaxNumRowEntries();
+#  endif
     }
 
 
@@ -950,9 +1001,17 @@ namespace LinearAlgebra
         out << *graph;
       else
         {
+#  if DEAL_II_TRILINOS_VERSION_GTE(14, 0, 0)
           for (unsigned int i = 0; i < graph->getLocalNumRows(); ++i)
+#  else
+          for (unsigned int i = 0; i < graph->getNodeNumRows(); ++i)
+#  endif
             {
+#  if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
               typename GraphType::local_inds_host_view_type indices;
+#  else
+              Teuchos::ArrayView<const int> indices;
+#  endif
               graph->getLocalRowView(i, indices);
               int num_entries = indices.size();
               for (int j = 0; j < num_entries; ++j)
@@ -975,7 +1034,11 @@ namespace LinearAlgebra
       for (dealii::types::signed_global_dof_index row = 0; row < local_size();
            ++row)
         {
+#  if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
           typename GraphType::local_inds_host_view_type indices;
+#  else
+          Teuchos::ArrayView<const int> indices;
+#  endif
           graph->getLocalRowView(row, indices);
           int num_entries = indices.size();
 

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.