Szymon - 5 months ago 29

C Question

Our professor said that you can't calculate a^{b} if a<0 using

`pow()`

`pow()`

I have searched through

`math.h`

So the question is how is

`pow()`

Answer

If you're curious how the `pow`

function might be implemented in practice, you can look at the source code. There is a kind of "knack" to searching through unfamiliar (and large) codebases to find the section you are looking for, and it's good to get some practice.

One implementation of the C library is glibc, which has mirrors on GitHub. I didn't find an official mirror, but an unofficial mirror is at https://github.com/lattera/glibc

We first look at the `math/w_pow.c`

file which has a promising name. It contains a function `__pow`

which calls `__ieee754_pow`

, which we can find in `sysdeps/ieee754/dbl-64/e_pow.c`

(remember that not all systems are IEEE-754, so it makes sense that the IEEE-754 math code is in its own directory).

It starts with a few special cases:

```
if (y == 1.0) return x;
if (y == 2.0) return x*x;
if (y == -1.0) return 1.0/x;
if (y == 0) return 1.0;
```

A little farther down you find a branch with a comment

```
/* if x<0 */
```

Which leads us to

```
return (k==1)?__ieee754_pow(-x,y):-__ieee754_pow(-x,y); /* if y even or odd */
```

So you can see, for negative `x`

and integer `y`

, the glibc version of `pow`

will compute `pow(-x,y)`

and then make the result negative if `y`

is odd.

This is not the only way to do things, but my guess is that this is common to many implementations. You can see that `pow`

is full of special cases. This is common in library math functions, which are supposed to work correctly with unfriendly inputs like denormals and infinity.

The `pow`

function is especially hard to read because it is heavily-optimized code which does bit-twiddling on floating-point numbers.

The C standard (n1548 ยง7.12.7.4) has this to say about `pow`

:

A domain error occurs if x is finite and negative and y is finite and not an integer value.

So, according to the C standard, negative `x`

*should* work.

There is also the matter of appendix F, which gives much tighter constraints on how `pow`

works on IEEE-754 / IEC-60559 systems.