From 16759037033f5de05fd7f739e9d2016e849f1513 Mon Sep 17 00:00:00 2001
From: Bruno Turcksin <bruno.turcksin@gmail.com>
Date: Fri, 20 May 2016 17:02:33 -0400
Subject: [PATCH] Fix a bug where size_type was used instead of std::size_t.

---
 doc/news/changes.h                            |  6 ++++++
 include/deal.II/lac/sparse_matrix.h           | 16 ++++++++--------
 include/deal.II/lac/sparse_matrix.templates.h |  8 ++++----
 include/deal.II/lac/sparsity_pattern.h        |  8 ++++----
 source/lac/sparsity_pattern.cc                |  2 +-
 5 files changed, 23 insertions(+), 17 deletions(-)

diff --git a/doc/news/changes.h b/doc/news/changes.h
index 08159b1464..1caca6a5ff 100644
--- a/doc/news/changes.h
+++ b/doc/news/changes.h
@@ -255,6 +255,12 @@ inconvenience this causes.
 <h3>Specific improvements</h3>
 
 <ol>
+ <li> Fixed: Fix a bug where the SparsityPattern could not have more than 4 
+ billions entries when using 32bit indices.
+ <br>
+ (Bruno Turcksin, 2016/05/22)
+  </li>
+
  <li> New: Added PArpackSolver::reinit() when dealing with BlockVectors.
  <br>
  (Alberto Sartori, 2016/05/19)
diff --git a/include/deal.II/lac/sparse_matrix.h b/include/deal.II/lac/sparse_matrix.h
index 50ea51f48a..21a4a52d5a 100644
--- a/include/deal.II/lac/sparse_matrix.h
+++ b/include/deal.II/lac/sparse_matrix.h
@@ -685,7 +685,7 @@ public:
    * 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 n_nonzero_elements () const;
+  std::size_t n_nonzero_elements () const;
 
   /**
    * Return the number of actually nonzero elements of this matrix. It is
@@ -696,7 +696,7 @@ public:
    * count all entries of the sparsity pattern but only the ones that are
    * nonzero (or whose absolute value is greater than threshold).
    */
-  size_type n_actually_nonzero_elements (const double threshold = 0.) const;
+  std::size_t n_actually_nonzero_elements (const double threshold = 0.) const;
 
   /**
    * Return a (constant) reference to the underlying sparsity pattern of this
@@ -750,8 +750,8 @@ public:
    */
   template <typename number2>
   void set (const std::vector<size_type> &indices,
-            const FullMatrix<number2>       &full_matrix,
-            const bool                       elide_zero_values = false);
+            const FullMatrix<number2>    &full_matrix,
+            const bool                    elide_zero_values = false);
 
   /**
    * Same function as before, but now including the possibility to use
@@ -1045,7 +1045,7 @@ public:
    * @dealiiOperationIsMultithreaded
    */
   template <class OutVector, class InVector>
-  void vmult (OutVector &dst,
+  void vmult (OutVector      &dst,
               const InVector &src) const;
 
   /**
@@ -1064,7 +1064,7 @@ public:
    * Source and destination must not be the same vector.
    */
   template <class OutVector, class InVector>
-  void Tvmult (OutVector &dst,
+  void Tvmult (OutVector      &dst,
                const InVector &src) const;
 
   /**
@@ -1084,7 +1084,7 @@ public:
    * @dealiiOperationIsMultithreaded
    */
   template <class OutVector, class InVector>
-  void vmult_add (OutVector &dst,
+  void vmult_add (OutVector      &dst,
                   const InVector &src) const;
 
   /**
@@ -1103,7 +1103,7 @@ public:
    * Source and destination must not be the same vector.
    */
   template <class OutVector, class InVector>
-  void Tvmult_add (OutVector &dst,
+  void Tvmult_add (OutVector      &dst,
                    const InVector &src) const;
 
   /**
diff --git a/include/deal.II/lac/sparse_matrix.templates.h b/include/deal.II/lac/sparse_matrix.templates.h
index f4c7060b07..abc199329f 100644
--- a/include/deal.II/lac/sparse_matrix.templates.h
+++ b/include/deal.II/lac/sparse_matrix.templates.h
@@ -196,7 +196,7 @@ SparseMatrix<number>::operator = (const double d)
   // operator=. The grain size is chosen to reflect the number of rows in
   // minimum_parallel_grain_size, weighted by the number of nonzero entries
   // per row on average.
-  const size_type matrix_size = cols->n_nonzero_elements();
+  const std::size_t matrix_size = cols->n_nonzero_elements();
   const size_type grain_size =
     internal::SparseMatrix::minimum_parallel_grain_size *
     (cols->n_nonzero_elements()+m()) / m();
@@ -298,7 +298,7 @@ SparseMatrix<number>::get_row_length (const size_type row) const
 
 
 template <typename number>
-typename SparseMatrix<number>::size_type
+std::size_t
 SparseMatrix<number>::n_nonzero_elements () const
 {
   Assert (cols != 0, ExcNotInitialized());
@@ -308,13 +308,13 @@ SparseMatrix<number>::n_nonzero_elements () const
 
 
 template <typename number>
-typename SparseMatrix<number>::size_type
+std::size_t
 SparseMatrix<number>::n_actually_nonzero_elements (const double threshold) const
 {
   Assert (cols != 0, ExcNotInitialized());
   Assert (threshold >= 0, ExcMessage ("Negative threshold!"));
   size_type nnz = 0;
-  const size_type nnz_alloc = n_nonzero_elements();
+  const std::size_t nnz_alloc = n_nonzero_elements();
   for (size_type i=0; i<nnz_alloc; ++i)
     if (std::abs(val[i]) > threshold)
       ++nnz;
diff --git a/include/deal.II/lac/sparsity_pattern.h b/include/deal.II/lac/sparsity_pattern.h
index 3403d7652d..a5645e05a2 100644
--- a/include/deal.II/lac/sparsity_pattern.h
+++ b/include/deal.II/lac/sparsity_pattern.h
@@ -762,7 +762,7 @@ public:
    * This function may only be called if the matrix struct is compressed. It
    * does not make too much sense otherwise anyway.
    */
-  size_type n_nonzero_elements () const;
+  std::size_t n_nonzero_elements () const;
 
   /**
    * Return whether the structure is compressed or not.
@@ -848,7 +848,7 @@ public:
    * of this function is <i>log(N)</i>.
    */
   std::pair<size_type, size_type>
-  matrix_position (const size_type global_index) const;
+  matrix_position (const std::size_t global_index) const;
 
   /**
    * Check if a value at a certain position may be non-zero.
@@ -1035,7 +1035,7 @@ private:
    * for the #rowstart array, i.e. it may be larger than the actually used
    * part of the array.
    */
-  size_type max_vec_len;
+  std::size_t max_vec_len;
 
   /**
    * Maximum number of elements per row. This is set to the value given to the
@@ -1400,7 +1400,7 @@ SparsityPattern::column_number (const size_type row,
 
 
 inline
-SparsityPattern::size_type
+std::size_t
 SparsityPattern::n_nonzero_elements () const
 {
   Assert ((rowstart!=0) && (colnums!=0), ExcEmptyObject());
diff --git a/source/lac/sparsity_pattern.cc b/source/lac/sparsity_pattern.cc
index 070a846457..79dbf8b04d 100644
--- a/source/lac/sparsity_pattern.cc
+++ b/source/lac/sparsity_pattern.cc
@@ -760,7 +760,7 @@ SparsityPattern::row_position (const size_type i, const size_type j) const
 
 
 std::pair<SparsityPattern::size_type, SparsityPattern::size_type>
-SparsityPattern::matrix_position (const size_type global_index) const
+SparsityPattern::matrix_position (const std::size_t global_index) const
 {
   Assert (compressed == true, ExcNotCompressed());
   Assert (global_index < n_nonzero_elements(),
-- 
2.39.5