]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Use get_index_vector() instead of fill_index_vector().
authorWolfgang Bangerth <bangerth@colostate.edu>
Thu, 26 Jan 2023 23:14:25 +0000 (16:14 -0700)
committerWolfgang Bangerth <bangerth@colostate.edu>
Thu, 26 Jan 2023 23:30:58 +0000 (16:30 -0700)
include/deal.II/lac/parpack_solver.h
include/deal.II/matrix_free/operators.h
source/base/index_set.cc
source/base/partitioner.cc
source/lac/petsc_parallel_vector.cc
source/lac/trilinos_sparse_matrix.cc

index c4ac8fa0758981868df51612ab7331d12ab9a30c..f111a07bf309b7351d7085070eeddfbb17f51ad0 100644 (file)
@@ -686,7 +686,7 @@ void
 PArpackSolver<VectorType>::internal_reinit(const IndexSet &locally_owned_dofs)
 {
   // store local indices to write to vectors
-  locally_owned_dofs.fill_index_vector(local_indices);
+  local_indices = locally_owned_dofs.get_index_vector();
 
   // scalars
   nloc = locally_owned_dofs.n_elements();
index cf0d1ecd95155d8a08bd47fb8d518014199c35a9..03de54f81c618c21da0d629a3e24e421fd9715c4 100644 (file)
@@ -1352,10 +1352,10 @@ namespace MatrixFreeOperators
           }
 
         // setup edge_constrained indices
-        std::vector<types::global_dof_index> interface_indices;
-        mg_constrained_dofs[j]
-          .get_refinement_edge_indices(level)
-          .fill_index_vector(interface_indices);
+        const std::vector<types::global_dof_index> interface_indices =
+          mg_constrained_dofs[j]
+            .get_refinement_edge_indices(level)
+            .get_index_vector();
         edge_constrained_indices[j].clear();
         edge_constrained_indices[j].reserve(interface_indices.size());
         edge_constrained_values[j].resize(interface_indices.size());
index d6aaab50aca75f79e50c0c24ff6415a12aed4a4c..ef0a7c05e7f6cfdd9e1eff95f2af0413db437976 100644 (file)
@@ -684,7 +684,7 @@ IndexSet::at(const size_type global_index) const
 
 
 
-std::vector<size_type>
+std::vector<IndexSet::size_type>
 IndexSet::get_index_vector() const
 {
   compress();
@@ -762,8 +762,7 @@ IndexSet::make_tpetra_map(const MPI_Comm &communicator,
     );
   else
     {
-      std::vector<size_type> indices;
-      fill_index_vector(indices);
+      const std::vector<size_type>         indices = get_index_vector();
       std::vector<types::global_dof_index> int_indices(indices.size());
       std::copy(indices.begin(), indices.end(), int_indices.begin());
       const Teuchos::ArrayView<types::global_dof_index> arr_view(int_indices);
@@ -831,13 +830,12 @@ IndexSet::make_trilinos_map(const MPI_Comm &communicator,
     );
   else
     {
-      std::vector<size_type> indices;
-      fill_index_vector(indices);
+      const std::vector<size_type> indices = get_index_vector();
       return Epetra_Map(
         TrilinosWrappers::types::int_type(-1),
         TrilinosWrappers::types::int_type(n_elements()),
         (n_elements() > 0 ?
-           reinterpret_cast<TrilinosWrappers::types::int_type *>(
+           reinterpret_cast<const TrilinosWrappers::types::int_type *>(
              indices.data()) :
            nullptr),
         0,
index 693042d549a90d47cc50fed5d852bb6a608f0ed0..c42d77c638d24f22c62a4d39ee9885f702592643 100644 (file)
@@ -365,13 +365,14 @@ namespace Utilities
       // exchange step and see whether the ghost indices sent to us by other
       // processes (ghost_indices) are the same as we hold locally
       // (ghost_indices_ref).
-      std::vector<types::global_dof_index> ghost_indices_ref;
-      ghost_indices_data.fill_index_vector(ghost_indices_ref);
+      const std::vector<types::global_dof_index> ghost_indices_ref =
+        ghost_indices_data.get_index_vector();
       AssertDimension(ghost_indices_ref.size(), n_ghost_indices());
       std::vector<types::global_dof_index> indices_to_send(n_import_indices());
       std::vector<types::global_dof_index> ghost_indices(n_ghost_indices());
-      std::vector<types::global_dof_index> my_indices;
-      locally_owned_range_data.fill_index_vector(my_indices);
+
+      const std::vector<types::global_dof_index> my_indices =
+        locally_owned_range_data.get_index_vector();
       std::vector<MPI_Request> requests;
       n_ghost_indices_in_larger_set = n_ghost_indices_data;
       export_to_ghosted_array_start(127,
index 79a41b2dd827cb4c6abc587e57e1a9e8c9fa15dc..d305051996c903898a754b0d94a28410a4379d78 100644 (file)
@@ -276,8 +276,7 @@ namespace PETScWrappers
       ghosted       = true;
       ghost_indices = ghostnodes;
 
-      std::vector<size_type> ghostindices;
-      ghostnodes.fill_index_vector(ghostindices);
+      const std::vector<size_type> ghostindices = ghostnodes.get_index_vector();
 
       const PetscInt *ptr =
         (ghostindices.size() > 0 ?
index 5cdc79a38a65ef8a907ca3494382336b342ac4a6..5d0ce1392f64babd28c7df884fbe84dc3b40f295 100644 (file)
@@ -648,14 +648,14 @@ namespace TrilinosWrappers
       // the owned rows. In that case, do not create the nonlocal graph and
       // fill the columns by demand
       const bool have_ghost_rows = [&]() {
-        std::vector<dealii::types::global_dof_index> indices;
-        relevant_rows.fill_index_vector(indices);
+        const std::vector<dealii::types::global_dof_index> indices =
+          relevant_rows.get_index_vector();
         Epetra_Map relevant_map(
           TrilinosWrappers::types::int_type(-1),
           TrilinosWrappers::types::int_type(relevant_rows.n_elements()),
           (indices.empty() ?
              nullptr :
-             reinterpret_cast<TrilinosWrappers::types::int_type *>(
+             reinterpret_cast<const TrilinosWrappers::types::int_type *>(
                indices.data())),
           0,
           row_space_map.Comm());

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.