I'm making program that prints their nth power for numbers from 1 to 20, using recursive and non-recursive function. I would like to give an option for user to enter real numbers. How to do that? Can it be done without pow function? Thank you for your help!
CodePudding user response:
If you want to raise a number to a real exponent x
a^x
you can compute it as
e^(x ln a)
so you can write a function like
double mypow( double base, double exponent )
{
return exp( exponent * log( base ) )
}
If you want to get deep into the weeds and create your own implementations of exp
and log
, you can use the Taylor series
e^x = 1 x ((x^2)/2!) ((x^3)/3!) ...
and
ln x = 2( u (u^3)/3 (u^5)/5 ... ), u = (x 1)/(x-1)
Personally, I would just use the pow
function.
CodePudding user response:
In making double keiski_pow(double num, double exponent)
consider different code for certain conditions.
// exponent == 0, num != 0
return 1
// exponent == 0, num == 0
Domain error
// num > 0
return exp(log(num)*exponent)
// or more accurately
return exp2(log2(num)*exponent)
// num == 0, exponent < 0
Infinity
// num < 0, exponent is an integer > 0.
// Integer test for `double`: trunc(exponent) == exponent
t = keiski_pow(num, exponent/2)
return t*t*(exponent%2 ? num : 1)
// or
t = keiski_pow(fabs(num), exponent)
return t*exponent%2 ? -1 : 1)
// num < 0, exponent is an integer < 0
return 1.0/keiski_pow(num, -exponent)
// num < 0, exponent is not an integer
Domain error
I think this covers most cases. (NaN are another story)
Note: recursion()
is very very slow for large exponent
. Try recursion(1.0, 1e15)
.
A much faster recursive solution:
double recursion_alt(double num, double exponent) {
// Code only valid when exponent has a whole number value
assert(exponent = trunc(exponent));
if (exponent > 0.0) {
double exponent_mod2 = fmod(exponent, 2.0);
double exponent_div2 = (exponent - exponent_mod2)/2.0;
double t = recursion_alt(num, exponent_div2);
return t*t*(exponent_mod2 ? num : 1.0);
}
if (exponent < 0.0) {
return 1.0/recursion_alt(num, -exponent);
}
return 1.0;
}