Why does this code return two different values doing the same?

Question

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

int main()
{
    double a;
    double b;
    double q0 = 0.5 * M_PI + 0.5 * -2.1500000405000002;
    double q1 = 0.5 * M_PI + 0.5 * 0.0000000000000000;
    double w0 = 0.5 * M_PI + 0.5 * -43000.0008100000050000;
    double w1 = 0.5 * M_PI + 0.5 * -0.0000000000000000;
    double m = 1;
    double g = 43000000.81;
    double l1 = 0.1;
    double l2 = 0.1;
    double h = 0.0001;

    a = ((-g / l1) * sin(q0) + (sin(q1 - q0) * (cos(q1 - q0) * (w0 * w0 + (g / l1) * cos(q0)) + l2 * (w1 * w1 / l1))) / (m + pow(sin(q1 - q0), 2)));
    a = h * a;

    b = h * ((-g / l1) * sin(q0) + (sin(q1 - q0) * (cos(q1 - q0) * (w0 * w0 + (g / l1) * cos(q0)) + l2 * (w1 * w1 / l1))) / (m + pow(sin(q1 - q0), 2)));

    printf("%.20lf ", a);
    printf("%.20lf", b);
    return 0;
}

I do the same calculations with a and b, just with the difference that I get the value of a in two steps, and the b in one.

My code returns: -629.47620126173774000000 -629.47620126173763000000

What is the reason of the difference between the two last decimals?

Actually you can not compare double numbers like strings. There are always floating point errors. — Krassi Em, Oct 11 '17 at 00:50
I got different result, http://rextester.com/ZSM95350. Well maybe your system use strange float representation. This is strange because your operations are in the same order so the result should be the same. What ever never test floating number equality. — Stargateur, Oct 11 '17 at 01:00

score 2 · Accepted Answer · answered Oct 11 '17 at 01:13

The C Standard (99 and 11) says:

The values of operations with floating operands and values subject to the usual arithmetic conversions and of floating constants are evaluated to a format whose range and precision may be greater than required by the type.

So in an expression such as h*(X+Y) as you have in the assignment to b, the implementation is allowed to use greater precision for the intermediate result of X+Y than can possibly be stored in a double, even though the type of the subexpression is still considered to be double. But in a=X+Y; a=h*a;, the first assignment forces the value to be one that actually can be stored in a double, causing a slightly different result.

Another possibility is that the compiler has done "floating point contraction". To quote the C Standard again,

A floating expression may be contracted, that is, evaluated as though it were an atomic operation, thereby omitting rounding errors implied by the source code and the expression evaluation method.

This would most likely happen if the processor has a single instruction that can do a floating point addition and then a multiplication in one step, and the compiler decided to use it.

Assuming one or both of these is the cause, your value b is likely a more accurate representation of the computation you specified (given that all the inputs were restricted to values which can be represented in a double).

The cppreference page about macro FLT_EVAL_METHOD discusses both of these issues in a little more detail. It may be interesting to find out your value of FLT_EVAL_METHOD and play with #pragma STDC FP_CONTRACT OFF.

`#include #include int main(void) { printf("%d\n", (int)FLT_EVAL_METHOD); }` — aschepler, Oct 11 '17 at 04:35
I don't know why, but I get an error: 'FLT_EVAL_METHOD' undeclared. I also tested the code using long double precision, I got -629.47620126173786 for both outputs — Rafael, Oct 11 '17 at 13:41
@Rafael If you're using gcc or clang, make sure to try it with the `-std=c99` or `-std=c11` flag. — aschepler, Oct 11 '17 at 23:24

score 0 · Answer 2 · answered Oct 11 '17 at 00:55

0

The answer is because in the float numbers calculations the equation a=bcde does not have to be equal to x = bc y = de and a.= xy.

It is because floating point arithmetic has a limited precision.

answered Oct 11 '17 at 00:55

0___________

60,014
4
34
74

And what way has better accuracy? Thanks for answering – Rafael Oct 11 '17 at 01:07

Why does this code return two different values doing the same?

2 Answers2