C - Exponentiation with non integer exponent

Question

Is there a neat way to achieve this without using the pow()-function? Say wanting to calculate 2^(2.5).

EDIT: Perhaps pow() is the way to go after all. I was hoping to create a function of my own using only the four common operations to solve it. (Reason being I like to do things manually)

Answer 1

There really should be a pow(float, float) version of that function. If not, you can express a**b by calculating exp( b * log (a) ), see wikipedia.

Answer 2

x^(2.5) = x * x * sqrt(x)

Add the multiplications to get the exponent. If follows that if you need near fractions of powers of 2, you can build it on top of sqrt. x^1/4 = sqrt(sqrt(x))

Euler's number is special and e^x can be calculated from a quickly converging factorial series. pow is implemented on top of it.

This is code from my book Basic Algorithms (meant to explain how to do the calculations, not as optimised production code).

/*
   power function
 */
double power(double x, double y)
{
  int i;
  double answer = 1;

  if(y < 0)
     return 1/power(x, -y);
  if(y == (int) y)
  {
    for(i=0;i<y;i++)
      answer *= x;
    return answer;
  }
  if(x == 0)
    return 0;
   assert(x > 0);
   return exponent(y * logarithm(x));
}

/*
  log function
 */
double logarithm(double x)
{
  int intpart = 0;
  double fractpart;
  double power;
  int i;

   while(x > 3.0/2.0)
   {
      intpart++;
      x /= 2.718282;
    }
   while(x < 1.0/2.0)
   {
      intpart--;
      x *= 2.718282;
    }

    x -= 1.0;

   i = 1;
   power = 1;
  fractpart = 0;
  do
  {
       power *= x;
       fractpart += power/i * ((i % 2) ? 1 : -1);
       i++;
    }
   while(power/i > 0.00000001 || power/i < -0.00000001);

   return intpart + fractpart;
}

Answer 3

• A ^ (B/C) is the same as CthRoot(A ^ B).
• A ^ (1/2) is the same as 2ndRoot(A ^ 1); Same as Sqrt(A ^ 1).
• A ^ (3/4) is the same as 4thRoot(A ^ 3).
• 2 ^ (2.5) is the same as 2 ^ (5 / 2) is the same as Sqrt(2 ^ 5).

Test Run: http://ideone.com/ILlT85

Simple Naive code:

#include <stdio.h>
#include <math.h>
#include <stdint.h>

#ifndef DBL_EPSILON
#define DBL_EPSILON 2.2204460492503131e-16
#endif

float power(float a, float b)
{
    int isNegative = b < 0.0;


    float res = 1;

    for (int i = 0; i < fabs(b); ++i)
    {
        res *= a;
    }

    return isNegative ? 1 / res : res;
}

double root(int n, double x)
{
    double d, r = 1;
    if (!x)
    {
        return 0;
    }

    if (n < 1 || (x < 0 && !(n&1)))
    {
        return 0.0 / 0.0;
    }

    do
    {
        d = (x / power(r, n - 1) - r) / n;
        r += d;
    } while (d >= DBL_EPSILON * 10 || d <= -DBL_EPSILON * 10);

    return r;
}


long gcd(long a, long b)
{
    return b == 0 ? a : gcd(b, a % b);
}

void frac(float value, long* numerator, long* denominator)
{
    double integral = floor(value);
    double frac = value - integral;
    const long precision = 1000000;

    long commonDenominator = gcd(round(frac * precision), precision);
    *numerator = round(frac * precision) / commonDenominator;
    *denominator = precision / commonDenominator;

    *numerator = *numerator + (integral * *denominator);
}

int main() {

    float base = 2;
    float exp = 2.5;


    printf("FIRST: %f\n", pow(base, exp));

    //OR

    //A ^ (B/C) is the same as CthRoot(A ^ B)
    long num = 0;
    long den = 0;
    frac(exp, &num, &den);

    printf("SECOND: %f\n", root(den, power(base, num)));


    base = 3;
    exp = 2.7;

    printf("THIRD: %f\n", pow(base, exp));

    //OR:

    num = 0;
    den = 0;
    frac(exp, &num, &den);

    printf("FOURTH: %f\n", root(den, power(base, num)));

    return 0;
}

Essentially, you can write a Root function that calculates the Nth-Root of a value and a power function. You'll also need a function that creates fractions from floats.

This is probably a lot more work than other solutions.

Answer 4

You could calculate it approximately, if this fits your needs: what you need is called Taylor binomial series, and that's how pow function could be internally implemented (actually, it usually uses better converging series). All that you need is determine the number of terms, which would allow you to fit the approximation error needed. Actually, the error of N-terms series expansion would be no more, than N+1th term (if it converges).

For more, have a look at this