将动态矢量从C返回到R.

Question

我在C中编写一些代码,从R动态调用.

此代码生成一个随机泊松过程的路径,直到所需的时间T.因此,在每次调用我的C函数时,返回的向量的长度将根据生成的随机数而不同.

我必须创建什么R数据结构？一个LISTSXP？另一个？

我怎样才能创建它,我怎样才能添加它？特别是如何将它还给R？

感谢帮助.

Answer 1

这真的取决于你想要用作临时结构的东西,因为最终你必须为结果分配一个向量.所以无论你将使用什么,都不是你会回来的.有几种可能性:

1. 使用Calloc+ Realloc+ Free允许您根据需要扩展临时内存.获得完整集后,分配结果向量并将其返回.
2. 如果你可以轻易高估大小,你可以在SETLENGTH返回之前过度分配结果向量并使用.但是,这有问题,因为结果将保持过度分配,直到稍后重复.
3. 您可以使用您所暗示的链接列表中的矢量块,即分配和保护一个您可以根据需要将矢量附加到尾部的一个旋转.最后,您可以分配返回向量并复制已分配块的内容.这比上述两者更复杂.

每个都有缺点和好处,所以它真的取决于你的应用程序选择最适合你的.

编辑:添加了使用pairlists方法的示例.我仍然会推荐这种Realloc方法,因为它更容易,但仍然如此:

#define BLOCK_SIZE xxx /* some reasonable size for increments - could be adaptive, too */ 
SEXP block; /* last vector block */
SEXP root = PROTECT(list1(block = allocVector(REALSXP, BLOCK_SIZE)));
SEXP tail = root;
double *values = REAL(block);
int count = 0, total = 0;
do { /* your code to generate values - if you want to add one 
        first try to add it to the existing block, otherwise allocate new one */
    if (count == BLOCK_SIZE) { /* add a new block when needed */
        tail = SETCDR(tail, list1(block = allocVector(REALSXP, BLOCK_SIZE)));
        values = REAL(block);
        total += count;
        count = 0;
    }
    values[count++] = next_value;
} while (...);
total += count;
/* when done, we need to create the result vector */
{
    SEXP res = allocVector(REALSXP, total);
    double *res_values = REAL(res);
    while (root != R_NilValue) {
        int size = (CDR(root) == R_NilValue) ? count : BLOCK_SIZE;
        memcpy(res_values, REAL(CAR(root)), sizeof(double) * size);
        res_values += size;
        root = CDR(root);
    }
    UNPROTECT(1);
    return res;
 }

Answer 2

如果您愿意从 C 切换到 C++，您可以免费获得附加的 Rcpp 层。以下是（仍然相当不完整）RcppExample 包页面中的示例：

#include <RcppClassic.h>
#include <cmath>

RcppExport SEXP newRcppVectorExample(SEXP vector) {
BEGIN_RCPP

    Rcpp::NumericVector orig(vector);                // keep a copy 
    Rcpp::NumericVector vec(orig.size());            // create vector same size

    // we could query size via
    //   int n = vec.size();
    // and loop over the vector, but using the STL is so much nicer
    // so we use a STL transform() algorithm on each element
    std::transform(orig.begin(), orig.end(), vec.begin(), ::sqrt);

    return Rcpp::List::create(Rcpp::Named( "result" ) = vec,
                              Rcpp::Named( "original" ) = orig) ;

END_RCPP
}

正如你所看到的，没有显式的内存分配、释放PROTECT/UNPROTECT等，你就得到了一个第一类 R 列表对象。

还有更多的例子，包含在其他 SO 问题中，例如这个。

编辑： 你并没有真正说出你的路径会做什么，但作为一个简单的说明，这里是使用 Rcpp 添加的 C++ 代码cumsum()，rpois()其行为就像在 R 中一样：

R> library(inline)
R> 
R> fun <- cxxfunction(signature(ns="integer", lambdas="numeric"),
+                    plugin="Rcpp",
+                    body='
+    int n = Rcpp::as<int>(ns);
+    double lambda = Rcpp::as<double>(lambdas);
+ 
+    Rcpp::RNGScope tmp;                  // make sure RNG behaves
+ 
+    Rcpp::NumericVector vec = cumsum( rpois( n, lambda ) );
+ 
+    return vec;
+ ')
R> set.seed(42)
R> fun(3, 0.3)
[1] 1 2 2
R> fun(4, 0.4)
[1] 1 1 1 2

作为证明，回到 R 中，如果我们设置种子，我们可以生成完全相同的数字：

R> set.seed(42)
R> cumsum(rpois(3, 0.3))
[1] 1 2 2
R> cumsum(rpois(4, 0.4))
[1] 1 1 1 2
R>