-1

I don't understand why the following code gives the wrong answer for Matrix 'u'. Matrix 'u' should equal the identity Matrix but only some of the values are correct. Can anyone help me understand why this might be happening?

for (k=0; k<3; k++) {

    int j;
    for (j=0; j<3; j++) {

        int h;
        for (h=0; h<3; h++) {


            u[k][j]+=(F[k][h])*(B[h][j]);

        }

    }

}

Matrix F
2 -2.2 0.6 
-0 0.4 -0.2 
-3 3.2 -0.6  
Matrix B
2 3 1 
3 3 2 
6 1 4 
Matrix u
1 -4.44089e-16 0 
0 1 0 
8.88178e-16 1.33227e-15 1 
  • 1
    Interesting "*code is correct but answer is wrong*"; then what do you think causes the wrong output? – abhishek_naik Jul 15 '16 at 16:31
  • I won't say that numbers like 8.88178e-16 are the same as zero. But they are really really tiny. So it may be that errors due to arithmetic precision may be in play – infixed Jul 15 '16 at 16:34
  • Welcome to Stack Overflow! Please **[edit]** your question with a [mcve] or [SSCCE (Short, Self Contained, Correct Example)](http://sscce.org) – NathanOliver Jul 15 '16 at 16:34
  • It's an accuracy issue. In a nutshell your decimals are being approximated. The easiest fix is print you matrices to fewer decimal places. Use %1.5f instead of %f. – Stephen Quan Jul 15 '16 at 16:35
  • Thanks, I'll look into that. – Richard Fisher Jul 15 '16 at 16:46
  • http://stackoverflow.com/questions/9577179/c-floating-point-precision addresses the same issue, but I'm not sure its a duplicate enough – infixed Jul 15 '16 at 16:51

2 Answers2

0

The reason your math isn't working out the way you might expect is because you are working with floating point numbers (Floating point numbers) and not real "Real" numbers. If you wish to do precise work with floats or doubles, then you need to understand floating point math.

If you really don't care about precision and just want to get stuff to look right, then a quick fix is to define a numerical zero and set everything smaller than it to actual zero. Something like 1e-6 could work. Some example code is the following.

for (int i = 0; i < 3; i++) 
   for (int j = 0; j < 3; j++) 
      u[i][j] = abs(u[i][j]) < 1e-6 ? 0 : u[i][j];
tlawren
  • 38
  • 4
  • That just "fixes" 0.0 but 1.0 has similar issues. And so does 0.573 actually, or any other number. You can't guess. – MSalters Jul 15 '16 at 20:45
0

Initialise u[k][j]=0 when you move on to calculate sum for different element of u. This should work fine.

for (k=0; k<3; k++)
{
    int j;
    for (j=0; j<3; j++) 
    {
        int u[k][j]=0;
        int h;
        for (h=0; h<3; h++)
        {
            u[k][j]+=(F[k][h])*(B[h][j]);
        }
    }
}
sagar_jeevan
  • 761
  • 5
  • 16
  • 34