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

Question

此文件位于 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);

Answer 1

fulltosparse.c 中似乎存在错误计算新的nzmax当稀疏矩阵需要额外的空间来存储非零时。

在每个非零行 i和列 j的 m行和 n列，分别属于 full源矩阵(k th 非零元素)，进行检查(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 的情况。 :

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

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

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