Hello Kai and Everybody!

A couple of comments about whether to use special cases for
pow (x, y) for things like y = 0.5 or y = integer...

On Tue, May 3, 2011 at 6:58 AM, Kai Tietz wrote:
> 2011/5/3 Peter Rockett <[email protected]>:
>> On 03/05/2011 10:11, Kai Tietz wrote:
>> ...
>> The reason is not introducing a discontinuity into pow() around 0.5.
>> Forcing the two functions to agree may seem like a good idea but in
>> reality, you are creating another much more insidious problem. You are
>> embedding a potentially disastrous property into pow().
>>
>>> It might be that I haven't seen the forbidden
>>> thing on the spec.  The sqrt special-case has at least one advantage
>>> over the log2/exp2 variant. First it is faster, as just one FPU
>>> command is used to calculate sqrt. Secondly, is provides better result
>>> and compatible one to gcc's internal used soft-floating-point library.
>>> As I don't see in spec that this special casing shouldn't be done, I
>>> think we will keep it on mingw-w64's side.
>> ...
> I am aware of this what you are telling me, but there is in pow
> function already special-casing for x^y for y with integer values in
> range for y >= INT_MIN and y <= INT_MAX. As here more precise
> calculation via powi is used. This makes that results of pow having
> already discontinuity.
> So pow functions doesn't have the requirement of having continuity.
> ...

I'm not aware that either the C standard or the IEEE-754 standard specifies
one behavior or the other.  (I think that they don't, but I could be wrong.)

I believe that is has been commonplace over the years to use these special
cases for pow.  (This is just the sense I have from personal recollection and
oral tradition.)  I think speed optimization has been the reason, rather than
improved accuracy or agreement with things like sqrt, but perhaps it's been
a mixture of those.  Of course, just because there is precedent, doesn't
mean it's the right way to go.

If I had to choose between pow (x, 0.5) agreeing with sqrt (x) and pow (x, y)
being monotone in y, I think I would choose the latter (i.e., choose not to
use the special cases).  But there are good arguments on both sides.

One solution that would be attractive to me would be to implement pow
"exactly" (i.e., bankers' rounding to the nearest representable value), and
then use all the special cases you want as a speed optimization.  Now
the special cases won't change the value of pow (assuming sqrt is properly
implemented), so the special cases won't introduce any non-monotonicity
or "discontinuities" into pow.

What I don't know is how hard or expensive this is.  I'm pretty sure it's
not too hard to calculate exp "exactly."  But I think that calculating log
is harder.  Perhaps James Beard knows whether it's practical to calculate
pow "exactly," and maybe even how to do it.

(Another approach, which I very much don't recommend, is to use pow
to calculate sqrt.  You get reduced speed and reduced accuracy for
sqrt, but at least you do get sqrt and pow to agree without introducing
non-monotonicity into pow.)

Lastly, I have to believe that this is a well-understood and "solved" problem
in the numerical-analysis community (of which I am not a part).  Perhaps
it would make sense to find out how the experts who have doing this for
fifty-plus years deal with this issue.

> Regards,
> Kai

And thanks to the mingw and mingw64 teams who have invested the
time and effort to actually write the code we're so breezily discussing.
My comments are in no way meant to be a criticism of those efforts.
I'm just offering up one man's opinion from the user community in the
hope that it might be helpful.

Bets regards.


K. Frank

------------------------------------------------------------------------------
WhatsUp Gold - Download Free Network Management Software
The most intuitive, comprehensive, and cost-effective network 
management toolset available today.  Delivers lowest initial 
acquisition cost and overall TCO of any competing solution.
http://p.sf.net/sfu/whatsupgold-sd
_______________________________________________
Mingw-w64-public mailing list
[email protected]
https://lists.sourceforge.net/lists/listinfo/mingw-w64-public

Reply via email to