representational error in C

Discussion in 'C Programming' started by newgoat, Mar 11, 2006.

  1. newgoat

    newgoat Guest

    Hello,

    I wrote two program to test the inherent problem binary system has
    with representing floating points values. Adding 0.1 to itself 1000
    times results in a number slightly larger than 100. But to my surprise,
    adding it to itself 100000 times results in a number smaller than
    10000.
    Can anybody tell me why? Thanks!

    #include <stdio.h>

    float add_dimes_fversion(int);
    double add_dimes_dversion(int);

    int main(void){
    int add_times1 = 1000;
    int add_times2 = 100000;

    /* Add 0.1 to itself for 1000 times using float for all data would
    result in 100.099045. */
    printf("\nAdding 0.1 for 1000 times using float data results in
    %.10f\n", add_dimes_fversion(add_times1));



    /* Add 0.1 to itself for 100000 times using float for all data would
    result in 9998.6562500000. */
    printf("\nAdding 0.1 for 100000 times using float data results in
    %.10f\n\n", add_dimes_fversion(add_times2));



    /* Add 0.1 to itself for 1000 times using float for all data would
    result in 100.09999999999858744104. */
    printf("\nAdding 0.1 for 1000 times using double data results in
    %.20lf\n",add_dimes_dversion(add_times1));


    /* Add 0.1 to itself for 100000 times using float for all data would
    result in 9998.65625000000000000000. */
    printf("\nAdding 0.1 for 100000 times using double data results in
    %.20lf\n\n", add_dimes_fversion(add_times2));

    printf("0.1 * 1000 = %lf, 0.1 * 100000 = %lf\n\n", 0.1 * 1000.0, 0.1
    * 100000);
    return 0;
    }

    /* Test with float data type */
    float add_dimes_fversion(int n){
    float sum = 0.0;
    int counter = 0;

    while(counter <= n){
    sum = sum + 0.1;
    counter++;
    }

    return sum;
    }

    /* Test with double data type */
    double add_dimes_dversion(int n){
    double sum = 0.0;
    int counter = 0;

    while(counter <= n){
    sum = sum + 0.1;
    counter++;
    }

    return sum;
    }
    newgoat, Mar 11, 2006
    #1
    1. Advertising

  2. In article <>,
    newgoat <> wrote:
    >I wrote two program to test the inherent problem binary system has
    >with representing floating points values. Adding 0.1 to itself 1000
    >times results in a number slightly larger than 100. But to my surprise,
    >adding it to itself 100000 times results in a number smaller than
    >10000.
    >Can anybody tell me why? Thanks!


    The code you included does not appear to contain anything
    specific to C99. If you are using C89 (likely because true C99
    is still fairly uncommon), then some of the details of floating
    point implementation are not nailed down, and are thus implementation
    dependant.

    The exact answers produced by the code you provided would
    depend upon the rounding behaviour in effect. There are
    a few common -different- rounding behaviours (and probably some
    obscure ones along the way). C89 does not provide any mechanism
    to query or control rounding behaviours (but your OS might.)

    One of the recommended rounding behaviours alternates
    between rounding up and rounding down, in a deliberate attempt
    to reduce accumulation of round-off errors.

    My guess, though, is that your system is probably not using
    decimal arithmetic, and is instead using a binary floating
    point representation. The representation of 100000 probably
    does not start with the same binary expansion as 100 does.
    This could result in the least-significant bit being
    "crowded off the end" as you accumulated values: if that
    least significant bit was the difference between underflow of
    100000 and overflow of that value, then Yes, the effect
    you saw could definitely happen.
    --
    "It is important to remember that when it comes to law, computers
    never make copies, only human beings make copies. Computers are given
    commands, not permission. Only people can be given permission."
    -- Brad Templeton
    Walter Roberson, Mar 11, 2006
    #2
    1. Advertising

  3. newgoat

    Me Guest

    newgoat wrote:
    > Hello,
    >
    > I wrote two program to test the inherent problem binary system has
    > with representing floating points values. Adding 0.1 to itself 1000
    > times results in a number slightly larger than 100. But to my surprise,
    > adding it to itself 100000 times results in a number smaller than
    > 10000.
    > Can anybody tell me why? Thanks!


    Floats work like this (using whole numbers for niceness):

    (use a fixed width font for this)
    ***** * * * * * * * * * ... *
    01234 6 8 10121416 20 24 28 64

    What I tried to represent here is that there are "goal posts" at
    2**(2n) [note: this is not what C does]

    0 == 0
    4 == 2**2
    16 = 2**4
    64 = 2**6
    256 = 2**8
    ....

    where the distance between 0-4 is 1 unit, 4-16 is 2 units, 16-64 is 4
    units, 64-256 is 8 units, etc. With the C standard, these "goal posts"
    are actually powers of FLT_RADIX.

    (I'm going to assume your computer uses IEEE-754 floats, so FLT_RADIX
    == 2) The problem with your example is that 0.1 isn't a power of
    FLT_RADIX (i.e. it's not one of the goal posts, it's a number between
    the two goal posts).

    There are three issues here:

    a. 0.1 isn't exactly representable. The standard doesn't say what to do
    in this case but odds are, your compiler will round it to the nearest
    number. In our whole number example, it would be like what do to when
    given number 18.

    b. 0.1 isn't a goal post. Goal posts are good. If we start from 0, we
    can keep adding by goal posts up to some number *exactly*. Our whole
    number example: 0 to 16 by adding 2, 0 to 64 by adding 4.

    c. If a number isn't exactly representable, the C standard doesn't say
    what to do. Most machines will round to nearest. For our whole number
    machine, lets use round to rand to approximate round to nearest because
    our whole number machine has a very small range compared to real floats
    to save typing.


    So with all this in mind lets replace 0.1 with 3.5 and 1000 with 4

    So 3.5 can either be 3 or 4. Lets pick 3.

    0+3 = 7, lets round to 6
    6+3 = 9, lets round to 10
    10+3 = 13, lets round to 14
    14+3 = 17, lets round to 16

    3.5 * 4 == 14 but because of all the rounding that took place, we got
    16.

    If your computer was a (non-existant?) one that has FLT_RADIX == 10,
    this would have worked. Since IEEE-754/FLT_RADIX==2 floats are pretty
    much a defacto standard for C (and have desirable numerical
    properties), there is a proposal to add a separate decimal float type
    to C which would have ran your example with no surprises.

    > #include <stdio.h>
    >
    > float add_dimes_fversion(int);
    > double add_dimes_dversion(int);
    >
    > int main(void){
    > int add_times1 = 1000;
    > int add_times2 = 100000;
    >
    > /* Add 0.1 to itself for 1000 times using float for all data would
    > result in 100.099045. */
    > printf("\nAdding 0.1 for 1000 times using float data results in
    > %.10f\n", add_dimes_fversion(add_times1));


    p.s. .9 is enough for IEEE float (see
    http://www.open-std.org/jtc1/sc22/wg21/docs/papers/2005/n1822.pdf for
    explanation)

    > /* Add 0.1 to itself for 100000 times using float for all data would
    > result in 9998.6562500000. */
    > printf("\nAdding 0.1 for 100000 times using float data results in
    > %.10f\n\n", add_dimes_fversion(add_times2));
    >
    >
    >
    > /* Add 0.1 to itself for 1000 times using float for all data would
    > result in 100.09999999999858744104. */
    > printf("\nAdding 0.1 for 1000 times using double data results in
    > %.20lf\n",add_dimes_dversion(add_times1));


    p.s. .17 is enough for IEEE double

    > /* Add 0.1 to itself for 100000 times using float for all data would
    > result in 9998.65625000000000000000. */
    > printf("\nAdding 0.1 for 100000 times using double data results in
    > %.20lf\n\n", add_dimes_fversion(add_times2));
    >
    > printf("0.1 * 1000 = %lf, 0.1 * 100000 = %lf\n\n", 0.1 * 1000.0, 0.1
    > * 100000);
    > return 0;
    > }
    >
    > /* Test with float data type */
    > float add_dimes_fversion(int n){
    > float sum = 0.0;
    > int counter = 0;
    >
    > while(counter <= n){
    > sum = sum + 0.1;


    Note: 0.1 is a double literal. To get a float literal use 0.1f. What
    you're code does here is:
    sum = (double)sum + 0.1;
    The standard allows a compiler to do this with float literals but using
    a double literal forces it to happen when it might not happen
    otherwise.

    > counter++;
    > }
    >
    > return sum;
    > }
    >
    > /* Test with double data type */
    > double add_dimes_dversion(int n){
    > double sum = 0.0;
    > int counter = 0;
    >
    > while(counter <= n){
    > sum = sum + 0.1;
    > counter++;
    > }
    >
    > return sum;
    > }
    Me, Mar 11, 2006
    #3
  4. newgoat

    Michael Mair Guest

    newgoat schrieb:
    > Hello,
    >
    > I wrote two program to test the inherent problem binary system has
    > with representing floating points values. Adding 0.1 to itself 1000
    > times results in a number slightly larger than 100. But to my surprise,
    > adding it to itself 100000 times results in a number smaller than
    > 10000.
    > Can anybody tell me why? Thanks!

    <snip!>
    > /* Test with float data type */
    > float add_dimes_fversion(int n){
    > float sum = 0.0;
    > int counter = 0;
    >
    > while(counter <= n){
    > sum = sum + 0.1;
    > counter++;
    > }
    >
    > return sum;
    > }


    For one thing, you are working with 1001 and 100001 summands,
    respectively; assuming 0.1F*10 < 1, this is not as surprising
    as it seems: Let us, for a moment, assume that 0.1F was equal
    to 0.09991 in the real world; even if we have no rounding errors,
    0.09991*1001 = 99.91+0.0991 = 100.0091 and
    0.09991*1000001 = 99910+0.0991 = 99910.0991

    Apart from that, the other replies explain most.
    Have a look at
    http://docs.sun.com/source/806-3568/ncg_goldberg.html
    if you are interested in more information.

    Cheers
    Michael
    --
    E-Mail: Mine is an /at/ gmx /dot/ de address.
    Michael Mair, Mar 11, 2006
    #4
  5. In article <>,
    Me <> wrote:
    >(I'm going to assume your computer uses IEEE-754 floats, so FLT_RADIX
    >== 2)


    >Since IEEE-754/FLT_RADIX==2 floats are pretty
    >much a defacto standard for C


    There are quite a few widespread IEEE-754 systems that use FLT_RADIX=16
    or FLT_RADIX=4, so I would say that FLT_RADIX=2 is *not* a defacto
    standard for C.

    --
    "law -- it's a commodity"
    -- Andrew Ryan (The Globe and Mail, 2005/11/26)
    Walter Roberson, Mar 11, 2006
    #5
  6. newgoat

    Chris Torek Guest

    In article <duuts1$ent$>
    Walter Roberson <-cnrc.gc.ca> wrote:
    >There are quite a few widespread IEEE-754 systems that use FLT_RADIX=16
    >or FLT_RADIX=4, so I would say that FLT_RADIX=2 is *not* a defacto
    >standard for C.


    Just out of curiosity: I know IBM S/370 (and S/360) architecture
    uses 16 ("hex float"); but who uses 4? (And: Why? 4 would lose
    the implied-1 advantage of base 2, and add the jitter seen in hex
    float, yet not increase the exponent range significantly as 16
    does.)
    --
    In-Real-Life: Chris Torek, Wind River Systems
    Salt Lake City, UT, USA (40°39.22'N, 111°50.29'W) +1 801 277 2603
    email: forget about it http://web.torek.net/torek/index.html
    Reading email is like searching for food in the garbage, thanks to spammers.
    Chris Torek, Mar 11, 2006
    #6
  7. newgoat

    Joe Wright Guest

    Walter Roberson wrote:
    > In article <>,
    > Me <> wrote:
    >> (I'm going to assume your computer uses IEEE-754 floats, so FLT_RADIX
    >> == 2)

    >
    >> Since IEEE-754/FLT_RADIX==2 floats are pretty
    >> much a defacto standard for C

    >
    > There are quite a few widespread IEEE-754 systems that use FLT_RADIX=16
    > or FLT_RADIX=4, so I would say that FLT_RADIX=2 is *not* a defacto
    > standard for C.
    >

    By 'quite a few' you mean presumably more than two. Can you name them
    for us please? Thanks.

    --
    Joe Wright
    "Everything should be made as simple as possible, but not simpler."
    --- Albert Einstein ---
    Joe Wright, Mar 11, 2006
    #7
  8. newgoat

    Ben Pfaff Guest

    -cnrc.gc.ca (Walter Roberson) writes:

    > There are quite a few widespread IEEE-754 systems that use FLT_RADIX=16
    > or FLT_RADIX=4, so I would say that FLT_RADIX=2 is *not* a defacto
    > standard for C.


    I was under the impression that IEEE 754 was a binary floating
    point standard, i.e. FLT_RADIX=2. Why would an IEEE 754 system
    claim FLT_RADIX of 4 or 16?
    --
    "Some programming practices beg for errors;
    this one is like calling an 800 number
    and having errors delivered to your door."
    --Steve McConnell
    Ben Pfaff, Mar 11, 2006
    #8
  9. newgoat

    Jordan Abel Guest

    On 2006-03-11, Walter Roberson <-cnrc.gc.ca> wrote:
    > In article <>,
    > Me <> wrote:
    >>(I'm going to assume your computer uses IEEE-754 floats, so FLT_RADIX
    >>== 2)

    >
    >>Since IEEE-754/FLT_RADIX==2 floats are pretty
    >>much a defacto standard for C

    >
    > There are quite a few widespread IEEE-754 systems that use FLT_RADIX=16
    > or FLT_RADIX=4, so I would say that FLT_RADIX=2 is *not* a defacto
    > standard for C.


    Just because some systems deviate from it doesn't mean it's not a de
    facto standard. [Is this really IEEE 754? I thought "IEEE 754" meant a
    system matching the 'ideal' of the Intel 387, which, i've read, is what
    the IEEE standard is based on.]
    Jordan Abel, Mar 11, 2006
    #9
  10. newgoat

    jacob navia Guest

    Michael Mair a écrit :
    > newgoat schrieb:
    >
    >> Hello,
    >>
    >> I wrote two program to test the inherent problem binary system has
    >> with representing floating points values. Adding 0.1 to itself 1000
    >> times results in a number slightly larger than 100. But to my surprise,
    >> adding it to itself 100000 times results in a number smaller than
    >> 10000.
    >> Can anybody tell me why? Thanks!

    >
    > <snip!>
    > > /* Test with float data type */
    > > float add_dimes_fversion(int n){
    > > float sum = 0.0;
    > > int counter = 0;
    > >
    > > while(counter <= n){
    > > sum = sum + 0.1;
    > > counter++;
    > > }
    > >
    > > return sum;
    > > }

    >
    > For one thing, you are working with 1001 and 100001 summands,
    > respectively; assuming 0.1F*10 < 1, this is not as surprising
    > as it seems: Let us, for a moment, assume that 0.1F was equal
    > to 0.09991 in the real world; even if we have no rounding errors,
    > 0.09991*1001 = 99.91+0.0991 = 100.0091 and
    > 0.09991*1000001 = 99910+0.0991 = 99910.0991
    >
    > Apart from that, the other replies explain most.
    > Have a look at
    > http://docs.sun.com/source/806-3568/ncg_goldberg.html
    > if you are interested in more information.
    >
    > Cheers
    > Michael


    Michael, you have a GREAT EYE!!!!!!!!!!

    How EASY is to oversee the dammed <= in that program. I wondered
    and wondered... WHAT THE HELL IS GOING ON... ?????

    jacob
    jacob navia, Mar 11, 2006
    #10
  11. newgoat

    Jack Klein Guest

    On 11 Mar 2006 19:42:52 GMT, Jordan Abel <> wrote
    in comp.lang.c:

    > On 2006-03-11, Walter Roberson <-cnrc.gc.ca> wrote:
    > > In article <>,
    > > Me <> wrote:
    > >>(I'm going to assume your computer uses IEEE-754 floats, so FLT_RADIX
    > >>== 2)

    > >
    > >>Since IEEE-754/FLT_RADIX==2 floats are pretty
    > >>much a defacto standard for C

    > >
    > > There are quite a few widespread IEEE-754 systems that use FLT_RADIX=16
    > > or FLT_RADIX=4, so I would say that FLT_RADIX=2 is *not* a defacto
    > > standard for C.

    >
    > Just because some systems deviate from it doesn't mean it's not a de
    > facto standard. [Is this really IEEE 754? I thought "IEEE 754" meant a
    > system matching the 'ideal' of the Intel 387, which, i've read, is what
    > the IEEE standard is based on.]


    Sigh.

    The original IEEE 754 and 854 formats were based in large part on the
    Intel 8087, which existed half a decade or so before the 387 did.

    But since Intel released the part before standardization was complete,
    the 8087 does not completely conform to the final standard in a few
    areas.

    --
    Jack Klein
    Home: http://JK-Technology.Com
    FAQs for
    comp.lang.c http://c-faq.com/
    comp.lang.c++ http://www.parashift.com/c -faq-lite/
    alt.comp.lang.learn.c-c++
    http://www.contrib.andrew.cmu.edu/~ajo/docs/FAQ-acllc.html
    Jack Klein, Mar 12, 2006
    #11
    1. Advertising

Want to reply to this thread or ask your own question?

It takes just 2 minutes to sign up (and it's free!). Just click the sign up button to choose a username and then you can ask your own questions on the forum.
Similar Threads
  1. hfk0
    Replies:
    2
    Views:
    21,660
  2. JavaQueries
    Replies:
    1
    Views:
    3,658
    John C. Bollinger
    Mar 1, 2005
  3. Balaji
    Replies:
    3
    Views:
    10,093
  4. Bishop
    Replies:
    1
    Views:
    774
    Bishop
    Feb 24, 2007
  5. juvi
    Replies:
    3
    Views:
    1,049
    Alexey Smirnov
    Jan 22, 2009
Loading...

Share This Page