使用OpenCV模拟Matlab的mldivide

Question

使用OpenCV模拟Matlab的mldivide

3

那个问题得到了mrdivide的解决。但是现在我遇到了mldivide的问题，它当前的实现方式如下：

cv::Mat mldivide(const cv::Mat& A, const cv::Mat& B ) 
{
    //return  b * A.inv();
    cv::Mat a;
    cv::Mat b;
    A.convertTo( a, CV_64FC1 );
    B.convertTo( b, CV_64FC1 );

    cv::Mat ret;
    cv::solve( a, b, ret, cv::DECOMP_NORMAL );

    cv::Mat ret2;
    ret.convertTo( ret2, A.type() );
    return ret2;
}

根据我的理解，mrdivide能够正常工作应该意味着mldivide也能正常工作，但是我无法得到与Matlab相同的结果。再次强调，这次我试图对非方阵[19x19] \ [19x200]进行操作。

- Goz

可能是重复的问题，参考如何在OpenCV中实现矩阵左除法。 - beaker

@beaker 如果你看一下，你会发现我已经“实现”了那个链接中解释的方法。 - Goz

那么你就需要向我们展示问题所在。我将删除重复的投票。 - beaker

@beaker..我不确定该如何向您展示。上述答案给出的结果与Matlab使用相同数据给出的结果完全不同。 - Goz

你已经用每个结果计算了 A*x = B 并发现其中一个是错误的吗？ - beaker

1个回答

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Amro · Accepted Answer

如我之前在你的另一个问题中提到的那样，我正在使用MATLAB以及mexopencv，这样我就可以轻松比较MATLAB和OpenCV的输出结果。

话虽如此，我无法重现你的问题：我生成了随机矩阵，并重复了N=100次比较。我正在运行MATLAB R2015a，mexopencv编译的是OpenCV 3.0.0：

N = 100;
r = zeros(N,2);
d = zeros(N,1);
for i=1:N
    % double precision, i.e CV_64F
    A = randn(19,19);
    B = randn(19,200);

    x1 = A\B;
    x2 = cv.solve(A,B);   % this a MEX function that calls cv::solve

    r(i,:) = [norm(A*x1-B), norm(A*x2-B)];
    d(i) = norm(x1-x2);
end

所有结果均一致，误差非常小，约为1e-11：

>> mean(r)
ans =
   1.0e-12 *
    0.2282    0.2698

>> mean(d)
ans =
   6.5457e-12

(顺便说一句，我也尝试过x2 = cv.solve(A,B,'IsNormal'，true);，它设置了cv::DECOMP_NORMAL标志，但结果也没有很大的区别)。

这使我相信，要么您的矩阵恰好突出了OpenCV求解器中的某个边缘情况，导致它未能给出正确的解决方案，要么更可能是您的代码中存在其他bug。

我建议先仔细检查数据的加载方式，特别是要注意矩阵的布局（显然MATLAB是列优先，而OpenCV是行优先）...

此外，您从未告诉我们有关您的矩阵的任何信息; 它们是否表现出某种特征，是否有任何对称性，它们是否大多数为零，它们的秩等。

在OpenCV中，默认求解器方法是LU分解，如果合适的话，您必须明确地更改它。另一方面，MATLAB将自动选择最适合矩阵A的方法，而LU只是可能的分解之一。

编辑（回应评论）

在MATLAB中使用SVD分解时，左右特征向量U和V的符号是任意的（这实际上来自DGESVD LAPACK例程）。为了获得一致的结果，一个惯例是要求每个特征向量的第一个元素具有某种符号，并根据需要将每个向量乘以+1或-1来翻转符号。我还建议查看eigenshuffle。

再次，这是我进行的测试，以确认我在MATLAB和OpenCV中获取SVD相似的结果：

N = 100;
r = zeros(N,2);
d = zeros(N,3);
for i=1:N
    % double precision, i.e CV_64F
    A = rand(19);

    % compute SVD in MATLAB, and apply sign convention
    [U1,S1,V1] = svd(A);
    sn = sign(U1(1,:));
    U1 = bsxfun(@times, sn, U1);
    V1 = bsxfun(@times, sn, V1);
    r(i,1) = norm(U1*S1*V1' - A);

    % compute SVD in OpenCV, and apply sign convention
    [S2,U2,V2] = cv.SVD.Compute(A);
    S2 = diag(S2);
    sn = sign(U2(1,:));
    U2 = bsxfun(@times, sn, U2);
    V2 = bsxfun(@times, sn', V2)';  % Note: V2 was transposed w.r.t V1
    r(i,2) = norm(U2*S2*V2' - A);

    % compare
    d(i,:) = [norm(V1-V2), norm(U1-U2), norm(S1-S2)];
end

再次测试结果非常相似，错误接近机器精度并且可以忽略不计：

>> mean(r)
ans =
   1.0e-13 *
    0.3381    0.1215

>> mean(d)
ans =
   1.0e-13 *
    0.3113    0.3009    0.0578

在OpenCV中有一件事我不太确定，但是MATLAB的svd函数返回的奇异值是按降序排序的（与eig函数不同），且特征向量的列也以相应顺序排列。

现在如果由于某些原因OpenCV中的奇异值不能保证排序，那么如果要将结果与MATLAB进行比较，就必须手动排序，如下所示：

% not needed in MATLAB
[U,S,V] = svd(A);
[S, ord] = sort(diag(S), 'descend');
S = diag(S);
U = U(:,ord)
V = V(:,ord);