因?yàn)榇邓哪芰Σ患雅敉砸却騻€草稿冷蚂，今天的吹水過程大概是：
1、牛頓迭代法的演繹過程
2痘番、牛頓迭代法求n次方根
3捉片、牛頓迭代法求n次方根改進(jìn)版
4、牛逼哄哄的invsqrt求平方根倒數(shù)

1汞舱、牛頓迭代法的演繹過程

乍一聽伍纫，好像很高大上，其實(shí)理解上沒什么難度兵拢。
牛頓迭代法就是在不斷迭代的過程中逼近曲線的根翻斟，可以看做，它就是用來求曲線的根的说铃。
所以它是怎么個迭代的過程访惜，先上個動圖：

牛頓迭代法的過程

我用吹水的姿勢描述一遍動圖的過程，先是有一個曲線腻扇，然后要找到曲線和x軸的交點(diǎn)债热，這個時候，腦袋一拍幼苛，選x1作為第一點(diǎn)窒篱，在這個點(diǎn)上垂直于x軸的線與曲線的交點(diǎn)，在這個交點(diǎn)做與曲線的切線，而切線與x軸的交點(diǎn)就是下一個迭代的點(diǎn)墙杯，不斷重復(fù)這個過程配并，直到找到曲線和x軸交點(diǎn)這個目標(biāo)值。值得注意的是高镐，這個方法是不斷逼近溉旋，所以得到的值可能會與目標(biāo)值存在誤差，而誤差小于一定范圍就可以認(rèn)為是目標(biāo)值了嫉髓。一個更加圖文并茂更加通俗易懂的解釋可以點(diǎn)擊 這個這個這個這個观腊。

2、牛頓迭代法求n次方根

那為什么說牛頓迭代法能用來求n次根呢算行，其實(shí)在1中的目標(biāo)值梧油，即曲線與x軸相交點(diǎn)的x值，假設(shè)你的方程為 f(x) = x2 - 64州邢，當(dāng)f(x) = 0時儡陨，x就是64的2次方根，以此類推量淌，把2換成你想要求的n次方迄委，就可以求出64的n次方根了。

明白了演繹過程类少，那我們就來探索它的代數(shù)實(shí)現(xiàn)。在上代碼前渔扎，先上個圖

牛頓迭代法某一點(diǎn)的求值圖解

每次迭代要求的點(diǎn)就是右綠色箭頭指向的這個點(diǎn)硫狞，在圖中可以看到黑色箭頭標(biāo)識的線段長度為-f(x)，f'(x)是對f(x)的求導(dǎo)結(jié)果晃痴，也就是切線的斜率残吩，所以綠色箭頭線段的長度為-f(x)/f'(x)，而左綠色箭頭的點(diǎn)為x倘核，所有右綠色箭頭指向的值為x加上綠色線段的長度泣侮，即x-(f(x)/f'(x))。

上面的代數(shù)實(shí)現(xiàn)完成后紧唱，就可以上代碼了活尊。嘻嘻嘻，最近在看Python版的SICP漏益，所以貼上Python的代碼蛹锰，這篇文章也是因?yàn)榭戳诉@本書的1.6章節(jié)寫的。def就是定義一個方法的意思绰疤。

def newton_nth_root_of_a(n, a, x = 1, e = 1e-13):
    def f(x):
        return x ** n - a
    def df(x):
        return n * x ** (n - 1)
    
    gap = abs(f(x))
    while gap > e:
        x = x - f(x)/df(x)
        gap = abs(f(x))
    return x
    
print(newton_nth_root_of_a(1, 64))  #64.0
print(newton_nth_root_of_a(2, 64))  #8.0
print(newton_nth_root_of_a(3, 64))  #4.0
print(newton_nth_root_of_a(4, 64))  #2.82842712474619
print(newton_nth_root_of_a(5, 64))  #2.29739670999407
print(newton_nth_root_of_a(6, 64))  #2.0

這里是在線可運(yùn)行版本 铜犬，不多不少，剛好十行，newton_nth_root_of_a方法定義中癣猾，n代表冪敛劝，a代表對a求根，x代表起始值纷宇，e代表誤差范圍夸盟，x和e已經(jīng)設(shè)置了缺省值。f(x)是曲線呐粘，df(x)是f(x)的求導(dǎo)結(jié)果满俗，在while循環(huán)中，x - f(x)/df(x)就是我們上面代數(shù)實(shí)現(xiàn)中的結(jié)果作岖，直到f(x)的值小于誤差唆垃，x就是最后的結(jié)果。

3痘儡、牛頓迭代法求n次方根改進(jìn)版

明白了1辕万、2其實(shí)整個過程就是這樣了，那為啥有3沉删，因?yàn)橐仓^皮接著吹水??渐尿。
2中的10行代碼很簡潔，但是假如我要不斷的求某些值的2次根矾瑰，就要不斷調(diào)用newton_nth_root_of_a(2, a)砖茸。那么我們來改進(jìn)下代碼，保持上面的代碼不變殴穴，增加下面的代碼：

def curry_nth_root(n):
    def curry_nth_root_of_a(a):
        return newton_nth_root_of_a(n, a)    
    return curry_nth_root_of_a

_2th_root = curry_nth_root(2)

print(_2th_root(64)) #8.0
print(_2th_root(49)) #7.000000000000001
print(_2th_root(36)) #6.0
print(_2th_root(25)) #5.0
print(_2th_root(16)) #4.0
print(_2th_root(9))  #3.0

這里是在線可運(yùn)行版本 凉夯，哇咔咔，在curry化之后采幌，就可以方便的求n次方根了劲够，不用每次都輸入n。

curry化其實(shí)也很簡單休傍，那么我們就來個極端點(diǎn)的抽象代碼：

def improve_guess(update, close, guess=1):
    '''返回猜想值征绎，update為更新方法，close為逼近方法磨取，guess為初始猜想值'''
    while not close(guess):
        guess = update(guess)
    return guess
    
def newton_update(f, df):
    '''返回牛頓迭代的計(jì)算函數(shù), f為原函數(shù)人柿，df是對其求導(dǎo)'''
    def update(x):
        return x - f(x) / df(x)
    return update
    
def newton_find_zero_of_f(f, df, tolerance=1e-13):
    '''找出f(x)=0時，x的值忙厌，tolerance為與目標(biāo)真實(shí)值的誤差'''
    def near_zero(x):
        return is_eq(f(x), 0, tolerance)
    return improve_guess(newton_update(f, df), near_zero)
    
def is_eq(x, y, tolerance):
    '''x與y的差值是否在誤差范圍內(nèi)'''
    return abs(x - y) < tolerance
    
def nth_root_of_a(n, a):
    def f(x):
        return x ** n - a
    def df(x):
        return n * x ** (n - 1)
    return newton_find_zero_of_f(f, df)
    
print(nth_root_of_a(1, 64))  #64.0
print(nth_root_of_a(2, 64))  #8.0
print(nth_root_of_a(3, 64))  #4.0
print(nth_root_of_a(4, 64))  #2.82842712474619
print(nth_root_of_a(5, 64))  #2.29739670999407
print(nth_root_of_a(6, 64))  #2.0

這里是在線可運(yùn)行版本 顷扩，把每一步都抽象出來，哈哈哈哈哈慰毅，這看上去很帥隘截，雖然實(shí)際上沒啥必要，最后也可以再curry化。這是個極端的例子婶芭，但是根據(jù)業(yè)務(wù)的需要东臀，我們可以把需要復(fù)用的步驟抽象出來，而抽象的程度也要根據(jù)業(yè)務(wù)開展犀农。

4惰赋、牛逼哄哄的invsqrt求平方根倒數(shù)

吹水時間又到了，究竟這個牛頓迭代法有個毛線用呵哨。那你就要google下invsqrt這個改變游戲史進(jìn)程的函數(shù)赁濒，可以看下知乎上的這個回答 有哪些算法驚艷到了你？孟害。這個函數(shù)到底是干嘛的拒炎，就是求某個數(shù)的2次方根的倒數(shù)。因?yàn)檫@個函數(shù)在3d引擎的中經(jīng)常使用挨务，所以提高這個函數(shù)的效率是至關(guān)重要的击你，這個函數(shù)比系統(tǒng)的1/sqrt(x)快了一半。雖然是倒數(shù)谎柄，但是也是能用牛頓迭代法的啊丁侄。為毛啊，因?yàn)檫@也是曲線啊朝巫，過程跟1中的圖相識鸿摇。

但是！Ｅ场；琛！請看代碼

float InvSqrt (float x){
    float xhalf = 0.5f*x;
    int i = *(int*)&x;
    i = 0x5f3759df - (i>>1);
    x = *(float*)&i;
    x = x*(1.5f - xhalf*x*x);
    return x;
}

第六行就是利用牛頓迭代法的值x - f(x)/f'(x)演算后的結(jié)果糙臼，在這里f(x) = 1/x2，則f'(x)=2/x3恩商。
第三行是把一個float指針強(qiáng)轉(zhuǎn)成int指針变逃，從而進(jìn)行右移運(yùn)算，因?yàn)閒loat類型不能進(jìn)行位運(yùn)算怠堪，從而得到x值指數(shù)部分/2的結(jié)果揽乱，額，這里我不打算展開說啊粟矿，因?yàn)槲艺f不清楚啊啊啊啊凰棉，我引用一段別人的解釋，詳情請移步這里：

所以陌粹，32位的浮點(diǎn)數(shù)用十進(jìn)制實(shí)數(shù)表示就是：M2^{E撒犀。開根然后倒數(shù)就是：M}(-1/2)2^(-E/2)。現(xiàn)在就十分清晰了。語句i> >1其工作就是將指數(shù)除以2或舞，實(shí)現(xiàn)2^{(E/2)的部分荆姆。而前面用一個常數(shù)減去它，目的就是得到M}(1/2)同時反轉(zhuǎn)所有指數(shù)的符號映凳。

至于0x5f3759df胆筒，這就相當(dāng)無解了。wiki的介紹中诈豌，這一句的注釋是這樣的：

    i  = 0x5f3759df - ( i >> 1 );               // what the fuck?（這他媽的是怎么回事仆救？）

關(guān)于這個值后來也有人專門去推理過，還寫成論文矫渔，但是彤蔽，我肯定是看不懂的！０稣丁铆惑！
因?yàn)檫@一步，取到一個非常合適的初始值送膳，只需要一次牛頓迭代就能知道目標(biāo)值啦员魏，一次！５Ｋ貉帧！

最后

沒有最后碌补，今天吹水結(jié)束Ｂ彩！Ｏ谜隆Ｕ蛟取！Ｍ嗫小：骨帧！