c - 如何解决matlab示例 "fulltosparse.c"的错误

该文件位于 C:\Program Files\MATLAB\R2013a\extern\examples\refbook 中。在 mex 之后，我使用了:

aa = [1 2 3 ; 4 5 6]
fulltosparse(aa)

第一次，该命令可能有效。但是多试几次 fulltosparse(aa)。你会发现它会崩溃。谁能告诉我为什么？

mex -largeArrayDims fulltosparse.c
aa = [1 2 3; 4  5 7];
fulltosparse(aa);
fulltosparse(aa);
fulltosparse(aa);
fulltosparse(aa);
fulltosparse(aa);
fulltosparse(aa);

最佳答案

当稀疏矩阵需要额外的空间来存储非零值时，在计算新的 fulltosparse.c 时，nzmax 中出现了一个错误。

在 i 源矩阵(第 j 个非零元素)的 m 行和 n 列的每个非零行 full 和列 k 处，分别进行检查( k>=nzmax )以确保在稀疏矩阵数据缓冲区( sr 、 si 和 irs )。如果没有足够的空间容纳更多元素，缓冲区将通过 mxRealloc 扩展，并为增加的 nzmax 非零元素提供足够的空间。

问题在于它如何计算新的 nzmax :

mwSize oldnzmax = nzmax;
percent_sparse += 0.1;
nzmax = (mwSize)ceil((double)m*(double)n*percent_sparse);

/* make sure nzmax increases atleast by 1 */
if (oldnzmax == nzmax)
    nzmax++;

该函数以初始 percent_sparse = 0.2 开始，对于 2x3 完整矩阵 aa 对应于 nzmax = ceil(6*0.2) = 2 ，并开始循环遍历行(最快)和列。这是出了什么问题:

在 k=2(MATLAB 索引中的第三个元素； i=0 , j=1 )，它需要第一次扩展缓冲区，上面的代码在重新分配之前运行。 oldnzmax 为 2。 percent_sparse 增加到 0.3，给出 nzmax=ceil(6*0.3)=2 。由于 oldnzmax == nzmax ，它只是增加( nzmax++ )并重新分配给 nzmax=3 。行。

在 k=3 (第四个元素； i=1 , j=1 )，它经过类似的路径，将 percent_sparse 增加到 0.4，计算 nzmax=ceil(6*0.4)=3 并将 nzmax 增加到 4。

当它到达 k=4 (第五个元素； i=0 , j=2 )时，循环迭代从 nzmax=4 和 percent_sparse=0.4 开始。至此，很明显 percent_sparse 没有跟上: (nzmax=4)/6 = 0.666 > percent_sparse=0.4 。现在，当 percent_sparse += 0.1; 只给出 0.5 并且 给出新的 nzmax=ceil(6*0.5)=3 时，这个错误就变得很明显了，它小于 oldnzmax=4 ! 。

灾难:使用 k=4 和 nzmax=3 ，它重新分配稀疏矩阵缓冲区( sr 、 si 和 irs )，并溢出两个缓冲区:

sr[k] = pr[i];  // k=4, sr is length 3
irs[k] = i;     //      irs is length 3

缓冲区的大小实际上减少了。但是，即使将测试更改为 oldnzmax >= nzmax ，缓冲区的大小仍然不会增加，因为 nzmax 是根据 percent_sparse 计算得出的，而 nzmax<oldnzmax 增加的速度不够快 。我认为有两个变化是必要的。首先，测试和增量都需要处理 percent_sparse 时的情况:

if (oldnzmax >= nzmax)
    nzmax = oldnzmax+1;

其次，为了更好的衡量，当条件为真并且 nzmax 递增而不是从 percent_sparse 计算时，应该正确更新 ojit_code :

if (oldnzmax >= nzmax)
{
    nzmax = oldnzmax+1;
    percent_sparse = (double)nzmax/((double)m*(double)n);
}

关于c - 如何解决matlab示例 "fulltosparse.c"的错误，我们在Stack Overflow上找到一个类似的问题：https://stackoverflow.com/questions/20530782/