4

So I was googling for a long time and i found almost nothing. I found some info about possible implementation of Math.Pow from this url, but they are inaccurate, for example this code

public static double PowerA(double a, double b)
{
    int tmp = (int)(BitConverter.DoubleToInt64Bits(a) >> 32);
    int tmp2 = (int)(b * (tmp - 1072632447) + 1072632447);
    return BitConverter.Int64BitsToDouble(((long)tmp2) << 32);
}
static void Main(string[] args)
{
    double x = 12.53, y = 16.45;
    Console.WriteLine(Math.Pow(x, y));
    Console.WriteLine(PowerA(x, y));
}

provides output:

1,15158266266297E+18
8,9966384455562E+17

So inaccurate...

I was thinking that it works like a sum of series but I don't know for certain.

Community
  • 1
  • 1
Alex Zhukovskiy
  • 9,565
  • 11
  • 75
  • 151
  • 3
    The obligatory [What Every Computer Scientist Should Know About Floating-Point Arithmetic](http://docs.oracle.com/cd/E19957-01/806-3568/ncg_goldberg.html). I'm sure there are bits you didn't know about floating point in that paper. – ppeterka Sep 25 '13 at 18:08
  • You can download code that should be equivalent to the actual code, if you want to. [SSCLI](http://www.microsoft.com/en-us/download/details.aspx?id=4917). It's not been updated for a while, but `Math.Pow` is fairly old. – Damien_The_Unbeliever Sep 25 '13 at 18:17
  • 11
    See [How is Math.Pow() implemented in .Net Framework?](http://stackoverflow.com/q/8870442/781792) – Tim S. Sep 25 '13 at 18:18
  • 1
    I would suggest that PowerA is the one that is inaccurate. PowerA(12.0, 5.0) returns the incorrect 227008.5 whereas Math.Pow(12.0, 5.0) returns the correct 248832, which is small enough to test. Plus a whole number double raised to a whole number double should return another whole number double, and not something ending with .5! – Rick Davin Sep 25 '13 at 18:41
  • 2
    @RickDavin I read that as what the OP was saying. The method he has doesn't correctly mirror `Math.Pow`, it is incorrect, and he's looking for an implementation that will correctly mirror `Math.Pow`. – Servy Sep 25 '13 at 19:38
  • to: Tim S - i've readed it before asking, but there's only explanation that it's native C++ code so this is why i published this question. to: ppeterka 66 - i've readed an article based on article from your URL. So the greater part of it i know, but i'l defenitly read the original. But how can it help me? Link to disassembler manual should be more practice :) – Alex Zhukovskiy Sep 25 '13 at 22:18
  • The linked page does state "In my tests it usually within an error margin of 5% to 12%, in extreme cases sometimes up to 25%." @RickDavin 's results are well within defined and documented error margins. – sisve Sep 29 '13 at 04:51
  • those are approximations not computations of pow which means that it is not accurate on whole RxR space just in some small portion of it. To compute pow you need to transform a from a^b to known base usually 2 so you compute 2^c, for more info see my answer – Spektre Sep 29 '13 at 08:16

1 Answers1

7

pow is usually evaluated by this formula:

x^y = exp2(y*log2(x))

Functions exp2(x),log2(x) are directly implemented in FPU. If you want to implement bignums then they can also be evaluated by basic operators with use of precomputed table of sqrt-powers like:

2^1/2, 2^1/4, 2^1/8, 2^1/16, 2^1/32 ...

to speed up the process

In case you need to handle also rooting for negative bases see this:

Spektre
  • 49,595
  • 11
  • 110
  • 380
  • functions exp2(x),log2(x) are directly implemented in FPU - it's the most important info that i didn't know. Tnx. – Alex Zhukovskiy Sep 29 '13 at 15:13
  • functions exp2(x),log2(x) are directly implemented in FPU - are you sure you don't mean *FP library*, since FPU is a hardware entity that varies enormously with CPU architecture (not explicit in question) and may NOT NECESSARILY implement high-order functions like this. – Eight-Bit Guru Mar 07 '15 at 12:16
  • @Eight-BitGuru I work on `x87` FPU's and there these functions are present natively , so yes I mean implemented in FPU not in FP lib. – Spektre Mar 07 '15 at 12:40
  • @Eight-BitGuru see instructions like `f2xm1,fyl2x,fyl2xp1` ... by "work on" I mean "I work with" of coarse ... the wrong translation just occur to me :) – Spektre Mar 07 '15 at 12:48
  • @Spektre I think you missed the core of my point - not all FPU architectures implement high-order functions of this type, so asserting that 'functions exp2(x),log2(x) are directly implemented in FPU' is making an assumption about the nature of the hardware which is not explicitly specified in the question. – Eight-Bit Guru Mar 07 '15 at 16:21
  • @Eight-BitGuru OP tags and text implies PC platform which contains x87 and the question is about pow not exp2,log2 so that is why I do not see your point. if exp2 and log2 are not presents then they can be computed by elementary functions. for example log can be done like this: first exponent is computed (on integers) then FP is converted to fixed point and binary search is applied on mantisa. this is for example sqrt implementation http://stackoverflow.com/a/28808095/2521214 in this way. – Spektre Mar 08 '15 at 09:02
  • @Eight-BitGuru exp2 can be computed directly by table of `2^(1/(2^i))` powers where i is integer>1 value by exploiting 2^(a0+a1+a2+...)=(2^a0)*(2^a1)*... btw mine implementation of log2 use this same table for speed up,... and there are many other ways for sure – Spektre Mar 08 '15 at 09:21
  • @Spektre: I suppose `log2x` is base 2 logarithm but what is `exp2x` that is implemented in FPU? – Jim May 09 '22 at 10:53
  • @Jim `exp2` is `pow(2,x)` its also implemented in FPU directly ... you know log,exp are inverse to each other ... `log2(exp2(x)) = x` – Spektre May 09 '22 at 14:17