Search code examples
c++assemblyc++builderinline-assemblyx87

Power function in inline 64-bits Assembly on C++ Builder for floating points base


I am using Embarcadero C++ Builder 12 and I am having a giant problem with my power function in Assembly, its called vpow(base, exp). It works perfectly with integers, but I need it to work with doubles. My target is a common PC processor x86-64. Here is my function:

double vpow(double base, int exp) {
   double result = 0.0;
   __asm {
      finit
      fld base
      fild exp
      fyl2x
      fld st(0)
      frndint
      fsub st(1), st(0)
      fxch st(1)
      f2xm1
      fld1
      fadd
      fscale
      fstp result
   }
   return result;
}

What am I doing wrong? I need it working for floating point bases. Like 2,5^3= 15,625. The exp, for now, will be integer. It's 64 bits assembly and it needs to run inline in the C++ Builder IDE.

I tried

double vpow(double base, int exp) {
   double result = 0.0;
   __asm {
      finit
      fld base
      fild exp
      fyl2x
      fld st(0)
      frndint
      fsub st(1), st(0)
      fxch st(1)
      f2xm1
      fld1
      fadd
      fscale
      fstp result
   }
   return result;
}

I just need it to work for doubles too. It works for integers already!


Solution

  • The fix is simple. Read the specification of FYL2X or use a debugger. It is a bit odd to be using x87 instructions today on x64 - they are much slower than AVX2 (although still default code generation on some compilers).

    FYL2X computes ST(1)*log2(ST(0)) you have loaded the arguments onto the stack in the wrong order! FINIT also removed since it isn't at all helpful.

    If it appeared to work OK for integer arguments it can only be because of inadequate testing.

    #include <stdio.h>
    
    double vpow(double base, int exp) {
    double result = 0.0;
    __asm {
           fild exp   // note order of arguments swapped here
           fld base
       
           fyl2x
           fld st(0)
           frndint
           fsub st(1), st(0)
           fxch st(1)
           f2xm1
           fld1
           fadd
           fscale
           fstp result
        }
        return result;
    }
    
    
    int main()
    {
       printf("%f\n", vpow(2.5, 2));
       for (int i=0; i<6; i++)
       printf("%18.9g\n", vpow(1.1, i));
    }
    

    Is a sufficient test case to validate that it works for doubles. Output is

    6.250000
                 1
               1.1
              1.21
             1.331
            1.4641
           1.61051