[Rcpp-commits] r3142 - in pkg/RcppEigen/inst: include unitTests
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Sat Jul 16 00:53:38 CEST 2011
Author: dmbates
Date: 2011-07-16 00:53:37 +0200 (Sat, 16 Jul 2011)
New Revision: 3142
Modified:
pkg/RcppEigen/inst/include/RcppEigenForward.h
pkg/RcppEigen/inst/include/RcppEigenWrap.h
pkg/RcppEigen/inst/unitTests/runit.RcppEigen.R
Log:
added as methods for Eigen::SparseMatrix and Eigen::MappedSparseMatrix classes and tests of same
Modified: pkg/RcppEigen/inst/include/RcppEigenForward.h
===================================================================
--- pkg/RcppEigen/inst/include/RcppEigenForward.h 2011-07-13 16:46:39 UTC (rev 3141)
+++ pkg/RcppEigen/inst/include/RcppEigenForward.h 2011-07-15 22:53:37 UTC (rev 3142)
@@ -33,19 +33,6 @@
namespace Rcpp {
/* support for wrap */
- // [romain] : don't need all this anymore
- // template<typename Derived> SEXP wrap(const Eigen::EigenBase<Derived>&);
- // template<typename T> SEXP wrap(const Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic>&);
- // template<typename T> SEXP wrap(const Eigen::Matrix<T, Eigen::Dynamic, 1>&);
- // template<typename T> SEXP wrap(const Eigen::Matrix<T, 1, Eigen::Dynamic>&);
- // template<typename T> SEXP wrap(const Eigen::Array<T, Eigen::Dynamic, Eigen::Dynamic>&);
- // template<typename T> SEXP wrap(const Eigen::Array<T, Eigen::Dynamic, 1>&);
- // template<typename T> SEXP wrap(const Eigen::Map<Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic> >&);
- // template<typename T> SEXP wrap(const Eigen::Map<Eigen::Matrix<T, Eigen::Dynamic, 1> >&);
- // template<typename T> SEXP wrap(const Eigen::Map<Eigen::Matrix<T, 1, Eigen::Dynamic> >&);
- // template<typename T> SEXP wrap(const Eigen::Map<Eigen::Array<T, Eigen::Dynamic, Eigen::Dynamic> >&);
- // template<typename T> SEXP wrap(const Eigen::Map<Eigen::Array<T, Eigen::Dynamic, 1> >&);
- // template<typename T> SEXP wrap(const Eigen::SparseMatrix<T>&);
template<typename T> SEXP wrap(const Eigen::Map<Eigen::SparseMatrix<T> >&);
namespace traits {
@@ -56,6 +43,8 @@
template<typename T> class Exporter< Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic> >;
template<typename T> class Exporter< Eigen::Matrix<T, Eigen::Dynamic, 1> >;
template<typename T> class Exporter< Eigen::Matrix<T, 1, Eigen::Dynamic> >;
+ template<typename T> class Exporter< Eigen::MappedSparseMatrix<T> >;
+ template<typename T> class Exporter< Eigen::SparseMatrix<T> >;
} // namespace traits
Modified: pkg/RcppEigen/inst/include/RcppEigenWrap.h
===================================================================
--- pkg/RcppEigen/inst/include/RcppEigenWrap.h 2011-07-13 16:46:39 UTC (rev 3141)
+++ pkg/RcppEigen/inst/include/RcppEigenWrap.h 2011-07-15 22:53:37 UTC (rev 3142)
@@ -100,79 +100,6 @@
} /* namespace RcppEigen */
- // /* wrap */
- // [romain] : no longer necessary
- // template <typename Derived>
- // SEXP wrap(const Eigen::EigenBase<Derived>& object) {
- // //FIXME: Check IsRowMajor and transpose if needed
- // ::Rcpp::RObject x = ::Rcpp::wrap(object.data(), object.data() + object.size());
- // if (object.ColsAtCompileTime == 1) return x; // represented as a vector
- // x.attr("dim") = ::Rcpp::Dimension(object.rows(), object.cols());
- // }
- //
- // template <typename T>
- // SEXP wrap(const Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic>& data) {
- // return RcppEigen::Eigen_wrap(data, Dimension(data.rows(), data.cols()));
- // }
- //
- // template <typename T>
- // SEXP wrap(const Eigen::Matrix<T, Eigen::Dynamic, 1>& object ){
- // return ::Rcpp::wrap(object.data(), object.data() + object.size());
- // }
- //
- // template <typename T>
- // SEXP wrap( const Eigen::Matrix<T, 1, Eigen::Dynamic>& data ){
- // return RcppEigen::Eigen_wrap(data, Dimension(1, data.size()));
- // }
- //
- // template <typename T>
- // SEXP wrap(const Eigen::Array<T, Eigen::Dynamic, Eigen::Dynamic>& data) {
- // return RcppEigen::Eigen_wrap(data, Dimension(data.rows(), data.cols()));
- // }
- //
- // template <typename T>
- // SEXP wrap(const Eigen::Array<T, Eigen::Dynamic, 1>& object ){
- // return ::Rcpp::wrap(object.data(), object.data() + object.size());
- // }
- //
- // template <typename T>
- // SEXP wrap(const Eigen::Map<Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic> >& data) {
- // return RcppEigen::Eigen_wrap(data, Dimension(data.rows(), data.cols()));
- // }
- //
- // template <typename T>
- // SEXP wrap(const Eigen::Map<Eigen::Matrix<T, Eigen::Dynamic, 1> >& object ){
- // return ::Rcpp::wrap(object.data(), object.data() + object.size());
- // }
- //
- // template <typename T>
- // SEXP wrap(const Eigen::Map<Eigen::Matrix<T, 1, Eigen::Dynamic> >& data ){
- // return RcppEigen::Eigen_wrap(data, Dimension(1, data.size()));
- // }
- //
- // template <typename T>
- // SEXP wrap(const Eigen::Map<Eigen::Array<T, Eigen::Dynamic, Eigen::Dynamic> >& data) {
- // return RcppEigen::Eigen_wrap(data, Dimension(data.rows(), data.cols()));
- // }
- //
- // template <typename T>
- // SEXP wrap(const Eigen::Map<Eigen::Array<T, Eigen::Dynamic, 1> >& object ){
- // return ::Rcpp::wrap(object.data(), object.data() + object.size());
- // }
- // template <typename T>
- // SEXP wrap(const Eigen::Map<Eigen::SparseMatrix<T> >& object ) {
- // int nnz = object.nonZeros(), p = object.outerSize();
- // Dimension dim(object.innerSize(), p);
- // const int *ip = object._innerIndexPtr(), *pp = object._outerIndexPtr();
- // const T *xp = object._valuePtr();
- // IntegerVector iv(ip, ip + nnz), pv(pp, pp + p + 1);
- // NumericVector xv(xp, xp + nnz);
- //
- // return ::Rcpp::wrap(List::create(_["Dim"] = dim,
- // _["i"] = iv,
- // _["p"] = pv,
- // _["x"] = xv));
- // }
/* support for Rcpp::as */
@@ -238,7 +165,61 @@
Exporter(SEXP x) :
MatrixExporter< Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic>, T >(x){}
};
-
+
+ template<typename T>
+ class Exporter<Eigen::MappedSparseMatrix<T> > {
+ public:
+ Exporter(SEXP x)
+ : d_x(x), d_dims(d_x.slot("Dim")), d_i(d_x.slot("i")), d_p(d_x.slot("p")) {
+ if (!d_x.is("CsparseMatrix"))
+ throw std::invalid_argument("Need S4 class CsparseMatrix for an mapped sparse matrix");
+ const int RTYPE = ::Rcpp::traits::r_sexptype_traits<T>::rtype ;
+ SEXP xx = d_x.slot("x");
+ if (TYPEOF(xx) != RTYPE)
+ throw std::invalid_argument("Wrong R type for mapped sparse matrix");
+ typedef typename ::Rcpp::traits::storage_type<RTYPE>::type STORAGE;
+ d_start = ::Rcpp::internal::r_vector_start<RTYPE,STORAGE>(xx);
+ }
+ Eigen::MappedSparseMatrix<T> get() {
+ return Eigen::MappedSparseMatrix<T>(d_dims[0], d_dims[1], d_p[d_dims[1]],
+ d_p.begin(), d_i.begin(), d_start);
+ }
+ protected:
+ S4 d_x;
+ T* d_start;
+ IntegerVector d_dims, d_i, d_p;
+ };
+
+ template<typename T>
+ class Exporter<Eigen::SparseMatrix<T> > {
+ public:
+ Exporter(SEXP x)
+ : d_x(x), d_dims(d_x.slot("Dim")), d_i(d_x.slot("i")), d_p(d_x.slot("p")) {
+ if (!d_x.is("CsparseMatrix"))
+ throw std::invalid_argument("Need S4 class CsparseMatrix for an mapped sparse matrix");
+ const int RTYPE = ::Rcpp::traits::r_sexptype_traits<T>::rtype ;
+ SEXP xx = d_x.slot("x");
+ if (TYPEOF(xx) != RTYPE) // should coerce instead - see Rcpp/inst/include/Rcpp/internal/export.h
+ throw std::invalid_argument("Wrong R type for sparse matrix");
+ typedef typename ::Rcpp::traits::storage_type<RTYPE>::type STORAGE;
+ d_start = ::Rcpp::internal::r_vector_start<RTYPE,STORAGE>(xx);
+ }
+ Eigen::SparseMatrix<T> get() {
+ Eigen::SparseMatrix<T> ans(d_dims[0], d_dims[1]);
+ ans.reserve(d_p[d_dims[1]]);
+ for(int j = 0; j < d_dims[1]; ++j) {
+ ans.startVec(j);
+ for (int k = d_p[j]; k < d_p[j + 1]; ++k) ans.insertBack(d_i[k], j) = d_start[k];
+ }
+ ans.finalize();
+ return ans;
+ }
+ protected:
+ S4 d_x;
+ T* d_start;
+ IntegerVector d_dims, d_i, d_p;
+ };
+
} // namespace traits
}
Modified: pkg/RcppEigen/inst/unitTests/runit.RcppEigen.R
===================================================================
--- pkg/RcppEigen/inst/unitTests/runit.RcppEigen.R 2011-07-13 16:46:39 UTC (rev 3141)
+++ pkg/RcppEigen/inst/unitTests/runit.RcppEigen.R 2011-07-15 22:53:37 UTC (rev 3142)
@@ -200,253 +200,51 @@
), msg = "as<MMat>" )
}
-if (FALSE) {
-
-test.as.Mat <- function(){
-
- fx <- cxxfunction( signature(input_ = "list" ) , '
+test.as.MSpMat <- function() {
+ suppressMessages(require("Matrix"))
+ data("KNex", package = "Matrix")
+ fx <- cxxfunction( signature(input_ = "list"), '
List input(input_) ;
- Eigen::MatrixXi m1 = input[0] ; /* implicit as */
- Eigen::MatrixXd m2 = input[1] ; /* implicit as */
- Eigen::Matrix<unsigned int, Eigen::Dynamic, Eigen::Dynamic> m3 = input[0] ; /* implicit as */
- Eigen::MatrixXf m4 = input[1] ; /* implicit as */
+ const Eigen::MappedSparseMatrix<double> m1 = input[0]; // maps share storage and do not allow conversion
- List res = List::create(m1.sum(), m2.sum(), m3.sum(), m4.sum());
+ List res = List::create(_["nnz"] = m1.nonZeros(),
+ _["nr"] = m1.rows(),
+ _["nc"] = m1.cols(),
+ _["inSz"] = m1.innerSize(),
+ _["outSz"] = m1.outerSize(),
+ _["sum"] = m1.sum());
return res ;
- ', plugin = "RcppEigen" )
- integer_mat <- matrix( as.integer(diag(4)), ncol = 4, nrow = 4 )
- numeric_mat <- diag(5)
- res <- fx( list( integer_mat, numeric_mat ) )
- checkEquals( unlist( res), c(4L, 5L, 4L, 5L ), msg = "as<Mat>" )
-}
-
-test.wrap.Glue <- function(){
-
- fx <- cxxfunction( , '
-
- arma::mat m1 = arma::eye<arma::mat>( 3, 3 ) ;
- arma::mat m2 = arma::eye<arma::mat>( 3, 3 ) ;
-
- List res ;
- res["mat+mat"] = m1 + m2 ;
- return res ;
-
', plugin = "RcppEigen" )
- res <- fx()
- checkEquals( res[[1]], 2*diag(3), msg = "wrap(Glue)" )
+ KNX <- KNex[[1]]
+ res <- fx(KNex)
+ checkEquals(unname(unlist(res)),
+ c(nnzero(KNX), nrow(KNX), ncol(KNX), nrow(KNX), ncol(KNX), sum(KNX at x)),
+ msg = "as<MSPMatrix>")
}
-test.wrap.Op <- function(){
-
- fx <- cxxfunction( , '
-
- arma::mat m1 = arma::eye<arma::mat>( 3, 3 ) ;
-
- List res ;
- res["- mat"] = - m1 ;
- return res ;
-
- ', plugin = "RcppEigen" )
- res <- fx()
- checkEquals( res[[1]], -1*diag(3), msg = "wrap(Op)" )
-}
-
-test.as.Col <- function(){
- fx <- cxxfunction( signature(input_ = "list" ) , '
-
+test.as.SpMat <- function() {
+ suppressMessages(require("Matrix"))
+ data("KNex", package = "Matrix")
+ fx <- cxxfunction( signature(input_ = "list"), '
List input(input_) ;
- arma::icolvec m1 = input[0] ; /* implicit as */
- arma::colvec m2 = input[1] ; /* implicit as */
- arma::ucolvec m3 = input[0] ; /* implicit as */
- arma::fcolvec m4 = input[1] ; /* implicit as */
+ const Eigen::SparseMatrix<double> m1 = input[0];
- List res = List::create(
- arma::accu( m1 ),
- arma::accu( m2 ),
- arma::accu( m3 ),
- arma::accu( m4 ) ) ;
+ List res = List::create(_["nnz"] = m1.nonZeros(),
+ _["nr"] = m1.rows(),
+ _["nc"] = m1.cols(),
+ _["inSz"] = m1.innerSize(),
+ _["outSz"] = m1.outerSize(),
+ _["sum"] = m1.sum());
return res ;
-
', plugin = "RcppEigen" )
- res <- fx( list( 1:10, as.numeric(1:10) ) )
- checkEquals( unlist( res ), rep(55.0, 4 ), msg = "as<Col>" )
+ KNX <- KNex[[1]]
+ res <- fx(KNex)
+ checkEquals(unname(unlist(res)),
+ c(nnzero(KNX), nrow(KNX), ncol(KNX), nrow(KNX), ncol(KNX), sum(KNX at x)),
+ msg = "as<MSPMatrix>")
}
-
-test.as.Row <- function(){
- fx <- cxxfunction( signature(input_ = "list" ) , '
-
- List input(input_) ;
- arma::irowvec m1 = input[0] ; /* implicit as */
- arma::rowvec m2 = input[1] ; /* implicit as */
- arma::urowvec m3 = input[0] ; /* implicit as */
- arma::frowvec m4 = input[1] ; /* implicit as */
-
- List res = List::create(
- arma::accu( m1 ),
- arma::accu( m2 ),
- arma::accu( m3 ),
- arma::accu( m4 ) ) ;
- return res ;
-
- ', plugin = "RcppEigen" )
-
- res <- fx( list( 1:10, as.numeric(1:10) ) )
- checkEquals( unlist( res ), rep(55.0, 4 ), msg = "as<Row>" )
-}
-
-test.cxmat <- function(){
-
- fx <- cxxfunction( signature() , '
-
- arma::cx_mat m1 = arma::eye<arma::cx_mat> ( 3, 3 ) ;
- arma::cx_fmat m2 = arma::eye<arma::cx_fmat>( 3, 3 ) ;
- return List::create( _["double"] = m1, _["float"] = m2 ) ;
-
- ', plugin = "RcppEigen" )
- checkEquals( fx(),
- list( double = (1+0i)*diag(3), float = (1+0i)*diag(3) ),
- msg = "support for complex matrices" )
-
-}
-
-test.mtOp <- function(){
-
- fx <- cxxfunction( signature() , '
-
- std::complex<double> x( 1.0, 2.0 ) ;
- arma::mat m1 = arma::eye<arma::mat> ( 3, 3 ) ;
-
- return wrap( x * m1 ) ;
-
- ', plugin = "RcppEigen" )
- checkEquals( fx(),
- (1+2i)*diag(3),
- msg = "support for mtOp" )
-
-}
-
-test.mtGlue <- function(){
-
- fx <- cxxfunction( signature() , '
-
- arma::imat m2 = arma::eye<arma::imat> ( 3, 3 ) ;
- arma::mat m1 = arma::eye<arma::mat> ( 3, 3 ) ;
-
- return wrap( m1 + m2 ) ;
-
- ', plugin = "RcppEigen" )
- checkEquals( fx(),
- 2.0 * diag(3) ,
- msg = "support for mtOp" )
-
-}
-
-
-test.sugar <- function(){
-
- fx <- cxxfunction( signature(x= "numeric") , '
- NumericVector xx(x) ;
- arma::mat m = forward( xx + xx ) ;
- return wrap( m ) ;
-
- ', plugin = "RcppEigen" )
- checkEquals( fx(1:10),
- matrix( 2*(1:10), nrow = 10 ) ,
- msg = "RcppEigen and sugar" )
-
-}
-
-test.sugar.cplx <- function(){
-
- fx <- cxxfunction( signature(x= "complex") , '
- ComplexVector xx(x) ;
- arma::cx_mat m = forward( exp( xx ) ) ;
-
- return wrap( m ) ;
-
- ', plugin = "RcppEigen" )
- x <- 1:10*(1+1i)
- checkEquals( fx(x),
- matrix( exp(x), nrow = 10 ) ,
- msg = "RcppEigen and sugar (complex)" )
-
-}
-
-test.armadillo.sugar.ctor <- function(){
-
- fx <- cxxfunction( signature(x= "numeric") , '
- NumericVector xx(x) ;
- arma::mat m = xx + xx ;
- arma::colvec co = xx ;
- arma::rowvec ro = xx ;
- return List::create(
- _["mat"] = m + m,
- _["rowvec"] = ro,
- _["colvec"] = co
- );
- ', plugin = "RcppEigen" )
- checkEquals( fx(1:10),
- list(
- mat = matrix( 4*(1:10), nrow = 10 ),
- rowvec = matrix( 1:10, nrow = 1 ),
- colvec = matrix( 1:10, ncol = 1 )
- )
- ,
- msg = "Mat( sugar expression )" )
-
-}
-
-
-test.armadillo.sugar.matrix.ctor <- function(){
-
- inc <- '
- double norm( double x, double y){
- return ::sqrt( x*x + y*y );
- }
- '
- fx <- cxxfunction( signature(x= "numeric") , '
- NumericVector xx(x) ;
- NumericVector yy = NumericVector::create( 1 ) ;
- arma::mat m = diag( xx ) ;
- arma::colvec co = outer( xx, yy, ::norm ) ;
- arma::rowvec ro = outer( yy, xx, ::norm ) ;
- return List::create(
- _["mat"] = m + m ,
- _["rowvec"] = ro,
- _["colvec"] = co
- );
- ', plugin = "RcppEigen", includes = inc )
- res <- fx(1:10)
- norm <- function(x, y) sqrt( x*x + y*y )
- checkEquals( res,
- list(
- mat = diag(2*(1:10)),
- rowvec = outer( 1, 1:10, norm ),
- colvec = outer( 1:10, 1, norm )
- ),
- msg = "Mat( sugar expression )" )
-
-}
-
-test.armadillo.rtti.check <- function() {
-
- inc <- '
- void blah(arma::mat& X) {
- X.set_size(5,5);
- }
- '
- src <- '
- arma::vec V;
- blah(V); // if blah() worked, we have a problem
- '
- fun <- cxxfunction(signature(), body=src, inc=inc, plugin = "RcppEigen")
-
- checkException(fun(), msg="RTTI check on matrix constructor exception")
-
-}
-}
More information about the Rcpp-commits
mailing list