Precision error on matrix multiplication
- by Wam
Hello all,
Coding a matrix multiplication in my program, I get precision errors (inaccurate results for large matrices).
Here's my code. The current object has data stored in a flattened array, row after row. Other matrix B has data stored in a flattened array, column after column (so I can use pointer arithmetic).
protected double[,] multiply (IMatrix B)
{
int columns = B.columns;
int rows = Rows;
int size = Columns;
double[,] result = new double[rows,columns];
for (int row = 0; row < rows; row++)
{
for (int col = 0; col < columns; col++)
{
unsafe
{
fixed (float* ptrThis = data)
fixed (float* ptrB = B.Data)
{
float* mePtr = ptrThis + row*rows;
float* bPtr = ptrB + col*columns;
double value = 0.0;
for (int i = 0; i < size; i++)
{
value += *(mePtr++) * *(bPtr++);
}
result[row, col] = value;
}
}
}
}
}
Actually, the code is a bit more complicated : I do the multiply thing for several chunks (so instead of having i from 0 to size, I go from localStart to localStop), then sum up the resulting matrices.
My problem : for a big matrix I get precision error :
NUnit.Framework.AssertionException: Error at (0,1)
expected: <6.4209571409444209E+18>
but was: <6.4207619776304906E+18>
Any idea ?
