返回从 C 到 R 的动态向量

Question

我正在用 C 编写一些代码，以便从 R 动态调用。

该代码生成一条随机 泊松 过程的路径，直到所需的时间 T。 因此，每次调用我的 C 函数时，返回向量的长度都会有所不同，具体取决于生成的随机数。

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

如何创建它以及如何附加它？特别是我如何将其返回给 R？

感谢您的帮助。

Answer 1

这实际上取决于您想要使用什么作为临时结构，因为最终您无论如何都必须为结果分配一个向量。因此，无论您将使用什么，都不会返回什么。有几种可能性：

1. 使用Calloc + Realloc + Free，它允许您根据需要扩展临时内存。一旦获得完整的集合，就可以分配结果向量并返回它。
2. 如果您可以轻松高估大小，则可以过度分配结果向量并在返回之前使用 SETLENGTH 。不过，这存在问题，因为结果将保持过度分配，直到稍后重复。
3. 您可以使用您暗示的向量块链表，即分配和保护一个对列表，根据需要将向量附加到尾部。最后，您分配返回向量并复制您分配的块的内容。这比上面两个更复杂。

它们各自都有缺点和优点，因此实际上取决于您的应用程序来选择最适合您的一种。

编辑：添加了使用配对列表方法的示例。我仍然推荐 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 列表对象。

还有更多示例，包括 其他问题，例如就像这个。

编辑：你并没有真正说出你的路径会做什么，但作为一个简单的说明，这里是使用 Rcpp 添加 cumsum() 和 的 C++ 代码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>