Square root metafunction?

Question

Is it possible to compute the square root of an integer with a metafunction with the following signature :

template<unsigned int N> inline double sqrt();

(or maybe using the constexpr keyword, I don't know what is the best). With that, sqrt<2>() would be replaced by 1.414... at compile-time.

What would be the best implementation for a such function ?

Answer 1

This may not be what you are looking for, but I wanted to make sure you realized that typically with optimization the compiler will calculate the result at compile time anyway. For example, if you have this code:

void g()
{
  f(sqrt(42));
}

With g++ 4.6.3 with optimization -O2, the resulting assembly code is:

   9 0000 83EC1C                subl    $28, %esp
  11 0003 DD050000              fldl    .LC0
  12 0009 DD1C24                fstpl   (%esp)
  13 000c E8FCFFFF              call    _Z1fd
  14 0011 83C41C                addl    $28, %esp
  16 0014 C3                    ret
  73                    .LC0:
  74 0000 6412264A              .long   1244009060
  75 0004 47EC1940              .long   1075440711

The sqrt function is never actually called, and the value is just stored as part of the program.

Therefore to create a function that technically meets your requirements, you simply would need:

template<unsigned int N> inline double meta_sqrt() { return sqrt(N); }

Answer 2

Eigen contains this meta_sqrt which uses binary search:

template<int Y,
         int InfX = 0,
         int SupX = ((Y==1) ? 1 : Y/2),
         bool Done = ((SupX-InfX)<=1 ? true : ((SupX*SupX <= Y) && ((SupX+1)*(SupX+1) > Y))) >
                                // use ?: instead of || just to shut up a stupid gcc 4.3 warning
class meta_sqrt
{
    enum {
  MidX = (InfX+SupX)/2,
  TakeInf = MidX*MidX > Y ? 1 : 0,
  NewInf = int(TakeInf) ? InfX : int(MidX),
  NewSup = int(TakeInf) ? int(MidX) : SupX
};
  public:
    enum { ret = meta_sqrt<Y,NewInf,NewSup>::ret };
};

template<int Y, int InfX, int SupX>
class meta_sqrt<Y, InfX, SupX, true>
{
    public:  enum { ret = (SupX*SupX <= Y) ? SupX : InfX };
};

Answer 3

The issue I see is that metaP effectively abuses enums into variables. The problem is that enums are treated internally as integers which precludes trying to get a floating point value out of them. However, you may be able to create your own floating point format that creates two results, an integer portion and an exponent. You'll still have to process this into a float, as Out = Sqrt<42>::mantissa * pow(10,Sqrt<42>::exponent);. Actually determining the values is left as an exercise for the reader, but you will probably have to scale the input upwards (by an even power of 10), computing the root, and storing the -power/2 you used earlier.

To compute sqrt<42>, you would first set the exponent enum to a suitable power, such as '-4' (the lower, the more decimals, but watch for overflow). You then multiply the input by '10^(-2*exponent)'. In this case you get 42*10^8 = 4200000000. You then take the root of this value getting '64807' as the final value. At runtime, you calculate the "val * 10 ^ exponent" = "64807 * 10 ^ -4" = 64807 * 0.0001 = 6.4807m and store it to a float.

The extra conversion work kinda defeats the purpose, but you can reduce it a bit by storing the exponent as 10^k (ie 10^4) then doing out=sqrt<x>::mantissa/sqrt<x>::exponent.

edit I just noticed that with the mantissa/exponent method, the choice of exponent is arbitrary as long as it's bigger than the integer portion of the final root. It can even be a constant, which eases the design of your meta functions. In the case of 42 for example you could select the 'exponent' to always be 6000. You then multiply the input by 6000^2, take the integer root of the product, then at runtime, divide the result by 6000 to get the root. Instead of treating the output as a*10^b, it uses the relation sqr(x*b^2)=sqr(x)*b. The math checks out:

• 42*6000*6000=1512000000
• sqr(1512000000)=38884
• 38884/6000=6.4806 (squared is 41.999)