Fix for non-symmetric matrix diagonal fill (Closes #619)

Author: coatless

ChangeLog

@@ -1,3 +1,11 @@
+2016-12-31  James J Balamuta  <[email protected]>
+
+        * inst/include/Rcpp/vector/Matrix.h: Fixed non-symmetric case of matrix
+        fills by switching to a loop based solution from iterator.
+        * inst/unitTests/runit.Matrix.R: Added unit tests for diagonally
+        filling matrices.
+        * inst/unitTests/cpp/Matrix.cpp: Idem
+
 2016-12-26  Dirk Eddelbuettel  <[email protected]>
 
         * R/Attributes.R: Added #nocov markers

inst/NEWS.Rd

@@ -11,6 +11,8 @@
       compiler versions (Jim Hester in \ghpr{598})
       \item Date and Datetime object and vector now have format methods and
       \code{operator<<} support (PR\ghpr{599})
+      \item Addressed improper diagonal fill for non-symmetric matrices
+      (James Balamuta in \ghpr{622} addressing \ghit{619})
     }
     \item Changes in Rcpp Sugar:
     \itemize{
@@ -21,7 +23,9 @@
     \itemize{
       \item Added Environment::find unit tests and an Environment::get(Symbol)
       test (James Balamuta in \ghpr{595} addressing issue \ghit{594}).
-    }	
+      \item Added diagonal matrix fill tests
+      (James Balamuta in \ghpr{622} addressing \ghit{619})
+    }
     \item Changes in Rcpp Documentation:
     \itemize{
       \item Exposed pointers macros were included in the Rcpp Extending vignette

inst/include/Rcpp/vector/Matrix.h

@@ -175,25 +175,25 @@ class Matrix : public Vector<RTYPE, StoragePolicy>, public MatrixBase<RTYPE, tru
 
     template <typename U>
     void fill_diag__dispatch( traits::false_type, const U& u) {
-          Shield<SEXP> elem( converter_type::get( u ) ) ;
-        int n = Matrix::ncol() ;
-        int offset = n +1 ;
-        iterator it( VECTOR::begin()) ;
-        for( int i=0; i<n; i++){
-            *it = ::Rf_duplicate( elem );
-            it += offset;
+        Shield<SEXP> elem( converter_type::get( u ) );
+
+        R_xlen_t bounds = ( Matrix::nrow() < Matrix::ncol() ) ?
+                             Matrix::nrow() : Matrix::ncol();
+        
+        for (R_xlen_t i = 0; i < bounds; ++i) {
+            (*this)(i, i) = elem;
         }
     }
 
     template <typename U>
     void fill_diag__dispatch( traits::true_type, const U& u) {
-          stored_type elem = converter_type::get( u ) ;
-        int n = Matrix::ncol() ;
-        int offset = n + 1 ;
-        iterator it( VECTOR::begin()) ;
-        for( int i=0; i<n; i++){
-            *it = elem ;
-            it += offset;
+        stored_type elem = converter_type::get( u );
+		
+        R_xlen_t bounds = ( Matrix::nrow() < Matrix::ncol() ) ?
+		                 Matrix::nrow() : Matrix::ncol();
+
+        for (R_xlen_t i = 0; i < bounds; ++i) {
+			(*this)(i, i) = elem;
         }
     }
 

inst/unitTests/cpp/Matrix.cpp

@@ -397,3 +397,16 @@ Rcpp::LogicalMatrix lgl_zeros(int n) {
     return Rcpp::LogicalMatrix::zeros(n);
 }
 
+// --- Diagonal Fill
+
+// [[Rcpp::export]]
+Rcpp::NumericMatrix num_diag_fill(Rcpp::NumericMatrix x, double diag_val) {
+    x.fill_diag(diag_val);
+    return x;
+}
+
+// [[Rcpp::export]]
+Rcpp::CharacterMatrix char_diag_fill(Rcpp::CharacterMatrix x, std::string diag_val) {
+    x.fill_diag(diag_val);
+    return x;
+}

inst/unitTests/runit.Matrix.R

@@ -338,5 +338,30 @@ if (.runThisTest) {
             "zeros - logical"
         )
     }
+    
+    
+    test.Matrix.diagfill <- function() {
+        
+        checkEquals(num_diag_fill(diag(1.0, 2, 4), 0.0),
+                    matrix(0.0, 2, 4),
+                    msg = "diagonal fill - case: n < p")
+        
+        checkEquals(num_diag_fill(diag(1.0, 4, 2), 0.0),
+                    matrix(0.0, 4, 2),
+                    msg = "diagonal fill - case: n > p")
+        
+        checkEquals(num_diag_fill(diag(1.0, 3, 3), 0.0),
+                    matrix(0.0, 3, 3),
+                    msg = "diagonal fill - case: n = p")
+        
+        m <- matrix("", 2, 4)
+        diag(m) <- letters[1:2]
+        
+        checkEquals(char_diag_fill(m, ""), 
+                    matrix("", 2, 4),
+                    msg = "diagonal fill - char")
+        
+        
+    }
 
 }