]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Enable a couple unit tests for LinearAlgebra::TpetraWrappers::Vector 16590/head
authorDaniel Arndt <arndtd@ornl.gov>
Fri, 2 Feb 2024 21:58:02 +0000 (16:58 -0500)
committerDaniel Arndt <arndtd@ornl.gov>
Sat, 3 Feb 2024 05:21:00 +0000 (00:21 -0500)
include/deal.II/lac/trilinos_tpetra_sparse_matrix.h
include/deal.II/lac/trilinos_tpetra_vector.h
tests/lac/linear_operator_04.cc
tests/lac/vector_reinit_02.cc

index 858115ca85fac56a181f89a3bccd2bb411555131..7e1eab0da52836f09e9e4c8d431f1d35229444c8 100644 (file)
@@ -598,6 +598,28 @@ namespace LinearAlgebra
       trilinos_rcp();
       /** @} */
 
+      /**
+       * @name Partitioners
+       */
+      /** @{ */
+
+      /**
+       * Return the partitioning of the domain space of this matrix, i.e., the
+       * partitioning of the vectors this matrix has to be multiplied with.
+       */
+      IndexSet
+      locally_owned_domain_indices() const;
+
+      /**
+       * Return the partitioning of the range space of this matrix, i.e., the
+       * partitioning of the vectors that are result from matrix-vector
+       * products.
+       */
+      IndexSet
+      locally_owned_range_indices() const;
+
+      /** @} */
+
       /**
        * @addtogroup Exceptions
        */
@@ -760,6 +782,24 @@ namespace LinearAlgebra
       return matrix;
     }
 
+
+
+    template <typename Number, typename MemorySpace>
+    inline IndexSet
+    SparseMatrix<Number, MemorySpace>::locally_owned_domain_indices() const
+    {
+      return IndexSet(matrix->getDomainMap());
+    }
+
+
+
+    template <typename Number, typename MemorySpace>
+    inline IndexSet
+    SparseMatrix<Number, MemorySpace>::locally_owned_range_indices() const
+    {
+      return IndexSet(matrix->getRangeMap());
+    }
+
   } // namespace TpetraWrappers
 
 } // namespace LinearAlgebra
index 1411617a8f68938810c2f37dcaf31c740079a13a..63f871badf87dca0704600cf9bf8d6690b6c8d1c 100644 (file)
@@ -1200,6 +1200,51 @@ namespace LinearAlgebra
 
 } // namespace LinearAlgebra
 
+
+namespace internal
+{
+  namespace LinearOperatorImplementation
+  {
+    template <typename>
+    class ReinitHelper;
+
+    /**
+     * A helper class used internally in linear_operator.h. Specialization for
+     * LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace>.
+     */
+    template <typename Number, typename MemorySpace>
+    class ReinitHelper<
+      LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace>>
+    {
+    public:
+      template <typename Matrix>
+      static void
+      reinit_range_vector(
+        const Matrix                                               &matrix,
+        LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace> &v,
+        bool omit_zeroing_entries)
+      {
+        v.reinit(matrix.locally_owned_range_indices(),
+                 matrix.get_mpi_communicator(),
+                 omit_zeroing_entries);
+      }
+
+      template <typename Matrix>
+      static void
+      reinit_domain_vector(
+        const Matrix                                               &matrix,
+        LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace> &v,
+        bool omit_zeroing_entries)
+      {
+        v.reinit(matrix.locally_owned_domain_indices(),
+                 matrix.get_mpi_communicator(),
+                 omit_zeroing_entries);
+      }
+    };
+
+  } // namespace LinearOperatorImplementation
+} /* namespace internal */
+
 /**
  * Declare dealii::LinearAlgebra::TpetraWrappers::Vector as distributed vector.
  */
index 598e4025fd9773577b248a0b20740c5d0e815e78..31e235a5a816f2ad1b8a2b3cfab884d50ecfdb51 100644 (file)
@@ -20,6 +20,8 @@
 #include <deal.II/lac/linear_operator.h>
 #include <deal.II/lac/trilinos_block_sparse_matrix.h>
 #include <deal.II/lac/trilinos_sparse_matrix.h>
+#include <deal.II/lac/trilinos_tpetra_sparse_matrix.h>
+#include <deal.II/lac/trilinos_tpetra_vector.h>
 #include <deal.II/lac/trilinos_vector.h>
 
 #include "../tests.h"
@@ -33,15 +35,27 @@ main(int argc, char *argv[])
   initlog();
   deallog << std::setprecision(10);
 
-  TrilinosWrappers::SparseMatrix a;
-
-  auto op_a  = linear_operator<TrilinosWrappers::MPI::Vector>(a);
-  auto op_a2 = linear_operator<TrilinosWrappers::MPI::Vector>(a);
-
-  TrilinosWrappers::BlockSparseMatrix b;
-
-  auto op_b  = linear_operator<TrilinosWrappers::MPI::BlockVector>(b);
-  auto op_b3 = linear_operator<TrilinosWrappers::MPI::BlockVector>(b);
+  {
+    TrilinosWrappers::SparseMatrix a;
+    auto op_a  = linear_operator<TrilinosWrappers::MPI::Vector>(a);
+    auto op_a2 = linear_operator<TrilinosWrappers::MPI::Vector>(a);
+  }
+
+#ifdef DEAL_II_TRILINOS_WITH_TPETRA
+  {
+    LinearAlgebra::TpetraWrappers::SparseMatrix<double> a;
+    auto                                                op_a =
+      linear_operator<LinearAlgebra::TpetraWrappers::Vector<double>>(a);
+    auto op_a2 =
+      linear_operator<LinearAlgebra::TpetraWrappers::Vector<double>>(a);
+  }
+#endif
+
+  {
+    TrilinosWrappers::BlockSparseMatrix b;
+    auto op_b  = linear_operator<TrilinosWrappers::MPI::BlockVector>(b);
+    auto op_b3 = linear_operator<TrilinosWrappers::MPI::BlockVector>(b);
+  }
 
   deallog << "OK" << std::endl;
 
index 3a2df37243254d850b621571b3f7037385d7007b..657e39c3ab991a2fcdff42bcf8439a40bba1405e 100644 (file)
@@ -21,6 +21,7 @@
 
 #include <deal.II/base/mpi.h>
 
+#include <deal.II/lac/trilinos_tpetra_vector.h>
 #include <deal.II/lac/trilinos_vector.h>
 #include <deal.II/lac/vector_memory.h>
 
@@ -113,4 +114,9 @@ main(int argc, char **argv)
 
   do_test<TrilinosWrappers::MPI::Vector>();
   do_test<TrilinosWrappers::MPI::Vector>();
+
+#ifdef DEAL_II_TRILINOS_WITH_TPETRA
+  do_test<LinearAlgebra::TpetraWrappers::Vector<double>>();
+  do_test<LinearAlgebra::TpetraWrappers::Vector<double>>();
+#endif
 }

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.