From e86613eafa016307c261974f771db9bb378a8c70 Mon Sep 17 00:00:00 2001
From: Marc Fehling <mafehling.git@gmail.com>
Date: Mon, 4 Apr 2022 20:26:49 -0600
Subject: [PATCH] Use `uint64` for `n_nonzero_elements`.

done: PETSc MatrixBase, Trilinos SparseMatrix, PETSc&Trilinos BlockSparseMatrix
todo: dealii SparseMatrix, all SparsityPattern objects

Should be backwards compatible with implicit type conversion.
---
 .../changes/incompatibilities/20220404Fehling |  5 +++++
 include/deal.II/base/types.h                  |  7 +++++-
 .../deal.II/lac/petsc_block_sparse_matrix.h   |  2 +-
 include/deal.II/lac/petsc_matrix_base.h       |  2 +-
 .../lac/trilinos_block_sparse_matrix.h        |  2 +-
 include/deal.II/lac/trilinos_index_access.h   | 22 ++++++-------------
 include/deal.II/lac/trilinos_sparse_matrix.h  | 14 ++++++------
 .../deal.II/lac/trilinos_sparsity_pattern.h   |  2 +-
 source/lac/petsc_matrix_base.cc               |  6 +++--
 .../lac/petsc_parallel_block_sparse_matrix.cc |  4 ++--
 source/lac/trilinos_block_sparse_matrix.cc    |  4 ++--
 source/lac/trilinos_sparsity_pattern.cc       |  6 ++---
 12 files changed, 40 insertions(+), 36 deletions(-)
 create mode 100644 doc/news/changes/incompatibilities/20220404Fehling

diff --git a/doc/news/changes/incompatibilities/20220404Fehling b/doc/news/changes/incompatibilities/20220404Fehling
new file mode 100644
index 0000000000..3dcfb29ad4
--- /dev/null
+++ b/doc/news/changes/incompatibilities/20220404Fehling
@@ -0,0 +1,5 @@
+Changed: All functions `n_nonzero_elements()`, which are members of the
+SparseMatrix and SparsityPattern family of classes, now always return
+`std::uint64_t` instead of types::global_dof_index.
+<br>
+(Marc Fehling, 2022/04/04)
diff --git a/include/deal.II/base/types.h b/include/deal.II/base/types.h
index ac262ad545..57d5a1a42d 100644
--- a/include/deal.II/base/types.h
+++ b/include/deal.II/base/types.h
@@ -167,11 +167,16 @@ namespace TrilinosWrappers
 {
   namespace types
   {
+    /**
+     * Declare type of 64 bit integer used in the Epetra package of Trilinos.
+     */
+    using int64_type = long long int;
+
 #ifdef DEAL_II_WITH_64BIT_INDICES
     /**
      * Declare type of integer used in the Epetra package of Trilinos.
      */
-    using int_type = long long int;
+    using int_type = int64_type;
 #else
     /**
      * Declare type of integer used in the Epetra package of Trilinos.
diff --git a/include/deal.II/lac/petsc_block_sparse_matrix.h b/include/deal.II/lac/petsc_block_sparse_matrix.h
index 60ac2c3e93..87b2499f42 100644
--- a/include/deal.II/lac/petsc_block_sparse_matrix.h
+++ b/include/deal.II/lac/petsc_block_sparse_matrix.h
@@ -255,7 +255,7 @@ namespace PETScWrappers
        * returns the number of entries in the sparsity pattern; if any of the
        * entries should happen to be zero, it is counted anyway.
        */
-      size_type
+      std::uint64_t
       n_nonzero_elements() const;
 
       /**
diff --git a/include/deal.II/lac/petsc_matrix_base.h b/include/deal.II/lac/petsc_matrix_base.h
index f0e94b0531..ac6c547374 100644
--- a/include/deal.II/lac/petsc_matrix_base.h
+++ b/include/deal.II/lac/petsc_matrix_base.h
@@ -653,7 +653,7 @@ namespace PETScWrappers
      * returns the number of entries in the sparsity pattern; if any of the
      * entries should happen to be zero, it is counted anyway.
      */
-    size_type
+    std::uint64_t
     n_nonzero_elements() const;
 
     /**
diff --git a/include/deal.II/lac/trilinos_block_sparse_matrix.h b/include/deal.II/lac/trilinos_block_sparse_matrix.h
index 4a165f5f69..440fe662ee 100644
--- a/include/deal.II/lac/trilinos_block_sparse_matrix.h
+++ b/include/deal.II/lac/trilinos_block_sparse_matrix.h
@@ -218,7 +218,7 @@ namespace TrilinosWrappers
      * Return the total number of nonzero elements of this matrix (summed
      * over all MPI processes).
      */
-    size_type
+    std::uint64_t
     n_nonzero_elements() const;
 
     /**
diff --git a/include/deal.II/lac/trilinos_index_access.h b/include/deal.II/lac/trilinos_index_access.h
index 73396a86bc..c6467298b5 100644
--- a/include/deal.II/lac/trilinos_index_access.h
+++ b/include/deal.II/lac/trilinos_index_access.h
@@ -32,18 +32,14 @@ DEAL_II_NAMESPACE_OPEN
 namespace TrilinosWrappers
 {
   /**
-   * A helper function that queries the size of an Epetra_BlockMap object
-   * and calls either the 32 or 64 bit function to get the number of global
-   * elements in the map.
+   * A helper function that queries the size of an Epetra_BlockMap object. We
+   * always call the 64 bit function to get the number of global elements in the
+   * map.
    */
-  inline TrilinosWrappers::types::int_type
+  inline TrilinosWrappers::types::int64_type
   n_global_elements(const Epetra_BlockMap &map)
   {
-#  ifdef DEAL_II_WITH_64BIT_INDICES
     return map.NumGlobalElements64();
-#  else
-    return map.NumGlobalElements();
-#  endif
   }
 
   /**
@@ -133,17 +129,13 @@ namespace TrilinosWrappers
   }
 
   /**
-   * A helper function that finds the number of global entries by calling
-   * either the 32 or 64 bit function.
+   * A helper function that finds the number of global entries. We always call
+   * the 64 bit function.
    */
-  inline TrilinosWrappers::types::int_type
+  inline TrilinosWrappers::types::int64_type
   n_global_entries(const Epetra_CrsGraph &graph)
   {
-#  ifdef DEAL_II_WITH_64BIT_INDICES
     return graph.NumGlobalEntries64();
-#  else
-    return graph.NumGlobalEntries();
-#  endif
   }
 
   /**
diff --git a/include/deal.II/lac/trilinos_sparse_matrix.h b/include/deal.II/lac/trilinos_sparse_matrix.h
index 8eec00f32d..88d6f536a6 100644
--- a/include/deal.II/lac/trilinos_sparse_matrix.h
+++ b/include/deal.II/lac/trilinos_sparse_matrix.h
@@ -947,7 +947,7 @@ namespace TrilinosWrappers
      * Return the total number of nonzero elements of this matrix (summed
      * over all MPI processes).
      */
-    size_type
+    std::uint64_t
     n_nonzero_elements() const;
 
     /**
@@ -3015,14 +3015,14 @@ namespace TrilinosWrappers
 
 
 
-  inline SparseMatrix::size_type
+  inline std::uint64_t
   SparseMatrix::n_nonzero_elements() const
   {
-#      ifndef DEAL_II_WITH_64BIT_INDICES
-    return matrix->NumGlobalNonzeros();
-#      else
-    return matrix->NumGlobalNonzeros64();
-#      endif
+    // Trilinos uses 64bit functions internally for attribute access, which
+    // return `long long`. They also offer 32bit variants that return `int`,
+    // however those call the 64bit version and convert the values to 32bit.
+    // There is no necessity in using the 32bit versions at all.
+    return static_cast<std::uint64_t>(matrix->NumGlobalNonzeros64());
   }
 
 
diff --git a/include/deal.II/lac/trilinos_sparsity_pattern.h b/include/deal.II/lac/trilinos_sparsity_pattern.h
index eebc438d45..5598833b23 100644
--- a/include/deal.II/lac/trilinos_sparsity_pattern.h
+++ b/include/deal.II/lac/trilinos_sparsity_pattern.h
@@ -712,7 +712,7 @@ namespace TrilinosWrappers
     /**
      * Return the number of nonzero elements of this sparsity pattern.
      */
-    size_type
+    std::uint64_t
     n_nonzero_elements() const;
 
     /**
diff --git a/source/lac/petsc_matrix_base.cc b/source/lac/petsc_matrix_base.cc
index 303b0ac0c2..aa05073f28 100644
--- a/source/lac/petsc_matrix_base.cc
+++ b/source/lac/petsc_matrix_base.cc
@@ -286,14 +286,16 @@ namespace PETScWrappers
 
 
 
-  MatrixBase::size_type
+  std::uint64_t
   MatrixBase::n_nonzero_elements() const
   {
     MatInfo              mat_info;
     const PetscErrorCode ierr = MatGetInfo(matrix, MAT_GLOBAL_SUM, &mat_info);
     AssertThrow(ierr == 0, ExcPETScError(ierr));
 
-    return static_cast<size_type>(mat_info.nz_used);
+    // MatInfo logs quantities as PetscLogDouble. So we need to cast it to match
+    // out interface.
+    return static_cast<std::uint64_t>(mat_info.nz_used);
   }
 
 
diff --git a/source/lac/petsc_parallel_block_sparse_matrix.cc b/source/lac/petsc_parallel_block_sparse_matrix.cc
index 27ca8cb3be..9a35d250a4 100644
--- a/source/lac/petsc_parallel_block_sparse_matrix.cc
+++ b/source/lac/petsc_parallel_block_sparse_matrix.cc
@@ -135,10 +135,10 @@ namespace PETScWrappers
       return index_sets;
     }
 
-    BlockSparseMatrix::size_type
+    std::uint64_t
     BlockSparseMatrix::n_nonzero_elements() const
     {
-      size_type n_nonzero = 0;
+      std::uint64_t n_nonzero = 0;
       for (size_type rows = 0; rows < this->n_block_rows(); ++rows)
         for (size_type cols = 0; cols < this->n_block_cols(); ++cols)
           n_nonzero += this->block(rows, cols).n_nonzero_elements();
diff --git a/source/lac/trilinos_block_sparse_matrix.cc b/source/lac/trilinos_block_sparse_matrix.cc
index bd6b591eb9..586b9bb8c6 100644
--- a/source/lac/trilinos_block_sparse_matrix.cc
+++ b/source/lac/trilinos_block_sparse_matrix.cc
@@ -232,10 +232,10 @@ namespace TrilinosWrappers
 
 
 
-  BlockSparseMatrix::size_type
+  std::uint64_t
   BlockSparseMatrix::n_nonzero_elements() const
   {
-    size_type n_nonzero = 0;
+    std::uint64_t n_nonzero = 0;
     for (size_type rows = 0; rows < this->n_block_rows(); ++rows)
       for (size_type cols = 0; cols < this->n_block_cols(); ++cols)
         n_nonzero += this->block(rows, cols).n_nonzero_elements();
diff --git a/source/lac/trilinos_sparsity_pattern.cc b/source/lac/trilinos_sparsity_pattern.cc
index 19710f75ad..7109b5a8e5 100644
--- a/source/lac/trilinos_sparsity_pattern.cc
+++ b/source/lac/trilinos_sparsity_pattern.cc
@@ -932,12 +932,12 @@ namespace TrilinosWrappers
 
 
 
-  SparsityPattern::size_type
+  std::uint64_t
   SparsityPattern::n_nonzero_elements() const
   {
-    TrilinosWrappers::types::int_type nnz = n_global_entries(*graph);
+    TrilinosWrappers::types::int64_type nnz = n_global_entries(*graph);
 
-    return static_cast<size_type>(nnz);
+    return static_cast<std::uint64_t>(nnz);
   }
 
 
-- 
2.39.5