用小数写自己的指数幂函数 [英] Writing your own exponential power function with decimals

问题描述

因此，我想使用某种算法在代码中编写一个函数，以计算任意幂的任何数字，包括小数.我使用JavaScript，它已经具有内置的pow函数:

So, I want to write a function in code using some sort of algorithm to calculate any number to any power, including decimals. I use JavaScript and it already has an inbuilt pow function:

Math.pow(2, 0.413) // 2^0.413 = 1.331451613236371, took under 1 second.

现在我要这样写自己的东西:

Now I want to write my own like this:

function pow(x, y) {
    // Algorithm
}

这是一个可以计算任何数字(x ^ 0.5)的平方根的函数，并且只有10个循环才非常准确:

This is a function that calculates the square root of any number (x^0.5), and it's very accurate with only 10 loops:

function sqrt(x, p) { // p = precision (accuracy)
    var a = 1;
    var b = x;

    while (p--) {
        a = (a + b) / 2
        b = x / a
    }

    return a
}

有没有简单的公式可以计算指数?

Is there any simple formula to calculate any exponential?

如果没有简单的东西，有没有困难的东西?

If there isn't a simple one, is there a hard one?

如果解决方案很慢，那么如何在一秒钟内估算JavaScript的功耗?

If the solution is slow, how can JavaScript's pow estimate under a single second?

推荐答案

这里有一个很好的正整数幂算法，它从处理一些简单的情况开始，然后使用循环测试指数的二进制位.例如，以二进制找到3^11 11就是1011，因此循环中的阶段是

Heres a nice algorithm for positive integer powers, it starts by dealing with some simple cases and then uses a loop testing the binary bits of the exponent. For example to find 3^11 11 in binary is 1011 so the stages in the loop are

• bitMask = 1011，evenPower = 3，结果= 3
• bitMask = 101，evenPower = 3 * 3 = 9，结果= 3 * 9 = 27
• bitMask = 10，evenPower = 9 * 9 = 81，结果= 27
• bitMask = 1，evenPower = 81 * 81 = 6561，结果= 27 * 6561 = 177147

这是每个循环的evenPower平方，如果最低位为1，则结果乘以evenPower.代码已从Patricia Shanahan的方法

That is the evenPower squares at each loop, and the result gets multiplied by the evenPower if the bottom bit is 1. The code has been lifted from Patricia Shanahan’s method http://mindprod.com/jgloss/power.html which in turn has its roots in Kunth and can be traced back to 200 BC in india.

/**
 * A fast routine for computing integer powers.
 * Code adapted from {@link <a href="http://mindprod.com/jgloss/power.html">efficient power</a>} by Patricia Shanahan [e-mail]
 * Almost identical to the method Knuth gives on page 462 of The Art of Computer Programming Volume 2 Seminumerical Algorithms.
 * @param l number to be taken to a power.
 * @param n power to take x to. 0 <= n <= Integer.MAX_VALUE
 * Negative numbers will be treated as unsigned positives.
 * @return x to the power n
 * 
 */
public static final double power(double l,int n)
{
    assert n>=0;

    double x=l;
    switch(n){
    case 0: x = 1.0; break;
    case 1: break;
    case 2: x *= x; break;
    case 3: x *= x*x; break;
    case 4: x *= x; x *= x; break;
    case 5: { double y = x*x; x *= y*y; } break;
    case 6: { double y = x*x; x = y*y*y; } break;
    case 7: { double y = x*x; x *= y*y*y; } break;
    case 8: x *= x; x *= x; x *= x; break;
    default:
    {
        int bitMask = n;
        double evenPower = x;
        double result;
        if ( (bitMask & 1) != 0 )
            result = x;
        else
            result = 1;
        bitMask >>>= 1;
        while ( bitMask != 0 ) {
            evenPower *= evenPower;
            if ( (bitMask & 1) != 0 )
                result *= evenPower;
            bitMask >>>= 1;
        } // end while
        x = result;
    }
    }
    return x;
}

对于真正的指数，您基本上需要找到exp和log的方法.您可以使用最简单的泰勒级数，但是有更好的方法.我们有

For a real exponent you will basically need ways of finding exp and log. You can use Taylor series which are the simplest to get but there are much better method. We have

exp(x) = 1 + x + x^2/2 + x^3/6 + x^4/24 + x^5/120 + x^6/6! + ...

ln(1+x) = x - x^2 /2 + x^3 /3 - x^4 / 4 + x^5 / 5 - x^6/6 + ... |x|<1

查找x ^ y注释ln(x^y) = y*ln(x).现在我们需要将参数设置在正确的范围内，以便可以使用幂级数.令x ＝ m * 2 ^ ex，选择尾数和指数，使得1/sqrt(2)＜＝ m ＜2. sqrt(2)和ln(m*2^ex) = ln(m) + ex*ln(2).令h = m-1并找到ln(1 + h).

To find x^y note ln(x^y) = y*ln(x). Now we need to get the argument in the right range so we can use our power series. Let x = m * 2^ex, the mantissa and exponent chosen so 1/sqrt(2)<= m < sqrt(2) and ln(m*2^ex) = ln(m) + ex*ln(2). Let h = m-1 and find ln(1+h).

使用java和floats是因为有一种简单的方法可以找到IEEE表示的内部结构(我已经使用float了，因为要处理的位数更少了)

Using java and floats as there is an easy way to find the internals of the IEEE representation (I've used float as there are fewer bits to cope with)

int b = Float.floatToIntBits(x);
int sign = (b & 0x80000000) == 0 ? 1 : -1;
int mattissa = b & 0x007fffff;
int ex = ((b & 0x7f800000) >> 23 ) - 127;

在javascript中，这对我们来说可能是最简单的 Number.toExponential 并分析结果.

in javascript it might be easiest to us Number.toExponential and parse the results.

接下来构造期望范围1/sqrt(2)＜的数字z. & sqrt(2)

Next construct a number z in the desired range 1/sqrt(2) < z < sqrt(2)

int bits = mattissa | 0x3f800000;
float z = Float.intBitsToFloat(bits);
if(z>root2) { 
    z = z/2;
    ++ex;
}

使用此函数使用泰勒级数查找1 + x的对数

Use this function to find the log of 1+x using a taylor series

static float ln1px(float x) {
    float x_2 = x*x; // powers of x
    float x_3 = x_2 * x;
    float x_4 = x_3 * x;
    float x_5 = x_4 * x;
    float x_6 = x_5 * x; 
    float res = x - x_2 /2 + x_3 /3 - x_4 / 4 + x_5 / 5 - x_6/6;
    return res;
}

这似乎对三个有效数字都很好，通常在x接近0时好得多.

this seems to be good to three significant figures, often much better when x is close to 0.

然后可以找到我们的数字x的日志

The log of our number x can then be found

float w = z - 1;
float ln_z = ln1px(w);
float ln_x = ln_z + ln2 * ex;
System.out.println("ln "+ln_x+"\t"+Math.log(x));

现在，如果我们写y = n + a，则为实际功效，其中n是整数，而a是小数.所以 x^y=x^(n+a) = x^n * x^a.使用此答案中的第一个算法找到x^n.先写x=m*2^ex然后写ln((m*2^ex)^a) = yln(m) + yex*ln(2)和

Now to the actual power if we write y = n + a where n is an integer and a is fractional. So x^y=x^(n+a) = x^n * x^a. use the first algorithm in this answer to find the x^n. Writing x=m*2^ex then ln((m*2^ex)^a) = yln(m) + yex*ln(2) and

x^a=exp(ln((m*2^ex)^a)) = exp(a * ln(m)) * exp(a * ln(2))^ex

两个指数项的值都较小，因此taylor级数应该很好.

the two exponential terms have fairly small values so the taylor series should be good.

指数函数的泰勒级数需要一个函数

We need one function for the taylor series of the exponential function

static float exp(float x) {
    float x_2 = x*x; // powers of x
    float x_3 = x_2 * x;
    float x_4 = x_3 * x;
    float x_5 = x_4 * x;
    float x_6 = x_5 * x; 
    float res = 1+ x + x_2 /2 + x_3 /6 + x_4 / 24 + x_5 / 120 + x_6/ 720;
    return res;
}

最后我们可以将各个部分放在一起

finally we can put the pieces together

// Get integer and fractional parts of y
int n = (int) Math.floor(y);
float a = y-n;

float x_n = power(x,n);         // x^n
float a_ln_m = a * ln_z;        // ln(m^a) = a ln(m)
float a_ln_2 = a * ln2;         // ln(2^a) = a ln(2)
float m_a = exp(a_ln_m);        // m^a = exp(a ln(m))
float _2_a = exp(a_ln_2);       // 2^a = exp(a ln(2))
float _2_a_ex = power(_2_a,ex); // (2^ex)^a = 2^(a*ex) = (2^a)^ex 
float x_a = m_a * _2_a_ex;      // x^a = m^a * 2^(a*ex)

float x_y = x_n * x_a;          // x^y = x^n * x^a

System.out.println("x^y "+x_y+"\t"+Math.pow(x,y));

那应该是完整的程序，您需要一些技巧来应对否定参数等.

That should be the complete program, you need some smarts to cope with negative arguments etc.

请注意，这并不是特别准确，因为我只使用了taylor系列的几个术语.其他SO问题具有更详细的答案我如何自己编写幂函数?

Note this is not particularly accurate as I've only used a few terms of the taylor series. Other SO questions have more detailed answers How can I write a power function myself?

这篇关于用小数写自己的指数幂函数的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！