11

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

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

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

如何创建它，以及如何附加到它？尤其是我怎样才能将它返回给 R？

感谢帮助。

4

2 回答 2

6

你想用什么作为临时结构真的取决于你，因为最终你必须为结果分配一个向量。因此，无论您将使用什么都不是您将返回的。有几种可能性：

  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;
 }
于 2012-01-10T00:46:57.293 回答
5

如果您愿意从 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>
于 2012-01-10T01:01:04.503 回答