4

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)

Donkey King
  • 87
  • 1
  • 8
  • That depend on what you are allowed to use. – zch Dec 30 '16 at 22:59
  • 3
    `sqrt(32)` – Or use `log` and `exp`. But why do you want to avoid `pow`? – Martin R Dec 30 '16 at 23:00
  • I want to avoid any given functions. For example when calculating integer exponents you can just make a loop and multiply the number with itself. Is there a similar way to do that with non integer exponents? – Donkey King Dec 30 '16 at 23:03
  • 1
    I would say using `pow` is the most neat solution you can implement using C but I'm happy to be corrected. – Marco Dec 30 '16 at 23:04
  • 2
    To put it another way, if `pow` wasn't already the neat way to achieve this, then somebody would have fixed it by now. – user3386109 Dec 30 '16 at 23:07
  • 2
    `pow()` is a _very challenging_ function to code accurately and has many special cases. Try coding `sqrt()` if you like to do things manually, then `log()`, then `exp()` for increasingly difficult functions. Save `pow()` for later - much later. – chux - Reinstate Monica Dec 31 '16 at 00:10

4 Answers4

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.

Brandon
  • 22,723
  • 11
  • 93
  • 186
  • 1
    Note that the question is tagged as C, not C++. That makes the tangential parts of your code (headers, calls prefixed with `std::`) inapposite, though readily fixable. What you're suggesting works OK while the exponent is a multiple of a half; it doesn't work for more general exponents. There are more efficient algorithms for the `power()` function. – Jonathan Leffler Dec 30 '16 at 23:50
  • Code has trouble with `power(base, negative_number)`. Often returns 1.0. – chux - Reinstate Monica Dec 31 '16 at 00:24
  • @chux; Fixed. x ^ -n = 1 / (x ^ n). – Brandon Dec 31 '16 at 03:46
  • Good improvement, UV. Note that `d >= DBL_EPSILON * 10` is a good _absolute_ test, yet FP limits should use a _relative_ test. `pow()` is a difficult function to code – chux - Reinstate Monica Dec 31 '16 at 05:13
0

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.

Ray
  • 686
  • 3
  • 6
  • Are you referring to `float powf(float, float)`, which is in standard C? In fact, the `pow()`, `powf()` and `powl()` functions all take two values of the same type and return another value of the same type — there are no versions in standard C that take integers for the second argument (or first argument, come to that). – Jonathan Leffler Dec 30 '16 at 23:46
  • With `a == -2, b == 4`, then `exp( b * log (a) )` fails rather than 16.0 – chux - Reinstate Monica Dec 31 '16 at 00:15
  • @chux: `(a <> 0.) ?((a > 0.) ?1. :(-1.)) * exp(b * log(fabs(a)) :0.` :) – alk Dec 31 '16 at 15:58
  • @alk Rather than `power(0.0, any_b)--> 0.0` per your [comment](http://stackoverflow.com/questions/41403730/c-exponentiation-with-non-integer-exponent/41403800?noredirect=1#comment70022571_41403800), `power(0.0, any_pos_b)--> 0.0, power(0.0, any_neg_b)--> inf`, leaving `power(0.0, 0.0)--> TBD`, (-0.0 may incur additional issues). As [commented](http://stackoverflow.com/questions/41403730/c-exponentiation-with-non-integer-exponent/41403800?noredirect=1#comment70011161_41403730), `pow()` is a challenging function - many corner cases. :) – chux - Reinstate Monica Dec 31 '16 at 17:17
  • thank you for taking care of the edge cases. I was a bit lazy by only giving the general approach, not dealing with edge cases :) – Ray Jan 23 '17 at 17:59
0

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;
}
Malcolm McLean
  • 6,258
  • 1
  • 17
  • 18
0

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

thodnev
  • 1,564
  • 16
  • 20