Float and double precision in c#

Having a bit of a problem here on a matrix multiplication code. I seem to lose precision on large matrices multiplications (my code runs fine on small matrices).

My loop is the following :

for (int j = 0; j < columns; j++)
{
    float[] column = otherMatrix.Column(j);
    for (int i = 0; i < rows; i++)
    {
        double s = 0;
        for (int k = 0; k < size; k++)
            s += this[i,k] * ((double) column[k]);
        result[i, j] = (float)s;
    }
}

As you can see, I force a (double) precision to make sure I don't lose precision when multiplying my two floats.

Looking at IL code, I can see two conv.r8 which make me think that IL code has this float-to-double precision conversion in it.

However, when running it and having a look at the disassembly (x86 machine), I see the following :

0000024e  fld         dword ptr [edx+eax*4+8] 
00000252  fmulp       st(1),st 
00000254  fadd        qword ptr [ebp-64h] 
00000257  fstp        qword ptr [ebp-20h] 

It make开发者_StackOverflows me think that JIT has thought that since I'm already multiplying floats, it shouldn't use double precision multiplication but single precision multiplication, giving me the errors I've been tracking.

Am I right ? Is there any way to force this double precision multiplication ?

Thanks

I think you're misinterpreting the assembly. I believe that FMULP always operates on the 80-bit registers. I would be surprised to see the JIT doing the wrong thing here.

I suggest you use my DoubleConverter to write out the precise values before and after the arithmetic. That way you should get a better idea of what's going on.

You might want to switch to decimal for better precision

First, is there any reason you can't use decimal? Decimals maintain the necessary precision of any number, and don't have the annoying "nearest binary representation of the mantissa" problems floats have that can result in x-y = z +- .0000000000...0000001 errors.

From MSDN:

The Decimal value type represents decimal numbers ranging from positive 79,228,162,514,264,337,593,543,950,335 to negative 79,228,162,514,264,337,593,543,950,335. The Decimal value type is appropriate for financial calculations requiring large numbers of significant integral and fractional digits and no round-off errors. The Decimal type does not eliminate the need for rounding. Rather, it minimizes errors due to rounding.

Failing this, try casting both sides to double before multiplying. float * double may result in a double, but as your double is a float in disguise being compared to another float, the compiler may be ignoring your desired precision, "knowing" that two floats won't make a double.