将 Rcpp 数字向量转换为 boost:ublas:vector

ste*_*980 2 c++ boost r rcpp

我正在尝试从rtype对象转换为ublasfrom boost。

使用我从Rcpp 开发列表中找到的一些关于 ublas 的代码,我能够返回ublas包装为rtype.

例如

// Converts from ublas to rtype

template <typename T>   
Rcpp::Vector< Rcpp::traits::r_sexptype_traits<T>::rtype >
ublas2rcpp( const boost::numeric::ublas::vector<T>& x ){
  return Rcpp::Vector< Rcpp::traits::r_sexptype_traits<T>::rtype >(
      x.begin(), x.end()
  ) ;
}
Run Code Online (Sandbox Code Playgroud)

为了模仿行为导入行为,我目前ublas 在相应长度的向量上使用循环,并将所有内容分配rtype给它。

从循环切换会提高性能吗?

// My attempt to convert from rtype to ublas

// so R can find the libraries
//[[Rcpp::depends(BH)]]
//[[Rcpp::plugins("cpp11")]]

#include <Rcpp.h>


#include <cstdlib>
#include <iostream>
#include <fstream>


#include <boost/numeric/odeint.hpp>
#include <boost/numeric/ublas/matrix.hpp>
#include <boost/numeric/ublas/io.hpp>

using namespace std;
using namespace Rcpp;
using namespace boost::numeric::ublas;

typedef boost::numeric::ublas::vector< double > vector_type;
typedef boost::numeric::ublas::matrix< double > matrix_type;

template <typename T>    /// this need to be fixed up, hopefully working for now though
Rcpp::Vector< Rcpp::traits::r_sexptype_traits<T>::rtype >
ublas2rcpp( const boost::numeric::ublas::vector<T>& x ){
  return Rcpp::Vector< Rcpp::traits::r_sexptype_traits<T>::rtype >(
      x.begin(), x.end()
  ) ;
}

// [[Rcpp::export]]
NumericVector main(NumericVector x1)
{


  int L =x1.length();


  vector_type x(L , 0 ); // initialize the vector to all zero

  for(int i=0;i<L;i++)
  {

    x(i) =  x1(i);

  }


return(ublas2rcpp(x));   

}
Run Code Online (Sandbox Code Playgroud)

coa*_*ess 5

首先,抱歉耽搁了。希望我所做的可以弥补它。

话虽如此,让我们这样做!

介绍

这个问题是正在尝试的两个方面:

  1. 从 R 转换为 C++ ( Rcpp::as<T>(obj))
  2. 从 C++ 转换为 R ( Rcpp::wrap(obj))

幸运的是,有一个名为Extending Rcpp的精彩 Rcpp 小插图可以解决自定义对象。遗憾的是,与其他文档相比,小插图的清晰度还有很多不足之处。(我可能会尝试做一个 PR 来改进它。)

因此,我将尝试通过一些评论来引导您完成这些步骤。请注意,所使用的方法是通过模板和部分专业化,最终会产生一些很好的自动魔术。

的解释

阶段 1 - 前向声明

在第一阶段,我们必须在参与之前声明我们希望使用的功能的意图Rcpp.h。为此,我们将加载不同的头文件并向Rcpp::traits命名空间添加一些定义。

原则上,当我们开始写入文件时,我们必须加载的第一个标头是RcppCommon.h而不是通常的Rcpp.h!! 如果我们没有在Rcpp.h通话前进行前向声明,我们将无法正确注册我们的分机。

然后,我们必须添加不同的插件标记,sourceCpp()以便在代码编译期间设置适当的标志。在插件之后,我们将包含我们想要使用的实际 boost 头文件。最后,我们必须在命名空间中添加两个特殊的 Rcpp 函数声明Rcpp::as<T>(obj)和。要启用多种类型,我们必须创建一个类,而不是更直接地调用.Rcpp::wrap(obj)Rcpp::traitsExportertemplate <> ClassName as( SEXP )

#include <RcppCommon.h>

// Flags for C++ compiler

// [[Rcpp::depends(BH)]]
// [[Rcpp::plugins("cpp11")]]

// Third party library includes that provide the template class of ublas
#include <boost/numeric/ublas/matrix_sparse.hpp>
#include <boost/numeric/ublas/matrix.hpp>

// Provide Forward Declarations
namespace Rcpp {

   namespace traits{

   // Setup non-intrusive extension via template specialization for
   // 'ublas' class boost::numeric::ublas

   // Support for wrap
   template <typename T> SEXP wrap(const boost::numeric::ublas::vector<T> & obj);

   // Support for as<T>
   template <typename T> class Exporter< boost::numeric::ublas::vector<T> >;

  }
}
Run Code Online (Sandbox Code Playgroud)

第 2 阶段 - 包括 Rcpp.h

仅仅为了申报进口订单而有一个阶段可能看起来很琐碎,但如果Rcpp.h在前向申报之前包含在内,Rcpp::traits则不会更新,我们进入深渊。

因此:

// >> Place <Rcpp.h> AFTER the forward declaration!!!! <<

#include <Rcpp.h>


// >> Place Definitions of Forward Declarations AFTER <Rcpp.h>!!!! <<
Run Code Online (Sandbox Code Playgroud)

第 3 阶段 - 实施扩展

现在,我们必须实际实现前向声明。特别是,唯一会有点问题的实现是as<>因为wrap()是直截了当的。

wrap()

为了实现,wrap()我们必须求助于 Rcpp 中的内置类型转换索引,称为Rcpp::traits::r_sexptype_traits<T>::rtype. 由此,我们能够获得一个int包含 的RTYPE,然后构造一个Rcpp::Vector。对于矩阵的构造,同样的想法也适用。

as()

对于as<>(),我们需要考虑将传入的模板。此外,我们typedef在Exporter类定义的正下方设置了 ,以便轻松定义OUT要在get()方法中使用的对象。除此之外,我们使用相同的技巧在 C++T类型和R类型之间来回移动。

为了完成as<>从 R 到 C++ 的直接移植,我不得不做一些肮脏的事情:我复制了向量内容。管理此输出的代码get()在Exporter类中给出。您可能希望花一些时间研究使用指针更改分配。我不太熟悉,ublas所以我没有看到解决指针传递的简单方法。

// Define template specializations for as<> and wrap

namespace Rcpp {

namespace traits{

// Defined wrap case
template <typename T> SEXP wrap(const boost::numeric::ublas::vector<T> & obj){
  const int RTYPE = Rcpp::traits::r_sexptype_traits<T>::rtype ;

  return Rcpp::Vector< RTYPE >(obj.begin(), obj.end());
};


// Defined as< > case
template<typename T>
class Exporter< boost::numeric::ublas::vector<T> > {
  typedef typename boost::numeric::ublas::vector<T> OUT ;

  // Convert the type to a valid rtype. 
  const static int RTYPE = Rcpp::traits::r_sexptype_traits< T >::rtype ;
  Rcpp::Vector<RTYPE> vec;

public:
  Exporter(SEXP x) : vec(x) {
    if (TYPEOF(x) != RTYPE)
      throw std::invalid_argument("Wrong R type for mapped 1D array");
  }
  OUT get() {

    // Need to figure out a way to perhaps do a pointer pass?
    OUT x(vec.size());

    std::copy(vec.begin(), vec.end(), x.begin()); // have to copy data

    return x;
  }
} ;


}
}
Run Code Online (Sandbox Code Playgroud)

第 4 阶段 - 测试

好的,让我们看看我们所做的工作是否得到了回报(剧透确实如此!剧透)。要进行检查,我们应该查看两个不同的区域:

  1. 功能内的跟踪诊断和;
  2. 自动测试。

两者都在下面给出。请注意,我选择将ublas设置缩短为:

// Here we define a shortcut to the boost ublas class to enable multiple ublas types via a template.
// ublas::vector<T> => ublas::vector<double>, ... , ublas::vector<int>
namespace ublas = ::boost::numeric::ublas;
Run Code Online (Sandbox Code Playgroud)

跟踪诊断

// [[Rcpp::export]]
void containment_test(Rcpp::NumericVector x1) {

  Rcpp::Rcout << "Converting from Rcpp::NumericVector to ublas::vector<double>" << std::endl;

  ublas::vector<double> x = Rcpp::as< ublas::vector<double> >(x1); // initialize the vector to all zero

  Rcpp::Rcout << "Running output test with ublas::vector<double>" << std::endl;

  for (unsigned i = 0; i < x.size (); ++ i)
    Rcpp::Rcout  << x(i) << std::endl;

  Rcpp::Rcout << "Converting from ublas::vector<double> to Rcpp::NumericVector" << std::endl;

  Rcpp::NumericVector test = Rcpp::wrap(x);

  Rcpp::Rcout << "Running output test with Rcpp::NumericVector" << std::endl;

  for (unsigned i = 0; i < test.size (); ++ i)
    Rcpp::Rcout  << test(i) << std::endl;

}
Run Code Online (Sandbox Code Playgroud)

测试电话:

containment_test(c(1,2,3,4))
Run Code Online (Sandbox Code Playgroud)

结果:

Converting from Rcpp::NumericVector to ublas::vector<double>
Running output test with ublas::vector<double>
1
2
3
4
Converting from ublas::vector<double> to Rcpp::NumericVector
Running output test with Rcpp::NumericVector
1
2
3
4
Run Code Online (Sandbox Code Playgroud)

该测试按预期执行。进入下一个测试!

自动测试

// [[Rcpp::export]]
ublas::vector<double> automagic_ublas_rcpp(ublas::vector<double> x1) {
  return x1;
}
Run Code Online (Sandbox Code Playgroud)

测试电话:

automagic_ublas_rcpp(c(1,2,3.2,1.2))
Run Code Online (Sandbox Code Playgroud)

结果:

[1] 1.0 2.0 3.2 1.2
Run Code Online (Sandbox Code Playgroud)

成功!

第 5 阶段 - 一起为Cntrl+C和Cntrl+V

下面是上面stage给出的代码块的组合。如果您复制并粘贴到该.cpp文件中,那么一切都应该工作。如果没有,请告诉我。

// -------------- Stage 1: Forward Declarations with `RcppCommon.h`

#include <RcppCommon.h>

// Flags for C++ compiler
// [[Rcpp::depends(BH)]]
// [[Rcpp::plugins("cpp11")]]

// Third party library includes that provide the template class of ublas
#include <boost/numeric/ublas/matrix_sparse.hpp>
#include <boost/numeric/ublas/matrix.hpp>


// Here we use ublas_vec to enable multiple ublas types via a template.
// ublas::vector<T> => ublas::vector<double>, ... , ublas::vector<int>
namespace ublas = ::boost::numeric::ublas;


// Provide Forward Declarations
namespace Rcpp {

  namespace traits{

    // Setup non-intrusive extension via template specialization for
    // 'ublas' class boost::numeric::ublas

    // Support for wrap
    template <typename T> SEXP wrap(const boost::numeric::ublas::vector<T> & obj);

    // Support for as<T>
    template <typename T> class Exporter< boost::numeric::ublas::vector<T> >;

  }
}

// -------------- Stage 2: Including Rcpp.h

// >> Place <Rcpp.h> AFTER the forward declaration!!!! <<

#include <Rcpp.h>


// >> Place Definitions of Forward Declarations AFTER <Rcpp.h>!!!! <<


// -------------- Stage 3: Implementation of Declarations

// Define template specializations for as<> and wrap

namespace Rcpp {

  namespace traits{

    // Defined wrap case
    template <typename T> SEXP wrap(const boost::numeric::ublas::vector<T> & obj){
      const int RTYPE = Rcpp::traits::r_sexptype_traits<T>::rtype ;

      return Rcpp::Vector< RTYPE >(obj.begin(), obj.end());
    };


    // Defined as< > case
    template<typename T>
    class Exporter< boost::numeric::ublas::vector<T> > {
      typedef typename boost::numeric::ublas::vector<T> OUT ;

      // Convert the type to a valid rtype. 
      const static int RTYPE = ::Rcpp::traits::r_sexptype_traits< T >::rtype ;
      Rcpp::Vector<RTYPE> vec;

    public:
      Exporter(SEXP x) : vec(x) {
        if (TYPEOF(x) != RTYPE)
          throw std::invalid_argument("Wrong R type for mapped 1D array");
      }
      OUT get() {

        // Need to figure out a way to perhaps do a pointer pass?
        OUT x(vec.size());

        std::copy(vec.begin(), vec.end(), x.begin()); // have to copy data

        return x;
      }
    } ;


  }
}

// -------------- Stage 4: Tests

// [[Rcpp::export]]
ublas::vector<double> automagic_ublas_rcpp(ublas::vector<double> x1) {
  return x1;
}


// [[Rcpp::export]]
void containment_test(Rcpp::NumericVector x1) {

  Rcpp::Rcout << "Converting from Rcpp::NumericVector to ublas::vector<double>" << std::endl;

  ublas::vector<double> x = Rcpp::as< ublas::vector<double> >(x1); // initialize the vector to all zero

  Rcpp::Rcout << "Running output test with ublas::vector<double>" << std::endl;

  for (unsigned i = 0; i < x.size (); ++ i)
    Rcpp::Rcout  << x(i) << std::endl;

  Rcpp::Rcout << "Converting from ublas::vector<double> to Rcpp::NumericVector" << std::endl;

  Rcpp::NumericVector test = Rcpp::wrap(x);

  Rcpp::Rcout << "Running output test with Rcpp::NumericVector" << std::endl;

  for (unsigned i = 0; i < test.size (); ++ i)
    Rcpp::Rcout  << test(i) << std::endl;

}
Run Code Online (Sandbox Code Playgroud)

闭幕致辞

哇...那太多了。希望以上提供了足够的理由,因为我相信您可能希望将1D向量扩展到 aublas::matrix等等。此外,等待应该是值得的,因为您现在拥有自动转换的魔力,Rcpp因此无需ublas2rcpp()在return()语句中调用。实际上,您可以简化将函数的返回类型指定为ublas::vector<double>!

在其他新闻中,我认为可能阻止您创建上面的备用模板函数(例如rcpp2ublas)的原因是无法推断出 C++T类型Rcpp::Vector是什么。就个人而言,在Rcpp GitHub Repo 中挖掘之后,我不确定这样的转换索引是否存在,或者是否应该由于使用而Shield使用它。深入挖掘转化是另一天的冒险。