]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Remove domain_partitioner and range_partitioner from TrilinosWrappers 9952/head
authorDaniel Arndt <arndtd@ornl.gov>
Fri, 24 Apr 2020 17:11:40 +0000 (13:11 -0400)
committerDaniel Arndt <arndtd@ornl.gov>
Fri, 24 Apr 2020 17:39:52 +0000 (13:39 -0400)
doc/news/changes/incompatibilities/20200424DanielArndt [new file with mode: 0644]
include/deal.II/lac/trilinos_block_sparse_matrix.h
include/deal.II/lac/trilinos_sparse_matrix.h
include/deal.II/lac/trilinos_vector.h
source/lac/trilinos_block_sparse_matrix.cc
source/lac/trilinos_solver.cc
source/lac/trilinos_sparse_matrix.cc

diff --git a/doc/news/changes/incompatibilities/20200424DanielArndt b/doc/news/changes/incompatibilities/20200424DanielArndt
new file mode 100644 (file)
index 0000000..4923c86
--- /dev/null
@@ -0,0 +1,5 @@
+Removed: The deprecated member functions domain_partitioner(),
+vector_partitioner() and range_partitioner() in TrilinosWrappers classes have
+been removed.
+<br>
+(Daniel Arndt, 2020/04/24)
index ab8d8a3f27d687e952bb770e6bebc02556ce5b0d..743f9560f0fa2d9688b1c876f3484bf60c9e6125 100644 (file)
@@ -240,32 +240,6 @@ namespace TrilinosWrappers
     MPI_Comm
     get_mpi_communicator() const;
 
-    /**
-     * Return a vector of the underlying Trilinos Epetra_Map that sets the
-     * partitioning of the domain space of this block matrix, i.e., the
-     * partitioning of the individual block vectors this matrix has to be
-     * multiplied with.
-     *
-     * @deprecated Use the methods of the individual matrices based on
-     * IndexSet arguments.
-     */
-    DEAL_II_DEPRECATED
-    std::vector<Epetra_Map>
-    domain_partitioner() const;
-
-    /**
-     * Return a vector of the underlying Trilinos Epetra_Map that sets the
-     * partitioning of the range space of this block matrix, i.e., the
-     * partitioning of the individual block vectors that are the result from
-     * matrix-vector products.
-     *
-     * @deprecated Use the methods of the individual matrices based on
-     * IndexSet arguments.
-     */
-    DEAL_II_DEPRECATED
-    std::vector<Epetra_Map>
-    range_partitioner() const;
-
     /**
      * Return the partitioning of the domain space for the individual blocks of
      * this matrix, i.e., the partitioning of the block vectors this matrix has
index 7c864504d3de957857f2ccd887fc7f0c7b16a5a8..a6f4709fe31fe3cb55800c76d09213e001f195c8 100644 (file)
@@ -1820,52 +1820,6 @@ namespace TrilinosWrappers
     const Epetra_CrsGraph &
     trilinos_sparsity_pattern() const;
 
-    /**
-     * Return a const reference to the underlying Trilinos Epetra_Map that
-     * sets the partitioning of the domain space of this matrix, i.e., the
-     * partitioning of the vectors this matrix has to be multiplied with.
-     *
-     * @deprecated Use locally_owned_domain_indices() instead.
-     */
-    DEAL_II_DEPRECATED
-    const Epetra_Map &
-    domain_partitioner() const;
-
-    /**
-     * Return a const reference to the underlying Trilinos Epetra_Map that
-     * sets the partitioning of the range space of this matrix, i.e., the
-     * partitioning of the vectors that are result from matrix-vector
-     * products.
-     *
-     * @deprecated Use locally_owned_range_indices() instead.
-     */
-    DEAL_II_DEPRECATED
-    const Epetra_Map &
-    range_partitioner() const;
-
-    /**
-     * Return a const reference to the underlying Trilinos Epetra_Map that
-     * sets the partitioning of the matrix rows. Equal to the partitioning of
-     * the range.
-     *
-     * @deprecated Use locally_owned_range_indices() instead.
-     */
-    DEAL_II_DEPRECATED
-    const Epetra_Map &
-    row_partitioner() const;
-
-    /**
-     * Return a const reference to the underlying Trilinos Epetra_Map that
-     * sets the partitioning of the matrix columns. This is in general not
-     * equal to the partitioner Epetra_Map for the domain because of overlap
-     * in the matrix.
-     *
-     * @deprecated Usually not necessary. If desired, access it via the
-     * Epetra_CrsMatrix.
-     */
-    DEAL_II_DEPRECATED
-    const Epetra_Map &
-    col_partitioner() const;
     //@}
 
     /**
index ecbe5f31e12c2ef620b5029a10eca8d1d0215318..72567d048611ce7549d0bb9206e605016d3c196f 100644 (file)
@@ -1199,16 +1199,6 @@ namespace TrilinosWrappers
       Epetra_FEVector &
       trilinos_vector();
 
-      /**
-       * Return a const reference to the underlying Trilinos Epetra_Map that
-       * sets the parallel partitioning of the vector.
-       *
-       * @deprecated Use trilinos_partitioner() instead.
-       */
-      DEAL_II_DEPRECATED
-      const Epetra_Map &
-      vector_partitioner() const;
-
       /**
        * Return a const reference to the underlying Trilinos Epetra_BlockMap
        * that sets the parallel partitioning of the vector.
@@ -2143,15 +2133,6 @@ namespace TrilinosWrappers
 
 
 
-    inline const Epetra_Map &
-    Vector::vector_partitioner() const
-    {
-      // TODO A dynamic_cast fails here. This is suspicious.
-      return static_cast<const Epetra_Map &>(vector->Map()); // NOLINT
-    }
-
-
-
     inline const Epetra_BlockMap &
     Vector::trilinos_partitioner() const
     {
index 6c0f5b55acc2f79d3d163ce403e5494c7999d2f0..495871a0bf2050b7b59cb7258b4039f77c6b4f56 100644 (file)
@@ -327,37 +327,6 @@ namespace TrilinosWrappers
 
 
 
-  std::vector<Epetra_Map>
-  BlockSparseMatrix::domain_partitioner() const
-  {
-    Assert(this->n_block_cols() != 0, ExcNotInitialized());
-    Assert(this->n_block_rows() != 0, ExcNotInitialized());
-
-    std::vector<Epetra_Map> domain_partitioner;
-    for (size_type c = 0; c < this->n_block_cols(); ++c)
-      domain_partitioner.push_back(
-        this->sub_objects[0][c]->domain_partitioner());
-
-    return domain_partitioner;
-  }
-
-
-
-  std::vector<Epetra_Map>
-  BlockSparseMatrix::range_partitioner() const
-  {
-    Assert(this->n_block_cols() != 0, ExcNotInitialized());
-    Assert(this->n_block_rows() != 0, ExcNotInitialized());
-
-    std::vector<Epetra_Map> range_partitioner;
-    for (size_type r = 0; r < this->n_block_rows(); ++r)
-      range_partitioner.push_back(this->sub_objects[r][0]->range_partitioner());
-
-    return range_partitioner;
-  }
-
-
-
   MPI_Comm
   BlockSparseMatrix::get_mpi_communicator() const
   {
index 82094b14198e59504f4b60effff756bf044dcad0..acaf840ed40a9f65738adf9b729327eb26077bab 100644 (file)
@@ -185,9 +185,9 @@ namespace TrilinosWrappers
     Assert(A.trilinos_matrix().Filled(),
            ExcMessage("Matrix is not compressed. Call compress() method."));
 
-    Epetra_Vector ep_x(View, A.domain_partitioner(), x.begin());
+    Epetra_Vector ep_x(View, A.trilinos_matrix().DomainMap(), x.begin());
     Epetra_Vector ep_b(View,
-                       A.range_partitioner(),
+                       A.trilinos_matrix().RangeMap(),
                        const_cast<double *>(b.begin()));
 
     // We need an Epetra_LinearProblem object to let the AztecOO solver know
@@ -229,12 +229,14 @@ namespace TrilinosWrappers
   {
     // In case we call the solver with deal.II vectors, we create views of the
     // vectors in Epetra format.
-    AssertDimension(x.local_size(), A.domain_partitioner().NumMyElements());
-    AssertDimension(b.local_size(), A.range_partitioner().NumMyElements());
+    AssertDimension(x.local_size(),
+                    A.trilinos_matrix().DomainMap().NumMyElements());
+    AssertDimension(b.local_size(),
+                    A.trilinos_matrix().RangeMap().NumMyElements());
 
-    Epetra_Vector ep_x(View, A.domain_partitioner(), x.begin());
+    Epetra_Vector ep_x(View, A.trilinos_matrix().DomainMap(), x.begin());
     Epetra_Vector ep_b(View,
-                       A.range_partitioner(),
+                       A.trilinos_matrix().RangeMap(),
                        const_cast<double *>(b.begin()));
 
     // We need an Epetra_LinearProblem object to let the AztecOO solver know
@@ -883,9 +885,9 @@ namespace TrilinosWrappers
     Assert(b.size() == A.m(), ExcDimensionMismatch(b.size(), A.m()));
     Assert(A.local_range().second == A.m(),
            ExcMessage("Can only work in serial when using deal.II vectors."));
-    Epetra_Vector ep_x(View, A.domain_partitioner(), x.begin());
+    Epetra_Vector ep_x(View, A.trilinos_matrix().DomainMap(), x.begin());
     Epetra_Vector ep_b(View,
-                       A.range_partitioner(),
+                       A.trilinos_matrix().RangeMap(),
                        const_cast<double *>(b.begin()));
 
     // We need an Epetra_LinearProblem object to let the Amesos solver know
@@ -904,11 +906,13 @@ namespace TrilinosWrappers
     dealii::LinearAlgebra::distributed::Vector<double> &      x,
     const dealii::LinearAlgebra::distributed::Vector<double> &b)
   {
-    AssertDimension(x.local_size(), A.domain_partitioner().NumMyElements());
-    AssertDimension(b.local_size(), A.range_partitioner().NumMyElements());
-    Epetra_Vector ep_x(View, A.domain_partitioner(), x.begin());
+    AssertDimension(x.local_size(),
+                    A.trilinos_matrix().DomainMap().NumMyElements());
+    AssertDimension(b.local_size(),
+                    A.trilinos_matrix().RangeMap().NumMyElements());
+    Epetra_Vector ep_x(View, A.trilinos_matrix().DomainMap(), x.begin());
     Epetra_Vector ep_b(View,
-                       A.range_partitioner(),
+                       A.trilinos_matrix().RangeMap(),
                        const_cast<double *>(b.begin()));
 
     // We need an Epetra_LinearProblem object to let the Amesos solver know
index bd04d29cc61bf115c4ab3787bc51ffc8c2157059..8a26060cefc04da4c1e70e40bfa5269276b52cef 100644 (file)
@@ -490,7 +490,7 @@ namespace TrilinosWrappers
     if (needs_deep_copy)
       {
         column_space_map =
-          std_cxx14::make_unique<Epetra_Map>(rhs.domain_partitioner());
+          std_cxx14::make_unique<Epetra_Map>(rhs.trilinos_matrix().DomainMap());
 
         // release memory before reallocation
         matrix = std_cxx14::make_unique<Epetra_FECrsMatrix>(*rhs.matrix);
@@ -951,8 +951,8 @@ namespace TrilinosWrappers
     if (this == &sparse_matrix)
       return;
 
-    column_space_map =
-      std_cxx14::make_unique<Epetra_Map>(sparse_matrix.domain_partitioner());
+    column_space_map = std_cxx14::make_unique<Epetra_Map>(
+      sparse_matrix.trilinos_matrix().DomainMap());
     matrix.reset();
     nonlocal_matrix_exporter.reset();
     matrix = std_cxx14::make_unique<Epetra_FECrsMatrix>(
@@ -2286,16 +2286,16 @@ namespace TrilinosWrappers
         {
           Assert(inputleft.n() == inputright.m(),
                  ExcDimensionMismatch(inputleft.n(), inputright.m()));
-          Assert(inputleft.domain_partitioner().SameAs(
-                   inputright.range_partitioner()),
+          Assert(inputleft.trilinos_matrix().DomainMap().SameAs(
+                   inputright.trilinos_matrix().RangeMap()),
                  ExcMessage("Parallel partitioning of A and B does not fit."));
         }
       else
         {
           Assert(inputleft.m() == inputright.m(),
                  ExcDimensionMismatch(inputleft.m(), inputright.m()));
-          Assert(inputleft.range_partitioner().SameAs(
-                   inputright.range_partitioner()),
+          Assert(inputleft.trilinos_matrix().RangeMap().SameAs(
+                   inputright.trilinos_matrix().RangeMap()),
                  ExcMessage("Parallel partitioning of A and B does not fit."));
         }
 
@@ -2319,8 +2319,8 @@ namespace TrilinosWrappers
           mod_B = Teuchos::rcp(
             new Epetra_CrsMatrix(Copy, inputright.trilinos_sparsity_pattern()),
             true);
-          mod_B->FillComplete(inputright.domain_partitioner(),
-                              inputright.range_partitioner());
+          mod_B->FillComplete(inputright.trilinos_matrix().DomainMap(),
+                              inputright.trilinos_matrix().RangeMap());
           Assert(inputright.local_range() == V.local_range(),
                  ExcMessage("Parallel distribution of matrix B and vector V "
                             "does not match."));
@@ -2438,38 +2438,6 @@ namespace TrilinosWrappers
 
 
 
-  const Epetra_Map &
-  SparseMatrix::domain_partitioner() const
-  {
-    return matrix->DomainMap();
-  }
-
-
-
-  const Epetra_Map &
-  SparseMatrix::range_partitioner() const
-  {
-    return matrix->RangeMap();
-  }
-
-
-
-  const Epetra_Map &
-  SparseMatrix::row_partitioner() const
-  {
-    return matrix->RowMap();
-  }
-
-
-
-  const Epetra_Map &
-  SparseMatrix::col_partitioner() const
-  {
-    return matrix->ColMap();
-  }
-
-
-
   MPI_Comm
   SparseMatrix::get_mpi_communicator() const
   {

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.