Efficient 10 power x algorithm

Question

Can anyone help me with finding an efficient code to find 10 power x?

 private int power(int base, int exp)
{
    int result = 1;
    while (exp != 0)
    {
        if ((exp & 1) == 1)
            result *= base;
        exp >>= 1;
        base *= base;
    }

    return result;
}

Source of code from here, but I am looking for a way where the input could be 3.14 (double). I also cannot use any library functions. The power can be a real number. So it is not just a simple integer algorithm where we can find by Exponentiation by Squaring.

Answer 1

Why dont you use Math.pow(double, double)

you can even check out the source.

---

Answer 2

Here's the pow method from a Java StrictMath implementation. Technically, StrictMath favors accuracy over performance, but in practice it's pretty fast as it's used everywhere in Java I copied and formatted the sources necessary for the pow function for you :) :

public class YourMathClass {
    public static double pow(double x, double y) {
        // Special cases first.
        if (y == 0)
            return 1;
        if (y == 1)
            return x;
        if (y == -1)
            return 1 / x;
        if (x != x || y != y)
            return Double.NaN;
        // When x < 0, yisint tells if y is not an integer (0), even(1),
        // or odd (2).
        int yisint = 0;
        if (x < 0 && floor(y) == y)
            yisint = (y % 2 == 0) ? 2 : 1;
        double ax = abs(x);
        double ay = abs(y);
        // More special cases, of y.
        if (ay == Double.POSITIVE_INFINITY) {
            if (ax == 1)
                return Double.NaN;
            if (ax > 1)
                return y > 0 ? y : 0;
            return y < 0 ? -y : 0;
        }
        if (y == 2)
            return x * x;
        if (y == 0.5)
            return sqrt(x);
        // More special cases, of x.
        if (x == 0 || ax == Double.POSITIVE_INFINITY || ax == 1) {
            if (y < 0)
                ax = 1 / ax;
            if (x < 0) {
                if (x == -1 && yisint == 0)
                    ax = Double.NaN;
                else if (yisint == 1)
                    ax = -ax;
            }
            return ax;
        }
        if (x < 0 && yisint == 0)
            return Double.NaN;
        // Now we can start!
        double t;
        double t1;
        double t2;
        double u;
        double v;
        double w;
        if (ay > TWO_31) {
            if (ay > TWO_64) // Automatic over/underflow.
                return ((ax < 1) ? y < 0 : y > 0) ? Double.POSITIVE_INFINITY : 0;
            // Over/underflow if x is not close to one.
            if (ax < 0.9999995231628418)
                return y < 0 ? Double.POSITIVE_INFINITY : 0;
            if (ax >= 1.0000009536743164)
                return y > 0 ? Double.POSITIVE_INFINITY : 0;
            // Now |1-x| is <= 2**-20, sufficient to compute
            // log(x) by x-x^2/2+x^3/3-x^4/4.
            t = x - 1;
            w = t * t * (0.5 - t * (1 / 3.0 - t * 0.25));
            u = INV_LN2_H * t;
            v = t * INV_LN2_L - w * INV_LN2;
            t1 = (float) (u + v);
            t2 = v - (t1 - u);
        } else {
            long bits = Double.doubleToLongBits(ax);
            int exp = (int) (bits >> 52);
            if (exp == 0) // Subnormal x.
            {
                ax *= TWO_54;
                bits = Double.doubleToLongBits(ax);
                exp = (int) (bits >> 52) - 54;
            }
            exp -= 1023; // Unbias exponent.
            ax = Double.longBitsToDouble((bits & 0x000fffffffffffffL) | 0x3ff0000000000000L);
            boolean k;
            if (ax < SQRT_1_5)  // |x|<sqrt(3/2).
                k = false;
            else if (ax < SQRT_3) // |x|<sqrt(3).
                k = true;
            else {
                k = false;
                ax *= 0.5;
                exp++;
            }
            // Compute s = s_h+s_l = (x-1)/(x+1) or (x-1.5)/(x+1.5).
            u = ax - (k ? 1.5 : 1);
            v = 1 / (ax + (k ? 1.5 : 1));
            double s = u * v;
            double s_h = (float) s;
            double t_h = (float) (ax + (k ? 1.5 : 1));
            double t_l = ax - (t_h - (k ? 1.5 : 1));
            double s_l = v * ((u - s_h * t_h) - s_h * t_l);
            // Compute log(ax).
            double s2 = s * s;
            double r = s_l * (s_h + s) + s2 * s2 * (L1 + s2 * (L2 + s2 * (L3 + s2 * (L4 + s2 * (L5 + s2 * L6)))));
            s2 = s_h * s_h;
            t_h = (float) (3.0 + s2 + r);
            t_l = r - (t_h - 3.0 - s2);
            // u+v = s*(1+...).
            u = s_h * t_h;
            v = s_l * t_h + t_l * s;
            // 2/(3log2)*(s+...).
            double p_h = (float) (u + v);
            double p_l = v - (p_h - u);
            double z_h = CP_H * p_h;
            double z_l = CP_L * p_h + p_l * CP + (k ? DP_L : 0);
            // log2(ax) = (s+..)*2/(3*log2) = exp + dp_h + z_h + z_l.
            t = exp;
            t1 = (float) (z_h + z_l + (k ? DP_H : 0) + t);
            t2 = z_l - (t1 - t - (k ? DP_H : 0) - z_h);
        }
        // Split up y into y1+y2 and compute (y1+y2)*(t1+t2).
        boolean negative = x < 0 && yisint == 1;
        double y1 = (float) y;
        double p_l = (y - y1) * t1 + y * t2;
        double p_h = y1 * t1;
        double z = p_l + p_h;
        if (z >= 1024) // Detect overflow.
        {
            if (z > 1024 || p_l + OVT > z - p_h)
                return negative ? Double.NEGATIVE_INFINITY : Double.POSITIVE_INFINITY;
        } else if (z <= -1075) // Detect underflow.
        {
            if (z < -1075 || p_l <= z - p_h)
                return negative ? -0.0 : 0;
        }
        // Compute 2**(p_h+p_l).
        int n = round((float) z);
        p_h -= n;
        t = (float) (p_l + p_h);
        u = t * LN2_H;
        v = (p_l - (t - p_h)) * LN2 + t * LN2_L;
        z = u + v;
        w = v - (z - u);
        t = z * z;
        t1 = z - t * (P1 + t * (P2 + t * (P3 + t * (P4 + t * P5))));
        double r = (z * t1) / (t1 - 2) - (w + z * w);
        z = scale(1 - (r - z), n);
        return negative ? -z : z;
    }

    public static int round(float f) {
        return (int) floor(f + 0.5f);
    }

    public static double floor(double a) {
        double x = abs(a);
        if (!(x < TWO_52) || (long) a == a)
            return a; // No fraction bits; includes NaN and infinity.
        if (x < 1)
            return a >= 0 ? 0 * a : -1; // Worry about signed zero.
        return a < 0 ? (long) a - 1.0 : (long) a; // Cast to long truncates.
    }

    public static double abs(double d) {
        return (d <= 0) ? 0 - d : d;
    }

    public static double sqrt(double x) {
        if (x < 0)
            return Double.NaN;
        if (x == 0 || !(x < Double.POSITIVE_INFINITY))
            return x;
        // Normalize x.
        long bits = Double.doubleToLongBits(x);
        int exp = (int) (bits >> 52);
        if (exp == 0) // Subnormal x.
        {
            x *= TWO_54;
            bits = Double.doubleToLongBits(x);
            exp = (int) (bits >> 52) - 54;
        }
        exp -= 1023; // Unbias exponent.
        bits = (bits & 0x000fffffffffffffL) | 0x0010000000000000L;
        if ((exp & 1) == 1) // Odd exp, double x to make it even.
            bits <<= 1;
        exp >>= 1;
        // Generate sqrt(x) bit by bit.
        bits <<= 1;
        long q = 0;
        long s = 0;
        long r = 0x0020000000000000L; // Move r right to left.
        while (r != 0) {
            long t = s + r;
            if (t <= bits) {
                s = t + r;
                bits -= t;
                q += r;
            }
            bits <<= 1;
            r >>= 1;
        }
        // Use floating add to round correctly.
        if (bits != 0)
            q += q & 1;
        return Double.longBitsToDouble((q >> 1) + ((exp + 1022L) << 52));
    }

    private static double scale(double x, int n) {
        //       if (Configuration.DEBUG && abs(n) >= 2048)
        //         throw new InternalError("Assertion failure");
        if (x == 0 || x == Double.NEGATIVE_INFINITY || !(x < Double.POSITIVE_INFINITY) || n == 0)
            return x;
        long bits = Double.doubleToLongBits(x);
        int exp = (int) (bits >> 52) & 0x7ff;
        if (exp == 0) // Subnormal x.
        {
            x *= TWO_54;
            exp = ((int) (Double.doubleToLongBits(x) >> 52) & 0x7ff) - 54;
        }
        exp += n;
        if (exp > 0x7fe) // Overflow.
            return Double.POSITIVE_INFINITY * x;
        if (exp > 0) // Normal.
            return Double.longBitsToDouble((bits & 0x800fffffffffffffL) | ((long) exp << 52));
        if (exp <= -54)
            return 0 * x; // Underflow.
        exp += 54; // Subnormal result.
        x = Double.longBitsToDouble((bits & 0x800fffffffffffffL) | ((long) exp << 52));
        return x * (1 / TWO_54);
    }

    private static final double TWO_31 = 0x80000000L, // Long bits 0x41e0000000000000L.
            TWO_52 = 0x10000000000000L, // Long bits 0x4330000000000000L.
            TWO_54 = 0x40000000000000L, // Long bits 0x4350000000000000L.
            TWO_64 = 1.8446744073709552e19; // Long bits 0x7fe0000000000000L.
    private static final double L1 = 0.5999999999999946, // Long bits 0x3fe3333333333303L.
            L2 = 0.4285714285785502, // Long bits 0x3fdb6db6db6fabffL.
            L3 = 0.33333332981837743, // Long bits 0x3fd55555518f264dL.
            L4 = 0.272728123808534, // Long bits 0x3fd17460a91d4101L.
            L5 = 0.23066074577556175, // Long bits 0x3fcd864a93c9db65L.
            L6 = 0.20697501780033842, // Long bits 0x3fca7e284a454eefL.
            P1 = 0.16666666666666602, // Long bits 0x3fc555555555553eL.
            P2 = -2.7777777777015593e-3, // Long bits 0xbf66c16c16bebd93L.
            P3 = 6.613756321437934e-5, // Long bits 0x3f11566aaf25de2cL.
            P4 = -1.6533902205465252e-6, // Long bits 0xbebbbd41c5d26bf1L.
            P5 = 4.1381367970572385e-8, // Long bits 0x3e66376972bea4d0L.
            DP_H = 0.5849624872207642, // Long bits 0x3fe2b80340000000L.
            DP_L = 1.350039202129749e-8, // Long bits 0x3e4cfdeb43cfd006L.
            OVT = 8.008566259537294e-17; // Long bits 0x3c971547652b82feL.
    private static final double SQRT_1_5 = 1.224744871391589, // Long bits 0x3ff3988e1409212eL.
            SQRT_3 = 1.7320508075688772, // Long bits 0x3ffbb67ae8584caaL.
            CP = 0.9617966939259756, // Long bits 0x3feec709dc3a03fdL.
            CP_H = 0.9617967009544373, // Long bits 0x3feec709e0000000L.
            CP_L = -7.028461650952758e-9, // Long bits 0xbe3e2fe0145b01f5L.
            LN2 = 0.6931471805599453, // Long bits 0x3fe62e42fefa39efL.
            LN2_H = 0.6931471803691238, // Long bits 0x3fe62e42fee00000L.
            LN2_L = 1.9082149292705877e-10, // Long bits 0x3dea39ef35793c76L.
            INV_LN2 = 1.4426950408889634, // Long bits 0x3ff71547652b82feL.
            INV_LN2_H = 1.4426950216293335, // Long bits 0x3ff7154760000000L.
            INV_LN2_L = 1.9259629911266175e-8; // Long bits 0x3e54ae0bf85ddf44L.
}

Now, if you're really into making things go really fast for large-scale computations (I mean hundreds of thousands of pow method calls), may I suggest this extremely fast pow function from the Apache Commons Math FastMath class?

According to your question, you don't want to have another library, so I copied the necessary parts of the FastMath and FastMathLiteralArrays in one file for you. It's too big to put on StackOverflow, so I put up on PasteAll:

http://www.pasteall.org/53238/java

WARNING: This file is huge in terms of code. Unless you are sure you are doing extremely long computations using pow, you shouldn't use this file because 1) It will use a ton of memory and 2) It may actually go slower if you are just going to use it once or twice.

Anyways, hope this answers all you pow efficiency needs :) (I also hope I'll get the bounty ;) )

Answer 3

Catering for double input without using library code is pointless - you either have to use libraries or plagiarise the source for library code - same thing.

However for integer arithmetic, not using libraries is plausible.

Answering your question as stated and titled:

private static int tenToPower(int power) {
    int result = 1;
    for (int p = 0; p < power; p++)
        result *= 10;
    return result;
}

Since your return type is int in your example code, I have assumed that negative powers, which would produce fractional results, are out of scope.

Answer 4

Instead of calculating 10^x you can calculate e^(x*ln(10)). So your question boil down to this one. I can offer a simple implementation by calculating the Taylor series:

double tenPow(double x) {
    // ln(10)
    double logTen = 2.3025850929940456840179914546843642076011014886287729;
    double sum = 0.;
    x *= logTen; // x * ln(10)
    double tmp = 1.;
    for (double i = 1.; tmp > 0.; i += 1.) {
        sum += tmp;
        tmp *= x;
        tmp /= i;
    }
    return sum;
}

EDIT1 As noted in the comment this method is not effective for large x. But we can break x into its integer and fractional parts, so x = xi + xf. e^(xi + xf) = e^xf * e^xi. e^xi can be effectively calculated using recursive exponentiation. And e^xf can be calculated using the Taylor series.

EDIT2 Here is my implementation of this algorithm:

public static double fastTenPow(double x) {
    if (x < 0.) {
        return 1. / fastTenPow(-x);
    }
    int intPart = (int) x;
    double fractPart = x - intPart;
    return fractTenPow(fractPart) * intTenPow(intPart);
}

private static double intTenPow(int x) { // copied it
    double res = 1.;
    double base = 10.;
    while (x != 0) {
        if ((x & 1) == 1) {
            res *= base;
        }
        x >>= 1;
        base *= base;
    }
    return res;
}

private static final double LOG_TEN = 2.3025850929940456840179914546843642076011014886287729;
private static final double EPS = 0.00000000000000001;

private static double fractTenPow(double x) {
    double sum = 0.;
    x *= LOG_TEN;
    double tmp = 1.;
    for (double i = 1.; tmp > EPS; i += 1.) {
        sum += tmp;
        tmp *= x;
        tmp /= i;
    }
    return sum;
}

Answer 5

If you want really fast and efficient power of 10, and only need integers, you can use a table (really old trick from the 90s), like this:

static final long powerOfTen[] = {1,10,100,1000,10000, ....} // you get the idea

and then

static long powerOfTen(long input, int power) {
    return input * powerOfTen[power];
}

If you need greater precision, you can also pre-calculate that, but keep in mind as the table gets bigger, you win less and less since it uses more memory/L1 cache. The point above is that because of its really small size the static method would be inlined, and the static array is so small that if fits a 64 byte L1 cache page, meaning that performance-wise, this would come down to an integer multiplication on registers.

Answer 6

can't you use a simple recursion?

 public static void main(String[] args) {
    System.out.println(power10x(4,1)); //replace ur X value with 4

}


public static int power10x(int a,int val){
     if(a==0)
        return val;
     else
         return power10x(a-1,val*10);

}