the pow() function
Author Message the pow() function

The pow() function.

When I use it with floats I get back a number that is very close to the
number I get on a calculator, but not quite.  For example:
---------------------------
float num;

num=pow(17.6, 4.5)
---------------------------

num is:  402538.16

the same calculation on a
calculator yields:    402538.11

the difference may seem small but it escalates into a bigger error as I
carry out more calculations on the result.

Is this normal?  Is there a way of correcting this?

Thanx

Gareth

Mon, 14 Oct 2002 03:00:00 GMT  the pow() function

Quote:
> The pow() function.

> When I use it with floats I get back a number that is very close to
the
> number I get on a calculator, but not quite.  For example:
> ---------------------------
> float num;

> num=pow(17.6, 4.5)

Hi Gareth,

First off, please note that the "pow()" function is defined as having
type "double" arguments and a type "double" return value. If you use
it with "float" variables you force your compiler to do automatic type
conversions. This applies only to your "num" variable. Floating point
constants, like 17.6, have type double by default.

Quote:
> num is:  402538.16

This might depend a bit on what precision you use for printing the
result.

Quote:
> the same calculation on a
> calculator yields:    402538.11

On my calculator it yields:  402538.1119

Quote:
> the difference may seem small but it escalates into a bigger error as
I
> carry out more calculations on the result.

> Is this normal?  Is there a way of correcting this?

Your main problem is very likely that some floating point values can not
be represented exactly in binary notation. Have a look at the answer to
this question in the comp.lang.c FAQ list:

14.1:   When I set a float variable to, say, 3.1, why is printf printing
it as 3.0999999?

--
Stephan
bringing the c.l.c campaign against grumpiness into the next millenium
fight the bug, read the FAQ: http://www.eskimo.com/~scs/C-faq/top.html

Sent via Deja.com http://www.deja.com/

Mon, 14 Oct 2002 03:00:00 GMT  the pow() function

Quote:

> float num;
> num=pow(17.6, 4.5)

> num is:  402538.16
> the same calculation on a
> calculator yields:    402538.11
> Is this normal?

On my PC (the repr. of floats are ieee 754 on PC, I think), this gives:
402538.111856 with double
402538.125000 with float
Note that pow takes and returns a double...

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

#ifdef TEST
#define REAL float
#else
#define REAL double
#endif

int
main (void) {
REAL num= pow (17.6, 4.5);
printf ("%f\n", num);
return 0;

Quote:
}

Mon, 14 Oct 2002 03:00:00 GMT  the pow() function

Quote:

>The pow() function.

>When I use it with floats I get back a number that is very close to the
>number I get on a calculator, but not quite.  For example:
>---------------------------
>float num;

>num=pow(17.6, 4.5)
>---------------------------

>num is:  402538.16

>the same calculation on a
>calculator yields:    402538.11

>the difference may seem small but it escalates into a bigger error as I
>carry out more calculations on the result.

To get higher accuracy, you can try type double instead of type float.

Quote:
>Is this normal?  Is there a way of correcting this?

Yes. No.

There is no way, to get an exact representation of the result
of a call to (say) pow for general arguments. The result has an
infinite number of "digits", so at some point there must be some
rounding or truncation. If you ask anybody, to write down the
exact result of sqrt(2) in decimal, he won't be able to do it.

BTW. The relative error in the C example is about 1e-7. This is
what you usually can expect from type float. With type double, you
can often expect a relative error of about 1e-16.

Mon, 14 Oct 2002 03:00:00 GMT  the pow() function

> The pow() function.
...
> num=pow(17.6, 4.5)
...
> num is:  402538.16
>
> the same calculation on a
> calculator yields:    402538.11

It looks like a poor implementation of pow.  The value is approximately:
402538.111855733998865661844358139
even when rounded to float the error is too large.  However, it might also
be a poor implementation of printf.

Note that the pow function is notably difficult to get right!
--
dik t. winter, cwi, kruislaan 413, 1098 sj  amsterdam, nederland, +31205924131
home: bovenover 215, 1025 jn  amsterdam, nederland; http://www.cwi.nl/~dik/

Mon, 14 Oct 2002 03:00:00 GMT  the pow() function

Quote:

> The pow() function.

> When I use it with floats I get back a number that is very close to the
> number I get on a calculator, but not quite.  For example:
> ---------------------------
> float num;

> num=pow(17.6, 4.5)
> ---------------------------

> num is:  402538.16

> the same calculation on a
> calculator yields:    402538.11

Since a float is guaranteed only 6 significant digits and you are
quibbling about the 8th digit, this is a non-issue.  Since there is
almost never a reson to prefer floats to doubles, change your code
appropiately.

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

int main(void)
{
float x;
double y;
x = y = pow(17.6, 4.5);
printf("The OP expected 402538.11 and obtaind 402538.16\n"
"In the following, there are 2 extra digit beyond the\n"
"Implementation's claimed significance,\n"
"just as in the OP's article.\n\n");
printf("The float result (%.*e) is between\n  %.*e and %.*e\n",
FLT_DIG + 1, x,
FLT_DIG + 1, x / (1 + FLT_EPSILON),
FLT_DIG + 1, x * (1 + FLT_EPSILON));
printf("The double result (%.*e) is between\n  %.*e and %.*e\n",
DBL_DIG + 1, y,
DBL_DIG + 1, y / (1 + DBL_EPSILON),
DBL_DIG + 1, y * (1 + DBL_EPSILON));
return 0;

Quote:
}

The OP expected 402538.11 and obtaind 402538.16
In the following, there are 2 extra digit beyond the
Implementation's claimed significance,
just as in the OP's article.

The float result (4.0253812e+05) is between
4.0253808e+05 and 4.0253817e+05
The double result (4.0253811185573414e+05) is between
4.0253811185573402e+05 and 4.0253811185573426e+05

Quote:

> the difference may seem small but it escalates into a bigger error as I
> carry out more calculations on the result.

> Is this normal?  Is there a way of correcting this?

Use a type with more than 6 digits of significance, obviously.

--

What one knows is, in youth, of little moment; they know enough who
know how to learn. - Henry Adams

A thick skin is a gift from God. - Konrad Adenauer
__________________________________________________________
Fight spam now!
Get your free anti-spam service: http://www.brightmail.com

Mon, 14 Oct 2002 03:00:00 GMT  the pow() function

Quote:

>  > The pow() function.
> ...
>  > num=pow(17.6, 4.5)
> ...
>  > num is:  402538.16

>  > the same calculation on a
>  > calculator yields:    402538.11
> It looks like a poor implementation of pow.  The value is approximately:
>      402538.111855733998865661844358139
> even when rounded to float the error is too large.

naw, FLT_DIG is only guaranteed to be 6, and it's usually not much better
(if any).  assuming your platform has 32-bit IEEE floats, just try:

volatile float a = 402538.111855733998865661844358139;
printf("%f\n", a);

on my system, i get 402538.125000.  that's just the way floats are.

Quote:
> However, it might also
> be a poor implementation of printf.

could be.  no way to know for sure i guess.

--
/"\                              m i k e    b u r r e l l

X        AGAINST HTML MAIL
/ \

Mon, 14 Oct 2002 03:00:00 GMT  the pow() function

Quote:

> > The pow() function.
>...
> > num=pow(17.6, 4.5)
>...
> > num is:  402538.16
>It looks like a poor implementation of pow.  The value is approximately:
>     402538.111855733998865661844358139

When I try with my multi-precision library, I get two different last
digits. Is this what you mean by approximately, or should I start to
worry?

Quote:
>even when rounded to float the error is too large.  However, it might also
>be a poor implementation of printf.

Is it correct, that your comment assumes IEEE-754 floats?
I have already given an answer, but after your post I got curious.
If I haven't done anything stupid, the three closest float values
to the correct result (assuming double parameters) are (exactly)
in float precision (23 bit + 1 implicit bit mantissa):

402538.09375, 402538.125, 402538.15625

And this would of course mean a very poor (double) pow implementation.
Even the naive exp(log(17.6)*4.5) would get this correct for float
precision.

But when I assume arguments in float precision (4.5 is exact,
but 17.6 will be stored as roughly 17.60000038), the result
indicated by the OP seems correct.

So I could speculate, that the original code was something like

float num, x=17.6, y=4.5;
num = pow(x,y);

and 402538.16 would be the correctly rounded and expected result.

--

Mon, 14 Oct 2002 03:00:00 GMT  the pow() function
On Thu, 27 Apr 2000 15:20:13 +0100, that hoopy frood "Gareth"

Quote:

>The pow() function.

>When I use it with floats I get back a number that is very close to the
>number I get on a calculator, but not quite.  For example:
>---------------------------
>float num;

>num=pow(17.6, 4.5)
>---------------------------

>num is:  402538.16

>the same calculation on a
>calculator yields:    402538.11

And Windows 98's Calculator yields: 402538.111855733998865661844358115

Quote:
>the difference may seem small but it escalates into a bigger error as I
>carry out more calculations on the result.

>Is this normal?  Is there a way of correcting this?

It's to be expected.  The problem is that computers can't store real
numbers to infinite precision.  Eventually, you have to round off and
be content with your fairly-accurate number.

You could try changing types from float to double to gain a little bit
of breathing space, but the floating point inaccuracy problem is
inescapable.

If your program *must* have perfect accuracy, you need to use some
sort of integral type.  How you choose to represent real numers as
integers is up to you.  :)

--
Chris Mears

ICQ: 36697123

C-FAQ: http://www.eskimo.com/~scs/C-faq/top.html

Tue, 15 Oct 2002 03:00:00 GMT  the pow() function

> >  > num is:  402538.16
...
> > It looks like a poor implementation of pow.  The value is approximately:
> >      402538.111855733998865661844358139
> > even when rounded to float the error is too large.
>
> naw, FLT_DIG is only guaranteed to be 6, and it's usually not much better
> (if any).  assuming your platform has 32-bit IEEE floats, just try:
>         volatile float a = 402538.111855733998865661844358139;
>         printf("%f\n", a);
> on my system, i get 402538.125000.  that's just the way floats are.

Yup, in HEX you get 6246A.2.  The proper value is between 6246A.1 and
6246A.2 (402538.0625 and 402538.125), both representable numbers as IEEE
floats.  The proper value starts with 6246A.1C.  Obviously both are
acceptable, but the second more so than the first.  However the OP got
402538.16, which would be 6246A.28FC..., the closest representable float
is 6246A.3, which is 402538.1875.  On this inpection I suspect the
implementation of printf of inaccuracies.

> naw, FLT_DIG is only guaranteed to be 6, and it's usually not much better
> (if any).
is wrong.  FLT_DIG is the number of decimal digits such that a decimal
f-p number with that number of digits can be rounded to float and back
again to decimal without change.  The actual precision of floating point
must be better to get this guarantee, and indeed is, and in some ranges
even much better.  Another question you might ask is how many decimal
digits do I need when I print out a number such that I read it back in
I get the same number.  This latter number is larger than FLT_DIG, for
IEEE float it is 8.  So when you want different floating point numbers
to print out differently you need to print out 8 digits, and a
well-behaved version of printf would indeed do this correctly.
--
dik t. winter, cwi, kruislaan 413, 1098 sj  amsterdam, nederland, +31205924131
home: bovenover 215, 1025 jn  amsterdam, nederland; http://www.cwi.nl/~dik/

Tue, 15 Oct 2002 03:00:00 GMT  the pow() function
I must be slipping.  I wrote the following, I note corrections:

> Yup, in HEX you get 6246A.2.  The proper value is between 6246A.1 and
> 6246A.2 (402538.0625 and 402538.125), both representable numbers as IEEE
> floats.

Must be 6246A.18 and 6246A.20 (402538.09375 and 402538.125).

>                                                     However the OP got
> 402538.16, which would be 6246A.28FC..., the closest representable float
> is 6246A.3, which is 402538.1875.

6246A.28 which is 402538.15625.
--
dik t. winter, cwi, kruislaan 413, 1098 sj  amsterdam, nederland, +31205924131
home: bovenover 215, 1025 jn  amsterdam, nederland; http://www.cwi.nl/~dik/

Tue, 15 Oct 2002 03:00:00 GMT  the pow() function

Quote:

> However, your comment:
>  > naw, FLT_DIG is only guaranteed to be 6, and it's usually not much better
>  > (if any).
> is wrong.  FLT_DIG is the number of decimal digits such that a decimal
> f-p number with that number of digits can be rounded to float and back
> again to decimal without change.  The actual precision of floating point
> must be better to get this guarantee, and indeed is, and in some ranges
> even much better.

well i don't know about 'wrong'.  people (when dealing with floats) usually
speak in base 10, so it's not a bad metric: all the nice degits in
314.31129381 don't help much (besides looking cool, i guess) when the number
is supposed to be 314.32 instead of 314.31.

--
/"\                              m i k e    b u r r e l l

X        AGAINST HTML MAIL
/ \

Tue, 15 Oct 2002 03:00:00 GMT  the pow() function

(snip)

Quote:
>And this would of course mean a very poor (double) pow implementation.
>Even the naive exp(log(17.6)*4.5) would get this correct for float
>precision.

I presume that is what pow() does.  There isn't much else it
can do.  For integral exponents it could do it differently,
but for non-intergral ones there isn't much else.  Note that if the
base and exponent only have float (24 bit) precision the pow() will
have somewhat less than 24 bit precision.

Quote:
>But when I assume arguments in float precision (4.5 is exact,
>but 17.6 will be stored as roughly 17.60000038), the result
>indicated by the OP seems correct.
>So I could speculate, that the original code was something like
>  float num, x=17.6, y=4.5;
>  num = pow(x,y);
>and 402538.16 would be the correctly rounded and expected result.

Yep.

-- glen

Sat, 19 Oct 2002 03:00:00 GMT

 Page 1 of 2 [ 20 post ] Go to page:  

Relevant Pages
 12. pow function