从 R 调用 C 时出现堆损坏错误，无法找到源问题

Question

更新 3：问题已解决，但我将代码保留在这里以供将来参考——我在下面发布了一个答案以及代码的最终状态，以防人们想看到最终产品。

更新 2：重构为使用 R_alloc 而不是 calloc 进行自动清理。不幸的是，问题仍然存在。

更新：如果我在 UNPROTECT(1) 之前添加此行：

Rprintf("%p %p %p", (void *)rans, (void *)fm, (void *)corrs);

然后函数执行，没有损坏的堆错误。也许有一个后台垃圾收集调用在执行完成之前破坏了其中一个指针，从而导致写入垃圾指针？重要的是要注意，如果我没有打印出所有三个指针地址，错误就会回来。

此外，我运行在 M1 Mac 上安装它并通过 R CMD SHLIB 使用 clang 进行编译，以防万一苹果芯片是罪魁祸首。

我在尝试调试这个问题时束手无策，我想我应该向 SO 寻求帮助。我正在用 C 编写一个函数来优化我的 R 代码的某些部分，并且在运行多次调用该函数时出现堆损坏错误。使用 .Call("trimCovar", ...) 接口从 R 调用函数 trimCovar()。

由于以下几个原因，我在调试时遇到很多困难：

我在 OSX，所以我不能使用 Valgrind
C 函数依赖于 R 的输入，所以我无法自行调试 C 代码
仅当在 R 函数中多次调用该函数时才会发生堆损坏（只是直接运行宁.Call一堆没有错误）
错误点不一致

我从两组向量开始，然后将它们压缩成一个频率矩阵，其中每一列是向量集中的一个位置，每一行是出现的特定字符。我在传入之前将它们连接成一个矩阵，因为它使预处理更容易。频率矩阵的玩具示例是：

INPUT:             
  v1_1 = 101
  v1_2 = 011

  v2_1 = 111
  v2_2 = 110

Frequency Matrix:

position: | 1_1 | 1_2 | 1_3 | 2_1 | 2_2 | 2_3 |

0:          0.5   0.5   0.0   0.0   0.0   0.5
1:          0.5   0.5   1.0   1.0   1.0   0.5

目标是在向量集中找到 NV 最高相关位置，我通过计算位置的成对 KL 散度来做到这一点。这些存储在按升序排序的链表中，最后我采用与第一个 NV 条目对应的位置。我拥有的 R 代码可以解析其他所有内容，所以我真的只需要最后一个位置向量（允许重复）。

该函数有 5 个参数：

fMAT：频率矩阵（RObject，因此作为平面向量读入）
fSP：矩阵中对应于第一个向量集的位置的列
sSP：与 fSP 相同，但用于第二个矢量集
NV : return
NR：fMAT 中的列数

错误 returned 是：

R(95564,0x104858580) malloc: Heap corruption detected, free list is damaged at 0x600000f10040
*** Incorrect guard value: 4626885667169763328
R(95564,0x104858580) malloc: *** set a breakpoint in malloc_error_break to debug

只有当我运行一个调用此函数 10 次以上的 R 函数时才会发生这种情况，所以我假设我只是遗漏了一两个小的悬挂指针，破坏了内存引用。我已经尝试运行ning 这个与 gc() 在每次调用后立即在 R 中调用，但它没有解决问题。我不太确定此时还能做什么，我试过使用 lldb 但我不太确定如何使用该程序。从运行宁大量的打印语句我确定它通常在主循环中崩溃（在下面的代码中识别），但它在崩溃时不一致。我也尝试过保存错误的输入——我可以毫无问题地单独重新运行它们，所以它一定是一些相对较小的东西，只出现在许多运行中。

如果有帮助，我们很乐意提供更多详细信息。代码列在底部。

这里唯一分配的是链表节点，我想在 returning 之前我已经 free() 了它们。我还仔细检查了输入值，所以我 99.99% 确定我在 firstSeqPos、secondSeqPos、ans 或 [=30= 上的引用永远不会越界].我还对围绕此的 R 代码进行了三重检查，可以自信地说这不是此错误的根源。

我已经很长时间没有用 C 编写代码了，所以我觉得我遗漏了一些明显的东西。如果我真的需要，我可以尝试获取一个 Linux 框到运行 valgrind，但如果有其他选择，我会更喜欢它。提前致谢！

代码：

#include <R.h>
#include <Rdefines.h>
#include <Rinternals.h>
#include <math.h>
#include <stdlib.h>
#include <stdbool.h>

typedef struct node {
  double data;
  int i1;
  int i2;
  struct node *next;
} node;

// Linked list
// data is the correlation value,
// i1 the position from first vector set,
// i2 the position from second vector set
node *makeNewNode(double data, int i1, int i2){
  node *newNode; 
  newNode = (node *)R_alloc(1, sizeof(node));
  newNode->data = data;
  newNode->i1 = i1;
  newNode->i2 = i2;
  newNode->next = NULL;
  return(newNode);
}

//insert link in sorted order (ascending)
void insertSorted(node **head, node *toInsert, int maxSize) {
  int ctr = 0;
  if ((*head) == NULL || (*head)->data >= toInsert->data){
    toInsert->next = *head;
    *head = toInsert;
  } else {
    node *temp = *head;
    while (temp->next != NULL && temp->next->data < toInsert->data){
      temp = temp->next;
      if (ctr == maxSize){
        // Performance optimization, if we aren't inserting in the first NR
        // positions then we can just skip since we only care about the NR 
        // lowest scores overall
        return;
      }
      ctr += 1;
    }
    toInsert->next = temp->next;
    temp->next = toInsert;
  }
}

// MAIN FUNCTION CALLED FROM R
// (This is the one that crashes)
SEXP trimCovar(SEXP fMAT, SEXP fSP, SEXP sSP, SEXP NV, SEXP NR){
  // Converting input SEXPs into C-compatible values
  int nv = asInteger(NV);
  int nr = asInteger(NR);
  int sp1l = length(fSP);
  int sp2l = length(sSP);

  int *firstSeqPos = INTEGER(coerceVector(fSP, INTSXP));
  int *secondSeqPos = INTEGER(coerceVector(sSP, INTSXP));
  double *fm = REAL(fMAT);
  int colv1, colv2;

  // Using a linked list for efficient insert
  node *corrs = NULL;
  int cv1, cv2;
  double p1, p2, score=0;

  // USUALLY FAILS IN THIS LOOP
  for ( int i=0; i<sp1l; i++ ){
    cv1 = firstSeqPos[i];
    colv1 = (cv1 - 1) * nr;
    for ( int j=0; j<sp2l; j++ ){
      cv2 = secondSeqPos[j];
      colv2 = (cv2 - 1) * nr;

      // KL Divergence
      score = 0;
      for ( int k=0; k<nr; k++){
        p1 = fm[colv1 + k];
        p2 = fm[colv2 + k];
        if (p1 != 0 && p2 != 0){
          score += p1 * log(p1 / p2);
        }
      }
      // Add result into LL
      node *newNode = makeNewNode(score, cv1, cv2);
      insertSorted(&corrs, newNode, nv);
    }
    R_CheckUserInterrupt();
  }
  SEXP ans;
  PROTECT(ans = allocVector(INTSXP, 2*nv));
  int *rans = INTEGER(ans); 
  int ctr=0;
  int pos1, pos2;
  node *ptr = corrs;

  for ( int i=0; i<nv; i++){
    rans[2*i] = ptr->i1;
    rans[2*i+1] = ptr->i2;
    ptr = ptr->next;
  }

  UNPROTECT(1);
  return(ans);
}

Answer 1

int *firstSeqPos = INTEGER(coerceVector(fSP, INTSXP));
int *secondSeqPos = INTEGER(coerceVector(sSP, INTSXP));

这样不好。对 coerceVector() 的 2 次调用返回的 SEXP 需要保护。然而，通常认为更好的做法是在进入 .Call 入口点之前在 R 级别执行此强制转换。请注意，如果 fSP 和 sSP 是整数矩阵，则无需将它们强制转换为整数，因为它们已被视为 C 级别的整数向量。这也避免了可能昂贵的副本（R 中的 as.integer() 和 C 中的 coerceVector() both 触发矩阵数据的完整副本）。

Answer 2

上面已经回答了这个问题，但我收到了一些人发来的消息，要求提供最终代码，因此我将把它作为答案包含在内，以保留原始问题。这里有一些优化（感谢 @hpages 的帮助和故障排除）：

原始代码失败，因为 coerceVector() 的输出没有受到 PROTECT() 的保护。我重构了 R 代码以在调用此 C 函数之前检查整数输入以避免此函数调用并提高内存效率（有关更多详细信息，请参阅已接受的答案）。
原始代码使用R_alloc()，这让R负责在函数调用结束时清理内存。但是，这会在函数运行期间引入大量内存开销，因为分配给未插入链表的节点的内存在函数调用结束之前不会被清除。
分配calloc()并不像在函数末尾切换调用free()那么简单，因为我们必须防范用户中断程序执行的情况。如果在函数结束之前抛出中断信号，我们将永远不会释放内存。

最终 C 代码：

#include <R.h>
#include <Rdefines.h>
#include <Rinternals.h>
#include <math.h>
#include <stdlib.h>
#include <stdbool.h>

typedef struct node {
  double data;
  int i1;
  int i2;
  struct node *next;
} node;

// Defining the head as a static so that we can access it globally
// Important for ensuring clean up in case of interrupt
static node *corrs = NULL;

// Function to clean up memory allocations in case of interrupt
void cleanupFxn(){
  node *ptr = corrs;
  // Free allocated memory in linked list
  while (corrs != NULL){
    ptr = corrs;
    corrs = corrs->next;
    free(ptr);
  }
}

node *makeNewNode(double data, int i1, int i2){
  node *newNode; 
  // very important to use calloc here so we have control of when we free it
  // R_alloc() memory won't be freed until after function finishes execution
  newNode = (node *)calloc(1, sizeof(node));
  newNode->data = data;
  newNode->i1 = i1;
  newNode->i2 = i2;
  newNode->next = NULL;
  return(newNode);
}

// insert link in sorted order
// returns a bool corresponding to if we inserted
bool insertSorted(node **head, node *toInsert, int maxSize) {
  int ctr = 0;
  if ((*head) == NULL || (*head)->data >= toInsert->data){
    toInsert->next = *head;
    *head = toInsert;
    return(true);
  } else {
    node *temp = *head;
    while (temp->next != NULL && temp->next->data < toInsert->data){
      temp = temp->next;
      if (ctr == maxSize){
        // Performance optimization, if we aren't inserting in the first NR
        // positions then we can just skip since we only care about the NR 
        // lowest scores overall. Saves a huge amount of time and memory.
        return(false);
      }
      ctr += 1;
    }
    toInsert->next = temp->next;
    temp->next = toInsert;
    return(true);
  }
}

SEXP trimCovar(SEXP fMAT, SEXP fSP, SEXP sSP, SEXP NV, SEXP NR){
  // Converting inputs into C-compatible forms
  int nv = asInteger(NV);
  int nr = asInteger(NR);
  int sp1l = length(fSP);
  int sp2l = length(sSP);

  // Note here we're not using coerceVector() anymore
  // typechecking done on R side
  int *firstSeqPos = INTEGER(fSP);
  int *secondSeqPos = INTEGER(sSP);
  double *fm = REAL(fMAT);
  int colv1, colv2;
  
  // Using a linked list for efficient insert
  corrs = NULL;
  int cv1, cv2;
  double p1, p2, score=0;
  bool success;
  for ( int i=0; i<sp1l; i++ ){
    cv1 = firstSeqPos[i];
    colv1 = (cv1 - 1) * nr;
    for ( int j=0; j<sp2l; j++ ){
      cv2 = secondSeqPos[j];
      colv2 = (cv2 - 1) * nr;

      score = 0;
      for ( int k=0; k<nr; k++){
        p1 = fm[colv1 + k];
        p2 = fm[colv2 + k];
        if (p1 != 0 && p2 != 0){
          score += p1 * log(p1 / p2);
        }
      }
      node *newNode = makeNewNode(score, cv1, cv2);
      success = insertSorted(&corrs, newNode, nv);
      // If we don't insert, free the associated memory
      // I'm checking for NULL here just out of an abundance of caution
      if (!success && newNode != NULL){
        free(newNode);
        newNode = NULL;
      }
    }
    R_CheckUserInterrupt();
  }


  SEXP ans;
  PROTECT(ans = allocVector(INTSXP, 2*nv));
  int *rans = INTEGER(ans); 
  node *ptr=corrs;

  for ( int i=0; i<nv; i++){
    rans[2*i] = ptr->i1;
    rans[2*i+1] = ptr->i2;
    ptr = ptr->next;
  }

  // Free allocated memory in linked list
  cleanupFxn();

  UNPROTECT(1);
  return(ans);
}

假设 C 文件名为 trimCovar.c，我们将使用 R CMD SHLIB trimCovar.c.

进行编译

R 运行这个函数的代码：

dyn.load("trimCovar.so")

# Wrapped into a function with on.exit(...) to ensure cleanup
# in the event the user or system interrupts execution early
CorrComp_C <- function(fm, fsp, ssp, nv, nr){
  # type checking to ensure input to C is integer vector
  # (could probably do more type checking here, mainly for illustration)
  stopifnot(is(fsp, 'integer'))
  stopifnot(is(ssp, 'integer'))
  on.exit(.C("cleanupFxn"))
  a <- .Call('trimCovar', fm, fsp, ssp, nv, nr)
  return(a)
}

从 R 调用 C 时出现堆损坏错误，无法找到源问题

Heap Corruption error when calling C from R, can't find the source issue

c

r