[Rcpp-commits] r3156 - in pkg/RcppEigen: . R inst/include inst/include/unsupported/Eigen inst/include/unsupported/Eigen/src/SparseExtra inst/skeleton inst/unitTests

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Thu Jul 28 18:01:20 CEST 2011


Author: dmbates
Date: 2011-07-28 18:01:20 +0200 (Thu, 28 Jul 2011)
New Revision: 3156

Modified:
   pkg/RcppEigen/DESCRIPTION
   pkg/RcppEigen/NAMESPACE
   pkg/RcppEigen/R/SHLIB.R
   pkg/RcppEigen/R/inline.R
   pkg/RcppEigen/inst/include/RcppEigen.h
   pkg/RcppEigen/inst/include/RcppEigenConfig.h
   pkg/RcppEigen/inst/include/RcppEigenForward.h
   pkg/RcppEigen/inst/include/RcppEigenWrap.h
   pkg/RcppEigen/inst/include/unsupported/Eigen/SparseExtra
   pkg/RcppEigen/inst/include/unsupported/Eigen/src/SparseExtra/CholmodSupport.h
   pkg/RcppEigen/inst/skeleton/Makevars
   pkg/RcppEigen/inst/unitTests/runit.sparse.R
Log:
Add support for Cholmod functions linked through the Matrix package.  Tests for same.  Wrap of an Eigen::SparseMatrix now returns an S4 dgCMatrix object.


Modified: pkg/RcppEigen/DESCRIPTION
===================================================================
--- pkg/RcppEigen/DESCRIPTION	2011-07-28 15:58:15 UTC (rev 3155)
+++ pkg/RcppEigen/DESCRIPTION	2011-07-28 16:01:20 UTC (rev 3156)
@@ -2,7 +2,7 @@
 Type: Package
 Title: Rcpp integration for the Eigen templated linear algebra library.
 Version: 0.1.2
-Date: 2011-07-xx
+Date: 2011-07-28
 Author: Douglas Bates, Romain Francois and Dirk Eddelbuettel
 Maintainer: Douglas Bates, Romain Francois and Dirk Eddelbuettel <RcppArmadillo-authors at r-enthusiasts.com>
 Description: R and Eigen integration using Rcpp.
@@ -24,5 +24,6 @@
 Depends: Rcpp (>= 0.9.5.1), methods
 LazyLoad: yes
 LinkingTo: Rcpp
+Imports: Matrix
 Suggests: inline, RUnit
 URL: http://eigen.tuxfamily.org

Modified: pkg/RcppEigen/NAMESPACE
===================================================================
--- pkg/RcppEigen/NAMESPACE	2011-07-28 15:58:15 UTC (rev 3155)
+++ pkg/RcppEigen/NAMESPACE	2011-07-28 16:01:20 UTC (rev 3156)
@@ -1,6 +1,14 @@
 useDynLib(RcppEigen)
 #import(Rcpp)
 
+importClassesFrom("Matrix"
+                  , ddiMatrix
+                  , dgCMatrix
+                  , dgeMatrix
+                  , dsCMatrix
+                  , dtCMatrix
+                  )
+
 #exportPattern("^[[:alpha:]]+")
 export(fastLm,
        fastLmPure,

Modified: pkg/RcppEigen/R/SHLIB.R
===================================================================
--- pkg/RcppEigen/R/SHLIB.R	2011-07-28 15:58:15 UTC (rev 3155)
+++ pkg/RcppEigen/R/SHLIB.R	2011-07-28 16:01:20 UTC (rev 3156)
@@ -18,6 +18,5 @@
 SHLIB <- Rcpp:::SHLIB.maker( 
     env = list( 
         PKG_LIBS = Rcpp:::RcppLdFlags(),
-        PKG_CPPFLAGS = Rcpp:::RcppCxxFlags()
    )
 )

Modified: pkg/RcppEigen/R/inline.R
===================================================================
--- pkg/RcppEigen/R/inline.R	2011-07-28 15:58:15 UTC (rev 3155)
+++ pkg/RcppEigen/R/inline.R	2011-07-28 16:01:20 UTC (rev 3156)
@@ -15,8 +15,10 @@
 ## You should have received a copy of the GNU General Public License
 ## along with RcppEigen.  If not, see <http://www.gnu.org/licenses/>.
 
-inlineCxxPlugin <- Rcpp:::Rcpp.plugin.maker(
-	include.before = "#include <RcppEigen.h>", 
-	package        = "RcppEigen"
-)
+inlineCxxPlugin <-
+    Rcpp:::Rcpp.plugin.maker(
+                             include.before = "#include <RcppEigen.h>", 
+                             package        = "RcppEigen"
+#                             , LinkingTo      = c("RcppEigen", "Rcpp")
+                             )
 

Modified: pkg/RcppEigen/inst/include/RcppEigen.h
===================================================================
--- pkg/RcppEigen/inst/include/RcppEigen.h	2011-07-28 15:58:15 UTC (rev 3155)
+++ pkg/RcppEigen/inst/include/RcppEigen.h	2011-07-28 16:01:20 UTC (rev 3156)
@@ -22,7 +22,6 @@
 #ifndef RcppEigen__RcppEigen__h
 #define RcppEigen__RcppEigen__h
 
-#define EIGEN_YES_I_KNOW_SPARSE_MODULE_IS_NOT_STABLE_YET
 #include <RcppEigenForward.h>
 #include <Rcpp.h>
 #include <RcppEigenWrap.h>

Modified: pkg/RcppEigen/inst/include/RcppEigenConfig.h
===================================================================
--- pkg/RcppEigen/inst/include/RcppEigenConfig.h	2011-07-28 15:58:15 UTC (rev 3155)
+++ pkg/RcppEigen/inst/include/RcppEigenConfig.h	2011-07-28 16:01:20 UTC (rev 3156)
@@ -22,5 +22,7 @@
 #ifndef RcppEigen__RcppEigenConfig__h
 #define RcppEigen__RcppEigenConfig__h
 
+#define EIGEN_YES_I_KNOW_SPARSE_MODULE_IS_NOT_STABLE_YET
+
 #endif
 

Modified: pkg/RcppEigen/inst/include/RcppEigenForward.h
===================================================================
--- pkg/RcppEigen/inst/include/RcppEigenForward.h	2011-07-28 15:58:15 UTC (rev 3155)
+++ pkg/RcppEigen/inst/include/RcppEigenForward.h	2011-07-28 16:01:20 UTC (rev 3156)
@@ -25,6 +25,7 @@
 #include <RcppCommon.h>
 #include <Rconfig.h>
 #include <RcppEigenConfig.h>
+#include <RcppEigenCholmod.h>
 #define EIGEN_PLAINOBJECTBASE_PLUGIN "PlainObjectBaseAddon.h"
 #include <Eigen/Dense>
 #include <unsupported/Eigen/SparseExtra> // also includes Eigen/Sparse

Modified: pkg/RcppEigen/inst/include/RcppEigenWrap.h
===================================================================
--- pkg/RcppEigen/inst/include/RcppEigenWrap.h	2011-07-28 15:58:15 UTC (rev 3155)
+++ pkg/RcppEigen/inst/include/RcppEigenWrap.h	2011-07-28 16:01:20 UTC (rev 3156)
@@ -50,21 +50,21 @@
         template <typename T> 
         SEXP eigen_wrap_plain_dense( const T& object, Rcpp::traits::false_type ){
             typedef typename T::Scalar Scalar ;
-            const int RTYPE = Rcpp::traits::r_sexptype_traits<Scalar>::rtype  ;  
+            const int  RTYPE = Rcpp::traits::r_sexptype_traits<Scalar>::rtype  ;  
             int          nnz = object.nonZeros(), p = object.outerSize();
 	        Dimension    dim(object.innerSize(), p);
 	        const int    *ip = object._innerIndexPtr(), *pp = object._outerIndexPtr();
-	        const Scalar      *xp = object._valuePtr();
+	        const Scalar *xp = object._valuePtr();
 	        IntegerVector iv(ip, ip + nnz), pv(pp, pp + p + 1);
 	        Vector<RTYPE> xv(xp, xp + nnz);
-	        
-	        return List::create(_["Dim"] = dim,
-	        								 _["i"]   = iv,
-	        								 _["p"]   = pv,
-	        								 _["x"]   = xv);
+	        S4           ans("dgCMatrix");
+			ans.slot("Dim")  = dim;
+			ans.slot("i")    = iv;
+			ans.slot("p")    = pv;
+			ans.slot("x")    = xv;
+			return  ans;
 	    } 
         
-        
         // plain object, so we can assume data() and size()
         template <typename T>
         inline SEXP eigen_wrap_is_plain( const T& obj, ::Rcpp::traits::true_type ){

Modified: pkg/RcppEigen/inst/include/unsupported/Eigen/SparseExtra
===================================================================
--- pkg/RcppEigen/inst/include/unsupported/Eigen/SparseExtra	2011-07-28 15:58:15 UTC (rev 3155)
+++ pkg/RcppEigen/inst/include/unsupported/Eigen/SparseExtra	2011-07-28 16:01:20 UTC (rev 3156)
@@ -57,7 +57,7 @@
 #include "src/SparseExtra/Solve.h"
 #include "src/SparseExtra/Amd.h"
 #include "src/SparseExtra/SimplicialCholesky.h"
-
+#include "src/SparseExtra/CholmodSupport.h"
 #include "src/SparseExtra/SparseLLT.h"
 #include "src/SparseExtra/SparseLDLTLegacy.h"
 #include "src/SparseExtra/SparseLU.h"

Modified: pkg/RcppEigen/inst/include/unsupported/Eigen/src/SparseExtra/CholmodSupport.h
===================================================================
--- pkg/RcppEigen/inst/include/unsupported/Eigen/src/SparseExtra/CholmodSupport.h	2011-07-28 15:58:15 UTC (rev 3155)
+++ pkg/RcppEigen/inst/include/unsupported/Eigen/src/SparseExtra/CholmodSupport.h	2011-07-28 16:01:20 UTC (rev 3156)
@@ -176,22 +176,22 @@
     CholmodDecomposition()
       : m_cholmodFactor(0), m_info(Success), m_isInitialized(false)
     {
-      cholmod_start(&m_cholmod);
+      M_R_cholmod_start(&m_cholmod);
       setMode(CholmodLDLt);
     }
 
     CholmodDecomposition(const MatrixType& matrix)
       : m_cholmodFactor(0), m_info(Success), m_isInitialized(false)
     {
-      cholmod_start(&m_cholmod);
+      M_R_cholmod_start(&m_cholmod);
       compute(matrix);
     }
 
     ~CholmodDecomposition()
     {
       if(m_cholmodFactor)
-        cholmod_free_factor(&m_cholmodFactor, &m_cholmod);
-      cholmod_finish(&m_cholmod);
+        M_cholmod_free_factor(&m_cholmodFactor, &m_cholmod);
+      M_cholmod_finish(&m_cholmod);
     }
     
     inline Index cols() const { return m_cholmodFactor->n; }
@@ -279,11 +279,11 @@
     {
       if(m_cholmodFactor)
       {
-        cholmod_free_factor(&m_cholmodFactor, &m_cholmod);
+        M_cholmod_free_factor(&m_cholmodFactor, &m_cholmod);
         m_cholmodFactor = 0;
       }
       cholmod_sparse A = viewAsCholmod(matrix.template selfadjointView<UpLo>());
-      m_cholmodFactor = cholmod_analyze(&A, &m_cholmod);
+      m_cholmodFactor = M_cholmod_analyze(&A, &m_cholmod);
       
       this->m_isInitialized = true;
       this->m_info = Success;
@@ -301,7 +301,7 @@
     {
       eigen_assert(m_analysisIsOk && "You must first call analyzePattern()");
       cholmod_sparse A = viewAsCholmod(matrix.template selfadjointView<UpLo>());
-      cholmod_factorize(&A, m_cholmodFactor, &m_cholmod);
+      M_cholmod_factorize(&A, m_cholmodFactor, &m_cholmod);
       
       this->m_info = Success;
       m_factorizationIsOk = true;
@@ -322,14 +322,14 @@
 
       // note: cd stands for Cholmod Dense
       cholmod_dense b_cd = viewAsCholmod(b.const_cast_derived());
-      cholmod_dense* x_cd = cholmod_solve(CHOLMOD_A, m_cholmodFactor, &b_cd, &m_cholmod);
+      cholmod_dense* x_cd = M_cholmod_solve(CHOLMOD_A, m_cholmodFactor, &b_cd, &m_cholmod);
       if(!x_cd)
       {
         this->m_info = NumericalIssue;
       }
       // TODO optimize this copy by swapping when possible (be carreful with alignment, etc.)
       dest = Matrix<Scalar,Dest::RowsAtCompileTime,Dest::ColsAtCompileTime>::Map(reinterpret_cast<Scalar*>(x_cd->x),b.rows(),b.cols());
-      cholmod_free_dense(&x_cd, &m_cholmod);
+      M_cholmod_free_dense(&x_cd, &m_cholmod);
     }
     
     /** \internal */
@@ -342,14 +342,14 @@
 
       // note: cs stands for Cholmod Sparse
       cholmod_sparse b_cs = viewAsCholmod(b);
-      cholmod_sparse* x_cs = cholmod_spsolve(CHOLMOD_A, m_cholmodFactor, &b_cs, &m_cholmod);
+      cholmod_sparse* x_cs = M_cholmod_spsolve(CHOLMOD_A, m_cholmodFactor, &b_cs, &m_cholmod);
       if(!x_cs)
       {
         this->m_info = NumericalIssue;
       }
       // TODO optimize this copy by swapping when possible (be carreful with alignment, etc.)
       dest = viewAsEigen<DestScalar,DestOptions,DestIndex>(*x_cs);
-      cholmod_free_sparse(&x_cs, &m_cholmod);
+      M_cholmod_free_sparse(&x_cs, &m_cholmod);
     }
     #endif // EIGEN_PARSED_BY_DOXYGEN
     

Modified: pkg/RcppEigen/inst/skeleton/Makevars
===================================================================
--- pkg/RcppEigen/inst/skeleton/Makevars	2011-07-28 15:58:15 UTC (rev 3155)
+++ pkg/RcppEigen/inst/skeleton/Makevars	2011-07-28 16:01:20 UTC (rev 3156)
@@ -1,3 +1,3 @@
 ## Use the R_HOME indirection to support installations of multiple R version
-PKG_LIBS = `$(R_HOME)/bin/Rscript -e "Rcpp:::LdFlags()"`
-
+PKG_LIBS = `$(R_HOME)/bin/Rscript -e "Rcpp:::LdFlags()"` 
+PKG_CPPFLAGS = `$(R_HOME)/bin/Rscript -e "RcppEigen:::RcppEigenCxxFlags()"`

Modified: pkg/RcppEigen/inst/unitTests/runit.sparse.R
===================================================================
--- pkg/RcppEigen/inst/unitTests/runit.sparse.R	2011-07-28 15:58:15 UTC (rev 3155)
+++ pkg/RcppEigen/inst/unitTests/runit.sparse.R	2011-07-28 16:01:20 UTC (rev 3156)
@@ -41,3 +41,42 @@
     colnames(rr) <- NULL
     checkEquals( res, rr, msg = "Sparsematrix wrap")
 }
+
+test.solveCholmod.R <- function() {
+    suppressMessages(require("Matrix", character.only=TRUE))
+    data("KNex", package = "Matrix")
+
+    fx <- cxxfunction( signature(input_ = "list"), '
+    using Eigen::VectorXd;
+    using Eigen::MatrixXd;
+    using Eigen::Lower;
+    using Eigen::Map;
+    using Eigen::MappedSparseMatrix;
+    using Eigen::SparseMatrix;
+    using Eigen::CholmodDecomposition;
+    using Eigen::CholmodAuto;
+    using Eigen::Success;
+
+    List input(input_) ;
+    const SparseMatrix<double>      m1 = input[0];
+    const Map<VectorXd>             v1 = input[1];
+    SparseMatrix<double>            m2(m1.cols(), m1.cols());
+
+    m2.selfadjointView<Lower>().rankUpdate(m1.adjoint());
+
+    CholmodDecomposition<SparseMatrix<double> > ff;
+    ff.setMode(CholmodAuto);
+    ff.compute(m2);
+    VectorXd                       res = ff.solve(m1.adjoint() * v1);
+    
+    return List::create(_["res"]   = res,
+                        _["rows"]  = ff.rows(),
+                        _["cols"]  = ff.cols());
+',
+                      plugin = "RcppEigen",
+                      inc='#include "RcppEigenStubs.cpp"')
+    rr <- fx(KNex)
+    checkEquals(rr[[1]], as.vector(solve(crossprod(KNex[[1]]),
+                                         crossprod(KNex[[1]], KNex[[2]]))),
+                                       "Cholmod solution")
+}



More information about the Rcpp-commits mailing list