math库的Python实现原理（根号运算）

本系列主要使用Python 实现，主要展示计算机基本数学运算是如何实现的，对于Python中math底层应该是c写的，所以直接在库里查不到源代码，其中所有内容均是查阅的资料，不一定是最高效的，但会尽量去找最高效的代码。

下面是根号运算，数学公式使用牛顿迭代法，其中迭代条件为精度，可根据需要修改，大量提高精度不会额外消耗过多时间。

def mysqrt(x):
    val = x
    last = 0.0
    while (abs(val - last) > 0.00000001):
        last = val
        val = (val + x / val) / 2
    return val

下面是多次根运算n=2时等同于上面，并且可以为小数，当n=1/a时，该方法等同于 x^a，同样支持小数

def mysqrt1(x, n):
    val = x
    last = 0.0
    while (abs(val - last) > 0.00000001):
        last = val
        val = ((n-1)*val + x / val**(n-1)) / n
    return val

测试

print(1.6**3)
print(mysqrt1(1.6, 1/3))
print(mypow(1.6, 3))
print(1.6**0.25)
print(mysqrt1(1.6, 4))

结果

4.096000000000001
4.096000000000004
4.096000000000004
1.1246826503806981

1.1246826503806981

整理后

Epsilon = 10e-16

def fab_h(x):
    '''
    求实数的绝对值
    :param x: R
    '''
    if x >= 0:
        return x
    else:
        return x * -1

def sqrt_h(x, n=2.0):
    '''
    n倍根号下x
    牛顿迭代法
    '''
    val = x
    last = 0.0
    if n == 2.0:
        while (fab_h(val - last) > Epsilon):
            last = val
            val = (val + x / val) / 2
        return val
    while (fab_h(val - last) > Epsilon):
        last = val
        val = ((n-1)*val + x / val**(n-1)) / n
    return val

总结：这是一种相当快的运算方法，只要迭代几次就可以达到相当的精度，比如上面的8位精度都可以在5次迭代内算出来。

不过可以看出来根号运算相当于次方运算的逆运算，x^1/4=x开4次方根，对于上面的运算来说实际上引用了x^n方法，只有根号运算是独立的，其他都可以用x^n方法替代。这里可以发现数学各种公式其实都是有关系的。

下面再提供一种C语音的根号实现，来源于雷神之锤3D引擎源代码，是现在最快的根号算法实现，由于C语言的指针特性，我没办法用python实现，就直接贴源代码了，最关键的是不知道常数0x5f375a86是如何想到的，不得不说，数学的顶端是玄学。

float InvSqrt(float x)
{
	float xhalf = 0.5f*x;
	int i = *(int*)&x; // get bits for floating VALUE 
	i = 0x5f375a86- (i>>1); // gives initial guess y0
	x = *(float*)&i; // convert bits BACK to float
	x = x*(1.5f-xhalf*x*x); // Newton step, repeating increases accuracy
	return x;
}