C 函数的 R-Wrapper 无法正常工作

Question

目前我正在向自己介绍 C 以扩展我的 R 函数。我用 C 编写了一个函数用于一些计算，它们工作正常。但是，一旦我在 R 本身中编写了一个包装器，似乎就会出现一些错误。将我的 C 函数视为 "colV"，将 "abc" 视为某个任意矩阵。

语句 (R) .Call("colV", abc, ncol(abc), nrow(abc)) 完全可以正常工作（每次，无论我使用它的频率如何），而

colV = function(x){
  nc = ncol(x)
  nr = nrow(x)
  .Call("colV", x, nc, nr)
}

第三次使用时给出错误结果：

> colV(abc)
[1] 1.4274933 0.6254796 0.6774042 1.7094617 0.7250386
> colV(abc)
[1] 1.4274933 0.6254796 0.6774042 1.7094617 0.7250386
> colV(abc)
[1] -8.370087e+22  6.254796e-01  6.774042e-01  1.709462e+00  7.250386e-01

如果我再次声明包装器，前两次运行正常，第三次尝试时出现相同的结果。请注意，除了外观之外，实际上只有第一个值发生了变化！

有人知道包装有什么问题吗？如前所述，如果我只使用 .Call 语句，则始终会返回正确的结果：

> .Call("colV", abc, ncol(abc), nrow(abc))
[1] 1.4274933 0.6254796 0.6774042 1.7094617 0.7250386
> .Call("colV", abc, ncol(abc), nrow(abc))
[1] 1.4274933 0.6254796 0.6774042 1.7094617 0.7250386
> .Call("colV", abc, ncol(abc), nrow(abc))
[1] 1.4274933 0.6254796 0.6774042 1.7094617 0.7250386
> .Call("colV", abc, ncol(abc), nrow(abc))
[1] 1.4274933 0.6254796 0.6774042 1.7094617 0.7250386
> .Call("colV", abc, ncol(abc), nrow(abc))
[1] 1.4274933 0.6254796 0.6774042 1.7094617 0.7250386
> .Call("colV", abc, ncol(abc), nrow(abc))
[1] 1.4274933 0.6254796 0.6774042 1.7094617 0.7250386

此外，我在使用 .C 接口时没有遇到过这个问题，但是由于速度和内存的原因，这对我来说不是解决方案。提前致谢艾琳

编辑：这是 C 代码：

#include <R.h>
#include <Rinternals.h>
#include <Rmath.h>


SEXP colV(SEXP y, SEXP n, SEXP r){
    int *nc = INTEGER(n);
    double *x = REAL(y);
    int d = length(y);
    int *nr = INTEGER(r);
    int i, j, z;
    //int d = nr * nc;
    double colMean[(*nc)];
    double xSq[(d)];
    double colMsq[(*nc)];
    double xSm[(*nc)];
    SEXP result;
    PROTECT(result = allocVector(REALSXP, (*nc)));
    memset(REAL(result), 0, (*nc) * sizeof(double));
    double *colVar = REAL(result);
    //PROTECT(colVar = NEW_DOUBLE(nc));
    int fr = ((*nr) - 1);


    for(z = 0; z < (d); z++){
        xSq[z] = x[z] * x[z];
    }

    for(i = 0; i < (*nc); i++){
        for(j = 0; j < (*nr); j++){
            colMean[i] += (x[(j + ((*nr) * i)) ]);
            xSm[i] += (xSq[(j + (*nr * i))]);
        }
        colMean[i] = (colMean[i] / (*nr));
        colMsq[i] = (*nr) * (colMean[i] * colMean[i]);
        colVar[i] = ((xSm[i] - colMsq[i]) / fr);
    }
    UNPROTECT(1);
    return(result);
}

Answer 1

编辑：正如评论中所说，这（下面）有效，但是对于另一个不同的 c-Function 会出现同样的问题。使用 .Call 直接提供预期的结果，而将 .Call 放在包装器中会出现错误。我还没有找到原因。

我现在找到了解决问题的方法，但是我不明白为什么会这样。如果有人可以解释这一点，请随意这样做！

改变循环中 colMean 部分的计算改变了一切：

...
for(i = 0; i < (*nc); i++){
        for(j = 0; j < (*nr); j++){
            colMean[i] += ((x[(j + ((*nr) * i)) ]) / (*nr));
            xSm[i] += (xSq[(j + (*nr * i))]);
        }
        //colMean[i] = (colMean[i] / (*nr));
        colMsq[i] = (*nr) * (colMean[i] * colMean[i]);
        colVar[i] = ((xSm[i] - colMsq[i]) / fr);
    }
...

因此将 nr 除法直接移动到加法语句中解决了问题。

祝大家day/evening愉快！

C 函数的 R-Wrapper 无法正常工作

R-Wrapper for C-Function does not work properly

r

c

function

call

wrapper