用Matlab为C++编写矩阵

Question

我有一个dim*dim的稀疏矩阵P作为指针通过

double* P

/* create the output matrix */
plhs[0] = mxCreateDoubleMatrix(dim,dim,mxREAL);

/* get a pointer to the real data in the output matrix*/
P = mxGetPr(plhs[0]);

我在一个mex文件中这样做,因为我需要很多for循环来填充P和c ++比matlab要快得多.

目前,dim = 22500并且c ++填充P需要大约2秒(matlab此任务花费50秒),并且在matlab中将矩阵标准化大约100秒,并且再次100秒来擦除matlab中的所有零列.我在matlab中使用以下代码执行此操作:

for i=1:size(P,1)
    if sum(P(i,:)) > 0
        sum(P(i,:))
        P(i,:)=(1/sum(P(i,:))).*P(i,:);
    end
end

% clear empty rows and colunms
P(~any(P,2),:)=[];
P(:,~any(P))=[];

我现在的问题是:我可以用c ++做到这一点吗？我尝试用以下方式在c ++中规范化P:

int i;
int j;
int sum;
int get_idx(int x, int y, int rows) {
   return x +y * rows;   
}
/* NORMALIZE */
for(i = 0; i <dim; i++)  {
    sum=0;
    for(j=0; j<dim;j++) {   
        sum = sum + P[get_idx(i,j,dim)];
    }
    if(sum > 0) {
         for(j=0; j<dim;j++) {   
           P[get_idx(i,j,p_rows)]=P[get_idx(i,j,dim)]*(1/sum);
        }
    }
}

但由于某些原因,这段代码似乎没有改变P,而且这在c ++中也需要大约85秒.还有更快的方式吗？此外,是否可以清除空行和列？

Answer 1

为什么选择C++？

在规范化之前清除空行/列 - 您不需要规范化空条目.

矢量化规范化:

s = sum(P, 2);
valid = s > 0;
P( valid,: ) = bsxfun(@rdivide, P(valid,:), s(valid) );

当当!

bsxfun非常有趣!

更新:关于行/列的减少.
经过短暂的调查后,我认为有一个~x3速度因素可以获得:

考虑以下三个选项:

1. P( ~any(P,2), :) = []; P( :, ~any(P,1) ) = [];
2. P( :, ~any(P,1) ) = []; P( ~any(P,2), :) = [];
3. P = P( any(P,2), any(P,1) );

测试这三种替代方案,你会发现第三种方法的速度要快〜x3,而第一种方式比第二种方式慢一点(但始终如一).

为什么？
如果你还记得,Matlab以列优先的方式将matices存储在内存中,因此在行之前消除了一些行,从而节省了一些复制和重新分配内存.

然而,第一个和第二个备选方案复制并重新分配内存两次:一次用于行,一次用于列,而第三种方法只用内存一次!