Matlab и c ++ матричный продукт: один и тот же вход, разные выходы

Я использую Matlab 2010 в моем Ubuntu 12.04 на компьютере с архитектурой x86 и g ++ 4.6.3.
Вот как я делаю производство и вклады:

   #include <Src/Tools/Math/Matrix_nxn.h>
#include <iostream>
using namespace std;
int main()
{
Matrix_nxn<double,4> A1,A2,Tb,aa;
A1[0][0] = 0.99958087959447828;   A1[0][1] = 1.7725781974830023e-18;A1[0][2] = 0.028949354900049871;  A1[0][3] = 0;
A1[1][0] = -0.028949354900049871; A1[1][1] = 6.1204654815537932e-17;A1[1][2] = 0.99958087959447828;   A1[1][3] = 0;
A1[2][0] = 0,           A1[2][1] = -1;            A1[2][2] = 6.1230317691118863e-17;A1[2][3] = 0.21129000000000001;
A1[3][0] = 0,           A1[3][1] = 0;         A1[3][2] = 0;             A1[3][3] = 1;

A2[0][0] = 0.90634806393366396;   A2[0][1] = -0.42253187690835708;A2[0][2] = 0;A2[0][3] = 0;
A2[1][0] = 0.42253187690835708;   A2[1][1] = 0.90634806393366396; A2[1][2] = 0;A2[1][3] = 0;
A2[2][0] = 0;           A2[2][1] = 0;           A2[2][2] = 1;A2[2][3] = 0;
A2[3][0] = 0;           A2[3][1] = 0;           A2[3][2] = 0;A2[3][3] = 1;

Tb[0][0] = 0.99956387949834924; Tb[0][1] = -0.00016363183229951183;   Tb[0][2] = -0.029530052943282908; Tb[0][3] = 0;
Tb[1][0] = 0;         Tb[1][1] = 0.99998464792303143; Tb[1][2] = -0.0055411116439683869;Tb[1][3] = 0;
Tb[2][0] = 0.029530506297888514;Tb[2][1] = 0.0055386950515785164; Tb[2][2] = 0.99954853411673616;   Tb[2][3] = 0;
Tb[3][0] = 0;         Tb[3][1] = 0;           Tb[3][2] = 0;             Tb[3][3] = 1;aa = Tb*A1*A2;
cout.precision(25);
cout <<aa[0][0]<<'   '<<aa[0][1]<<'  '<<aa[0][2]<<'  '<<aa[0][3]<<endl
<<aa[1][0]<<'   '<<aa[1][1]<<'  '<<aa[1][2]<<'  '<<aa[1][3]<<endl
<<aa[2][0]<<'   '<<aa[2][1]<<'  '<<aa[2][2]<<'  '<<aa[2][3]<<endl
<<aa[3][0]<<'   '<<aa[3][1]<<'  '<<aa[3][2]<<'  '<<aa[3][3]<<endl;
}

и это определение operator*:

Matrix_nxn<T, N> res;
size_t i, j, k;
for (i = 0; i < N; ++i)
{
for (j = 0; j < N; ++j)
{
for (k = 0; k < N; ++k)
{
res[i][j] += m1[i][k] * m2[k][j];
}
if (MVTools::isNearInf(res[i][j]))
{
if (MVTools::isNearPosInf(res[i][j]))
throw MVException(MVException::PosInfValue);
else
throw MVException(MVException::NegInfValue);
}
}
}
return res;

Странно то, что я делаю одинаковые матрицы с одинаковыми значениями внутри Matlab и получаю разные результаты. Вот код Matlab:

Tb = [0.99956387949834924,-0.00016363183229951183,-0.029530052943282908,0;0,0.99998464792303143,-0.0055411116439683869,0;0.029530506297888514,0.0055386950515785164,0.99954853411673616,0;0,0,0,1];
A1 = [0.99958087959447828,1.7725781974830023e-18,0.028949354900049871,0;-0.028949354900049871,6.1204654815537932e-17,0.99958087959447828,0;0,-1,6.1230317691118863e-17,0.21129000000000001;0,0,0,1];
A2 = [0.90634806393366396,-0.42253187690835708,0,0;0.42253187690835708,0.90634806393366396,0,0;0,0,1,0;0,0,0,1];
aa = Tb*A1*A2;
aa - aaa

ans =

1.0e-16 *

0  -0.555111512312578                   0                   0
0                   0                   0                   0
0                   0                   0                   0
0                   0                   0                   0

в то время как aaa является выходом реализации c ++.
Я знаю, что ошибки так мало, но я хочу знать, что вызывает проблему! Я хочу отладить много кода, и мне нужна нулевая разница для хорошей отладки.

0

Решение

Причина различного значения (каким бы незначительным оно ни было) в том, что алгоритмы, используемые matlab, и вы не тот же самый.

Ваш алгоритм прост O(N^3) умножение матриц. Существуют специальные алгоритмы для матриц небольшого размера, которые должны быть эффективно рассчитаны, а также сложные алгоритмы, которые имеют лучшую асимптотику, чем O(N^3),

Если вам интересно, смотрите:

2

Другие решения

Я вижу, вы ожидаете 25-значную точность от кода C ++. Это очень маловероятно, используя double тип. Вы можете получить лучшую точность, используя long double но, возможно, не так много, как 25 цифр.

Увидеть: Какова точность long double в C ++?

1

По вопросам рекламы [email protected]