]> https://gitweb.dealii.org/ - dealii.git/commitdiff
avoid using Trilinos col/row_partitioner 1213/head
authorTimo Heister <timo.heister@gmail.com>
Tue, 28 Jul 2015 18:49:40 +0000 (14:49 -0400)
committerTimo Heister <timo.heister@gmail.com>
Tue, 28 Jul 2015 18:50:36 +0000 (14:50 -0400)
include/deal.II/lac/trilinos_sparse_matrix.h
source/lac/trilinos_sparse_matrix.cc
source/lac/trilinos_vector.cc

index 9d7a7a99503f684d073f60cb9966e900dc62ebd8..b9157880b04ed08f001420652333c12735b6f0aa 100644 (file)
@@ -2449,7 +2449,7 @@ namespace TrilinosWrappers
           {
             int ierr;
             ierr = matrix->GlobalAssemble(*column_space_map,
-                                          row_partitioner(), false);
+                                          matrix->RowMap(), false);
 
             Assert (ierr == 0, ExcTrilinosError(ierr));
             (void)ierr; // removes -Wunused-but-set-variable in optimized mode
index 66bc09313b0b771a230beee0f1df491ae21823ae..63f27345bdadb77979b2d2bd2814c2e0f3b10fb8 100644 (file)
@@ -1399,7 +1399,7 @@ namespace TrilinosWrappers
     // one is when the pattern is already fixed. In the former case, we add
     // the possibility to insert new values, and in the second we just replace
     // data.
-    if (row_partitioner().MyGID(static_cast<TrilinosWrappers::types::int_type>(row)) == true)
+    if (matrix->RowMap().MyGID(static_cast<TrilinosWrappers::types::int_type>(row)) == true)
       {
         if (matrix->Filled() == false)
           {
@@ -1589,7 +1589,7 @@ namespace TrilinosWrappers
     // If the calling processor owns the row to which we want to add values, we
     // can directly call the Epetra_CrsMatrix input function, which is much
     // faster than the Epetra_FECrsMatrix function.
-    if (row_partitioner().MyGID(static_cast<TrilinosWrappers::types::int_type>(row)) == true)
+    if (matrix->RowMap().MyGID(static_cast<TrilinosWrappers::types::int_type>(row)) == true)
       {
         ierr = matrix->Epetra_CrsMatrix::SumIntoGlobalValues(row, n_columns,
                                                              col_value_ptr,
@@ -1632,18 +1632,18 @@ namespace TrilinosWrappers
         std::cout << "------------------------------------------"
                   << std::endl;
         std::cout << "Got error " << ierr << " in row " << row
-                  << " of proc " << row_partitioner().Comm().MyPID()
+                  << " of proc " << matrix->RowMap().Comm().MyPID()
                   << " when trying to add the columns:" << std::endl;
         for (TrilinosWrappers::types::int_type i=0; i<n_columns; ++i)
           std::cout << col_index_ptr[i] << " ";
         std::cout << std::endl << std::endl;
         std::cout << "Matrix row "
-                  << (row_partitioner().MyGID(static_cast<TrilinosWrappers::types::int_type>(row)) == false ? "(nonlocal part)" : "")
+                  << (matrix->RowMap().MyGID(static_cast<TrilinosWrappers::types::int_type>(row)) == false ? "(nonlocal part)" : "")
                   << " has the following indices:" << std::endl;
         std::vector<TrilinosWrappers::types::int_type> indices;
         const Epetra_CrsGraph *graph =
           (nonlocal_matrix.get() != 0 &&
-           row_partitioner().MyGID(static_cast<TrilinosWrappers::types::int_type>(row)) == false) ?
+           matrix->RowMap().MyGID(static_cast<TrilinosWrappers::types::int_type>(row)) == false) ?
           &nonlocal_matrix->Graph() : &matrix->Graph();
 
         indices.resize(graph->NumGlobalIndices(static_cast<TrilinosWrappers::types::int_type>(row)));
@@ -1974,7 +1974,7 @@ namespace TrilinosWrappers
   TrilinosScalar
   SparseMatrix::matrix_norm_square (const VectorBase &v) const
   {
-    Assert (row_partitioner().SameAs(domain_partitioner()),
+    Assert (matrix->RowMap().SameAs(matrix->ColMap()),
             ExcNotQuadratic());
 
     VectorBase temp_vector;
@@ -1990,7 +1990,7 @@ namespace TrilinosWrappers
   SparseMatrix::matrix_scalar_product (const VectorBase &u,
                                        const VectorBase &v) const
   {
-    Assert (row_partitioner().SameAs(domain_partitioner()),
+    Assert (matrix->RowMap().SameAs(matrix->ColMap()),
             ExcNotQuadratic());
 
     VectorBase temp_vector;
@@ -2119,14 +2119,14 @@ namespace TrilinosWrappers
                                                                      indices);
               Assert (num_entries >= 0, ExcInternalError());
 #ifndef DEAL_II_WITH_64BIT_INDICES
-              const size_type GID = inputleft.row_partitioner().GID(i);
+              const size_type GID = inputleft.trilinos_matrix().RowMap().GID(i);
               for (TrilinosWrappers::types::int_type j=0; j<num_entries; ++j)
-                sparsity_transposed.add (inputleft.col_partitioner().GID(indices[j]),
+                sparsity_transposed.add (inputleft.trilinos_matrix().ColMap().GID(indices[j]),
                                          GID);
 #else
-              const size_type GID = inputleft.row_partitioner().GID64(i);
+              const size_type GID = inputleft.trilinos_matrix().RowMap().GID64(i);
               for (TrilinosWrappers::types::int_type j=0; j<num_entries; ++j)
-                sparsity_transposed.add (inputleft.col_partitioner().GID64(indices[j]),
+                sparsity_transposed.add (inputleft.trilinos_matrix().ColMap().GID64(indices[j]),
                                          GID);
 #endif
             }
@@ -2141,14 +2141,14 @@ namespace TrilinosWrappers
                                                            values, indices);
               Assert (num_entries >= 0, ExcInternalError());
 #ifndef DEAL_II_WITH_64BIT_INDICES
-              const size_type GID = inputleft.row_partitioner().GID(i);
+              const size_type GID = inputleft.trilinos_matrix().RowMap().GID(i);
               for (TrilinosWrappers::types::int_type j=0; j<num_entries; ++j)
-                transposed_mat.set (inputleft.col_partitioner().GID(indices[j]),
+                transposed_mat.set (inputleft.trilinos_matrix().ColMap().GID(indices[j]),
                                     GID, values[j]);
 #else
-              const size_type GID = inputleft.row_partitioner().GID64(i);
+              const size_type GID = inputleft.matrix->RowMap().GID64(i);
               for (TrilinosWrappers::types::int_type j=0; j<num_entries; ++j)
-                transposed_mat.set (inputleft.col_partitioner().GID64(indices[j]),
+                transposed_mat.set (inputleft.trilinos_matrix().ColMap().GID64(indices[j]),
                                     GID, values[j]);
 #endif
             }
index 6889333b4a35f25be0161657aa215334f0846f2b..a6bafe28f97955138a91f300603f3da3aee5424a 100644 (file)
@@ -474,10 +474,11 @@ namespace TrilinosWrappers
               ExcMessage ("The input vector has overlapping data, "
                           "which is not allowed."));
 
-      if (vector->Map().SameAs(m.col_partitioner()) == false)
+      if (vector->Map().SameAs(m.trilinos_matrix().ColMap()) == false)
         {
-          Epetra_Map map = m.col_partitioner();
-          vector.reset (new Epetra_FEVector(map));
+          vector.reset (new Epetra_FEVector(
+                          m.trilinos_matrix().ColMap()
+                        ));
         }
 
       Epetra_Import data_exchange (vector->Map(), v.vector->Map());

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.