犰狳：矩阵乘法损失的精度

Question

我在执行矩阵乘法时似乎遇到了精度损失，并且想知道如何防止这种情况。例如，假设壮举和β是适当的尺寸，

Y = feat*beta.rows(0,N); 

我一起工作的数字是相当小的值，多数数字是小于1E-3所以有可能是什么我想实现是不可能的。我还应该注意到这是一个调用C++函数的MATLAB函数，因此涉及到MEX编译器。当他们到达时，我确实检查了mex函数中的数字，他们是正确的，只有在上面的这一行后，我才得到令人难以置信的错误答案。

编辑：我决定这不会给程序的完整背景带来伤害。这是我到目前为止。我有一条线，精度的损失来自一个评论。

编辑2：下面是一些有关矩阵的例子。
Feat_2d是5x4608

 0   0   0   0   0 
    0   0   0   0   0 
    0   0   0   0 4.8146 
    0   0 19.0266   0   0 
    0   0   0   0   0 

Beta_2d是4609x4。我把最后一排的壮举*的乘法Beta_2d

-7.1486e-05 -1.6801e-04 1.0970e-05 3.7837e-04 
-8.7524e-05 1.8275e-04 -6.7857e-04 2.6267e-04 
-9.1812e-05 -6.5495e-05 -1.7687e-03 -3.2168e-04 
    0e+00  0e+00  0e+00  0e+00 
-4.5089e-04 -5.6013e-05 1.4841e-04 2.4912e-04 

Y =

6.8995e-310  0e+00 4.7430e-322 1.7802e-306 
6.8995e-310  0e+00 4.9407e-322 1.4463e-307 
0e+00  0e+00  0e+00 1.4463e-307 
0e+00  0e+00 1.2352e-322 1.2016e-306 
6.8996e-310 6.8996e-310 6.8995e-310 1.7802e-306 

下面是从EDIT1

#include <mex.h> 
#include <iostream> 
#include <armadillo> 

using namespace arma; 

void predict_bbox_reg(double *beta, int beta_dim[2], double *t_inv, int tinv_dim[2], double mu, double *feat, int feat_dim[2], double *boxes, int box_dim[2]){ 


    //convert pointers 
    //beta 
    arma::mat beta_2d = mat(beta_dim[0], beta_dim[1]); 
    for(int i = 0; i<beta_dim[0]; i++){ 
     for(int j = 0; j<beta_dim[1]; j++){ 
      beta_2d(i, j) = beta[i+beta_dim[0]*j]; 
     } 
    } 

    //t_inv 
    arma::mat tinv_2d = mat(tinv_dim[0], tinv_dim[1]); 
    for(int i = 0; i<tinv_dim[0]; i++){ 
     for(int j = 0; j<tinv_dim[1]; j++){ 
      tinv_2d(i, j) = t_inv[i+tinv_dim[0]*j]; 
     } 
    } 
    //feadoublet_2d 
    arma::mat feat_2d = mat(feat_dim[0], feat_dim[1]); 
    for(int i = 0; i<feat_dim[0]; i++){ 
     for(int j = 0; j<feat_dim[1]; j++){ 
      feat_2d(i, j) = feat[i+feat_dim[0]*j]; 
     } 
    } 

    //boxes 
    arma::mat box_2d = mat(box_dim[0], box_dim[1]); 
    for(int i = 0; i<box_dim[0]; i++){ 
     for(int j = 0; j<box_dim[1]; j++){ 
      box_2d(i, j) = boxes[i+box_dim[0]*j]; 
     } 
    } 


    arma::mat Y = mat(feat_dim[0], beta_dim[1]); 

    Y = feat_2d*beta_2d.rows(0,beta_dim[0]-2);// this is the precision loss 

    arma::mat y1 = beta_2d.row(beta_2d.n_rows-1); 
    Y.each_row() += y1.row(0); 

    //RETURNS SOMETHING 


}   

void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]) { 

    int M = mxGetM(prhs[0]); 
    int N = mxGetN(prhs[0]); 

    int beta_dim[2] = {M,N}; 
    double *beta = mxGetPr(prhs[0]); 


    M = mxGetM(prhs[1]); 
    N = mxGetN(prhs[1]); 
    int tinv_dim[2] = {M,N}; 
    double *t_inv = mxGetPr(prhs[1]); 


    double mu = *mxGetPr(prhs[2]);  

    M = mxGetM(prhs[3]); 
    N = mxGetN(prhs[3]); 

    int feat_dim[2] = {M,N}; 
    double *feat = mxGetPr(prhs[3]); 

    M = mxGetM(prhs[4]); 
    N = mxGetN(prhs[4]); 

    int box_dim[2] = {M,N}; 
    double *ex_boxes = mxGetPr(prhs[4]); 

    predict_bbox_reg(beta, beta_dim, t_inv, tinv_dim, 
      mu, feat, feat_dim, ex_boxes, box_dim); 

    //RETURNS results to matlab 
} 

Answer 1

的代码，我看不到任何直接失误你的代码，但我怀疑犰狳有精度问题。我怀疑这可能是您指向无用数据转换的东西。我会使用犰狳在.../mex_interface/armaMex.hpp文件中提供的功能。

矩阵乘法（mult_test.cpp）的一个简单的例子：

#include <armadillo> 
#include "armaMex.hpp" 

void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]) 
{  
    // Convert to Armadillo 
    mat x1 = conv_to<mat>::from(armaGetPr(prhs[0],true)); 
    mat x2 = conv_to<mat>::from(armaGetPr(prhs[1],true)); 
    mat y(x1.n_rows,x2.n_cols); 

    // Do your stuff here: 
    y = x1*x2; 

    // Convert back to Matlab 
    plhs[0] = armaCreateMxMatrix(y.n_rows, y.n_cols, mxDOUBLE_CLASS, mxREAL); 
    armaSetPr(plhs[0], conv_to<mat>::from(y)); 
    return; 
} 

与.m文件

X1 = [1 2 3 ; 4 5 6]; 
X2 = 1e-7*rand(3,2); 

Y = mult_test(X1,X2); 

disp('Matlab:') 
disp(X1*X2) 
disp('Armadillo:') 
disp(Y) 

给输出

Matlab: 
1.0e-06 * 

0.240798243020273 0.410716970485213 
0.559953800808707 0.974915983937571 

Armadillo: 
1.0e-06 * 

0.240798243020273 0.410716970485213 
0.559953800808707 0.974915983937571