本系列主要使用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;
}