bug #449: add SparseMatrix::conservativeResize feature

Eigen/src/SparseCore/SparseMatrix.h

@@ -286,7 +286,8 @@ class SparseMatrix
       {
         std::size_t totalReserveSize = 0;
         // turn the matrix into non-compressed mode
-        m_innerNonZeros = new Index[m_outerSize];
+        m_innerNonZeros = static_cast<Index*>(std::malloc(m_outerSize * sizeof(Index)));
+        if (!m_innerNonZeros) internal::throw_std_bad_alloc();
         
         // temporarily use m_innerSizes to hold the new starting points.
         Index* newOuterIndex = m_innerNonZeros;
@@ -318,7 +319,9 @@ class SparseMatrix
       }
       else
       {
-        Index* newOuterIndex = new Index[m_outerSize+1];
+        Index* newOuterIndex = static_cast<Index*>(std::malloc((m_outerSize+1)*sizeof(Index)));
+        if (!newOuterIndex) internal::throw_std_bad_alloc();
+        
         Index count = 0;
         for(Index j=0; j<m_outerSize; ++j)
         {
@@ -345,7 +348,7 @@ class SparseMatrix
         }
         
         std::swap(m_outerIndex, newOuterIndex);
-        delete[] newOuterIndex;
+        std::free(newOuterIndex);
       }
       
     }
@@ -460,7 +463,7 @@ class SparseMatrix
         m_outerIndex[j+1] = m_outerIndex[j] + m_innerNonZeros[j];
         oldStart = nextOldStart;
       }
-      delete[] m_innerNonZeros;
+      std::free(m_innerNonZeros);
       m_innerNonZeros = 0;
       m_data.resize(m_outerIndex[m_outerSize]);
       m_data.squeeze();
@@ -506,6 +509,64 @@ class SparseMatrix
       m_data.resize(k,0);
     }
 
+    /** Resizes the matrix to a \a rows x \a cols matrix leaving old values untouched.
+      * \sa resizeNonZeros(Index), reserve(), setZero()
+      */
+    void conservativeResize(Index rows, Index cols) 
+    {
+        // No change
+        if (this->rows() == rows && this->cols() == cols) return;
+
+        Index innerChange = IsRowMajor ? cols - this->cols() : rows - this->rows();
+        Index outerChange = IsRowMajor ? rows - this->rows() : cols - this->cols();
+        Index newInnerSize = IsRowMajor ? cols : rows;
+
+        // Deals with inner non zeros
+        if (m_innerNonZeros)
+        {
+          // Resize m_innerNonZeros
+          Index *newInnerNonZeros = static_cast<Index*>(std::realloc(m_innerNonZeros, (m_outerSize + outerChange) * sizeof(Index)));
+          if (!newInnerNonZeros) internal::throw_std_bad_alloc();
+          m_innerNonZeros = newInnerNonZeros;
+          
+          for(Index i=m_outerSize; i<m_outerSize+outerChange; i++)          
+            m_innerNonZeros[i] = 0;
+        } 
+        else if (innerChange < 0) 
+        {
+          // Inner size decreased: allocate a new m_innerNonZeros
+          m_innerNonZeros = static_cast<Index*>(std::malloc((m_outerSize+outerChange+1) * sizeof(Index)));
+          if (!m_innerNonZeros) internal::throw_std_bad_alloc();
+          for(Index i = 0; i < m_outerSize; i++)
+            m_innerNonZeros[i] = m_outerIndex[i+1] - m_outerIndex[i];
+        }
+        
+        // Change the m_innerNonZeros in case of a decrease of inner size
+        if (m_innerNonZeros && innerChange < 0) {
+              for(Index i = 0; i<m_outerSize; i++)
+              {
+                Index &n = m_innerNonZeros[i];
+                Index start = m_outerIndex[i];
+                while (n >= 0 && m_data.index(start+n-1) >= newInnerSize) --n; 
+              }
+        }
+        
+        m_innerSize = newInnerSize;
+
+        // Re-allocate outer index structure if necessary
+        if (outerChange == 0)
+          return;
+            
+        Index *newOuterIndex = static_cast<Index*>(std::realloc(m_outerIndex, (m_outerSize + outerChange + 1) * sizeof(Index)));
+        if (!newOuterIndex) internal::throw_std_bad_alloc();
+        m_outerIndex = newOuterIndex;
+        Index last = m_outerSize == 0 ? 0 : m_outerIndex[m_outerSize];
+        for(Index i=m_outerSize; i<m_outerSize+outerChange+1; i++)          
+            m_outerIndex[i] = last; 
+        m_outerSize += outerChange;
+        
+    }
+    
     /** Resizes the matrix to a \a rows x \a cols matrix and initializes it to zero.
       * \sa resizeNonZeros(Index), reserve(), setZero()
       */
@@ -516,13 +577,15 @@ class SparseMatrix
       m_data.clear();
       if (m_outerSize != outerSize || m_outerSize==0)
       {
-        delete[] m_outerIndex;
+        std::free(m_outerIndex);
-        m_outerIndex = new Index [outerSize+1];
+        m_outerIndex = static_cast<Index*>(std::malloc((outerSize + 1) * sizeof(Index)));
+        if (!m_outerIndex) internal::throw_std_bad_alloc();
+        
         m_outerSize = outerSize;
       }
       if(m_innerNonZeros)
       {
-        delete[] m_innerNonZeros;
+        std::free(m_innerNonZeros);
         m_innerNonZeros = 0;
       }
       memset(m_outerIndex, 0, (m_outerSize+1)*sizeof(Index));
@@ -714,8 +777,8 @@ class SparseMatrix
     /** Destructor */
     inline ~SparseMatrix()
     {
-      delete[] m_outerIndex;
+      std::free(m_outerIndex);
-      delete[] m_innerNonZeros;
+      std::free(m_innerNonZeros);
     }
 
 #ifndef EIGEN_PARSED_BY_DOXYGEN
@@ -735,7 +798,7 @@ protected:
       resize(other.rows(), other.cols());
       if(m_innerNonZeros)
       {
-        delete[] m_innerNonZeros;
+        std::free(m_innerNonZeros);
         m_innerNonZeros = 0;
       }
     }
@@ -1105,7 +1168,7 @@ void SparseMatrix<Scalar,_Options,_Index>::sumupDuplicates()
   m_outerIndex[m_outerSize] = count;
 
   // turn the matrix into compressed form
-  delete[] m_innerNonZeros;
+  std::free(m_innerNonZeros);
   m_innerNonZeros = 0;
   m_data.resize(m_outerIndex[m_outerSize]);
 }

test/sparse_basic.cpp

@@ -379,6 +379,41 @@ template<typename SparseMatrixType> void sparse_basic(const SparseMatrixType& re
     initSparse<Scalar>(density, refMat2, m2);
     VERIFY_IS_APPROX(m2.diagonal(), refMat2.diagonal().eval());
   }
+  
+  // test conservative resize
+  {
+      std::vector< std::pair<int,int> > inc;
+      inc.push_back(std::pair<int,int>(-3,-2));
+      inc.push_back(std::pair<int,int>(0,0));
+      inc.push_back(std::pair<int,int>(3,2));
+      inc.push_back(std::pair<int,int>(3,0));
+      inc.push_back(std::pair<int,int>(0,3));
+      
+      for(size_t i = 0; i< inc.size(); i++) {
+        int incRows = inc[i].first;
+        int incCols = inc[i].second;
+        SparseMatrixType m1(rows, cols);
+        DenseMatrix refMat1 = DenseMatrix::Zero(rows, cols);
+        initSparse<Scalar>(density, refMat1, m1);
+        
+        m1.conservativeResize(rows+incRows, cols+incCols);
+        refMat1.conservativeResize(rows+incRows, cols+incCols);
+        if (incRows > 0) refMat1.bottomRows(incRows).setZero();
+        if (incCols > 0) refMat1.rightCols(incCols).setZero();
+        
+        VERIFY_IS_APPROX(m1, refMat1);
+        
+        // Insert new values
+        if (incRows > 0) 
+          m1.insert(refMat1.rows()-1, 0) = refMat1(refMat1.rows()-1, 0) = 1;
+        if (incCols > 0) 
+          m1.insert(0, refMat1.cols()-1) = refMat1(0, refMat1.cols()-1) = 1;
+          
+        VERIFY_IS_APPROX(m1, refMat1);
+          
+          
+      }
+  }
 }
 
 void test_sparse_basic()