乐趣区

关于算法:面试题精选求根号2简单高级算法你肯定不会

前两天逛 github 看到一道很简略的面试题——如何不必库函数疾速求出 $\sqrt2$ 的值,准确到小数点后 10 位!第一反馈这不很简略嘛,大学数据结构课讲二分查找的时候老师还用这个做过示例。但转念一想,能作为大厂的面试题,背地相对没有那么简略,于是我 google 了下,后果找到了更奇妙的数学方法,甚至发现了一件奇闻趣事…… 一道简简单单的面试题,不仅能考查到候选人的编程能力,还能间接考查到候选人的数学素养,难怪很多大厂都会问这个。。。

回到正题,求 $\sqrt2$ 到底有多少种解法,咱们由简入难一步步来看下咱们是如何让计算机更快计算 sqrt 的。

二分查找

首先是略微具备点编程和数据结构根底的人都能想到的二分查找,这里我不再具体解说思路,但还是要编码测试下,次要是测试须要迭代多少次能力达到小数点后 10 位的精度。

public class BSearch {
    static double sqrt2 = 1.4142135624;
    static double delta = 0.0000000001;
    public static void main(String[] args) {
        double l = 1.0;
        double r = 2.0;
        int cnt = 0;
        while (l < r) {double mid = (l + r) / 2;
            if (Math.abs(l - sqrt2) < delta) {break;}
            if (mid * mid > 2.0) {r = mid;} else {l = mid;}
            cnt++;
        }
        System.out.println("通过" + cnt + "次迭代后得" + l);
    }
}
通过 34 次迭代后得 1.4142135623260401

记住这个迭代次数 34。

牛顿迭代

数学学得好的同学必定晓得牛顿迭代法,它是求解线性方程近似解的办法,因为有些方程无奈疾速求出准确解,只能尽可能去迫近。

回到咱们求解 $\sqrt2$ 上,咱们能够结构出多项式方程 f(x),使得 $f(x)=0$ 的一个正根是 $\sqrt2$,最简略的就是 $f(x) = x^2 – 2= 0$,而后咱们就能够使用牛顿迭代求解它的根了。

$$
x_{1}=x_{0}-\frac{f\left(x_{0}\right)}{f^{\prime}\left(x_{0}\right)} = x_0 – \frac{x_0^2-2}{2x_0}
$$

为了更直观些,咱们图例展现下迭代两次的过程。

上图中彩色的曲线是 $f(x)=x^2-2$,咱们最终想要的是它和 x 轴的交点 X,也就是 $\sqrt2$ 的具体值。A 点的坐标是(2,2),绿色的线是 $f(x)=x^2-2$ 在点 A 处的切线,洋红色线是过 A 点的垂线。咱们选的牛顿迭代初始点就是 B,设 B 点的横坐标是 $x_0$,按求导的定义,AC 点的斜率就是 $f^{\prime}(x)$,AB 的长度就是 $f(x)$,晓得了 AB 的长度,AC 的斜率,BC 的长度也就很好求了,就是 $\frac{f\left(x_{0}\right)}{f^{\prime}\left(x_{0}\right)}$,所以 C 点的横坐标就是 $x_1 = x_0 – \frac{f\left(x_{0}\right)}{f^{\prime}\left(x_{0}\right)}$。

怎么样,迭代一次之后就很迫近咱们真正想要的 X 了,当初咱们在 C 点把同样的事件再做一遍,如下图。

再一次迭代后失去了点 E,点 E 比点 C 更加迫近 X 点,我把上图中点 E 部分放大如下。

能够看到点 E 曾经相当靠近 $\sqrt2$ 了,这只是两次迭代的成果,写个代码来看下到底须要迭代多少次就能达到准确到小数点后 10 位的精度。

具体代码如下,这里我取 x 的初始值为 2.0 因为 $\sqrt2$ 不可能大于 2,咱们晓得这点就能够取个近似值,缩小迭代次数。

public class NewtonMethod {
    static double sqrt2 = 1.4142135624;
    static double delta = 0.0000000001;
    public static void main(String[] args) {
        double x = 2.0;
        int cnt = 0;
        while (Math.abs(x - sqrt2) > delta){x = x - (x*x - 2)/(2*x);
            cnt++;
        }
        System.out.println("通过" + cnt + "次迭代后得" + x);
    }
}
通过 4 次迭代后得 1.4142135623746899  

只须要 4 次迭代就能达到二分 34 次迭代的成果了,的确显著快多了。这里再补充一点,实际上 x 的初始值能够取任意负数,然而会影响到性能,我尝试取 1 亿,最终须要 30 次迭代,不过还是比二分快。

为了主观比照下牛顿迭代和二分的性能差别,这里我还是用 JMH 压测下,后果如下,压测后果仅供参考。

Benchmark           Mode  Cnt         Score   Error  Units
Test.NewtonMethod  thrpt    2  96292165.813          ops/s
Test.bSearch       thrpt    2  11856462.059          ops/s

到这里文章进度条只有一半,你是不是感觉我上面要讲比牛顿迭代更好更快的办法?理论我目前没有找到比牛顿迭代又好又快的算法了,然而我找到一个相干的故事,以及它引出的以就义精度换取速度求 $\frac{1}{\sqrt{x}}$ 的神奇算法,当然它也能够用来求 $\sqrt2$。

神奇的数字 0x5f3759df

这首先要从一个诡异的常数说起—— __0x5f3759df__,提到这个常数还得提到一个游戏,Quake-III Arena (雷神之锤 3)。

这是 90 年代的经典游戏之一。该系列的游戏岂但画面和内容不错,而且即便计算机配置低,也能极其流畅地运行。这款游戏起初在 GPL 协定下开源,github 地址 https://github.com/id-Software/Quake-III-Arena,大家从中发现其中一个可能疾速求解 $\frac{1}{\sqrt{x}}$ 的函数,代码如下。

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;
}

看完代码的我

其实这个算法就是牛顿迭代单次的近似解法,具体证实请看卡马克疾速平方根倒数算法,它能以几十倍的速度劣势秒杀其余算法,要晓得几十年前的 CPU 速度可远不迭当初的,速度就是绝对优势。这段代码看起来神奇,它的起源也很离奇。

开始大家都认为这个算法是游戏的开发者 Carmack 发现的,但起初考察发现,该算法在这之前就在计算机图形学的硬件与软件畛域中有所利用,如 SGI 和 3dfx 就曾在产品中利用此算法,所以至今都无人知晓这个算法是谁创造的。

但传奇并没有完结。Lomont 计算出后果当前十分称心,于是拿本人计算出的起始值和卡马克的神秘数字做较量,看看谁的数字可能更快更准确求得平方根。后果是卡马克赢了 … 谁也不晓得卡马克是怎么找到这个数字的。最初 Lomont 怒了,采纳暴力办法一个数字一个数字试过去,终于找到一个比卡马克数字要好上那么一丁点的数字,尽管实际上这两个数字所产生的后果十分近似,这个暴力得出的数字是 0x5f375a86。Lomont 为此写下一篇论文 Fast Inverse Square Root。From 百度百科

来看看这个算法的理论运行成果怎么样吧,上面是我批改后的 java 代码,下面的 C 代码中有些操作 java 中并不反对,所以须要做些改变。

public class CarmackMethod {public static void main(String[] args) {
        float x = 2.0f;
        float xhalf = 0.5f*x;
//        int i = *(int*)&x;
        int i= Float.floatToRawIntBits(x);
        i = 0x5f3759df - (i >> 1);
//        x = *(float*)&i;
        x = Float.intBitsToFloat(i);
        x = x*(1.5f - xhalf*x*x);
        System.out.println(1.0/x);
    }
}
1.4145671304934657

当然这里精度就不行了,只准确到了 0.001。

这里我把下面 3 种算法的精度都升高到 0.001 而后用 JMH 做个不严格的性能测试。这里说不严格是因为我只做 $\sqrt2$ 的测试,而且用的是 java 实现的,而且像是 CarmackMethod 的实现,可能因为 java 和 c 的运行机制的不同,性能会受很大影响,上面这个后果 __仅供娱乐__,看看就好。

Benchmark            Mode  Cnt          Score   Error  Units
Test.CarmackMethod  thrpt    2  111286742.330          ops/s
Test.bSearch        thrpt    2   58705907.393          ops/s
Test.NewtonMethod   thrpt    2  110232331.339          ops/s

各种编程语言是如何实现 sqrt?

下面说了 3 个解法,那你是不是也很好奇目前各种编程语言库函数里对 sqrt 是如何实现的?我也很好奇,于是咱们帮你们翻了 jdk 源码,发现它基本就不须要本人实现 sqrt,毕竟 sqrt 在计算机领域是有个比拟罕用的计算,所以支流的 CPU 架构都提供了对 sqrt 的硬件反对,只须要一条汇编指令就能够了,在 x86 架构下 sqrt 能够间接用上面这条汇编指令。

  sqrtsd %1, %0":"=x"(res) :"xm" (x)

在 Risc- v 中能够能够用 fsqrt.sfsqrt.d指令,Rics- v 中文手册

硬件能够在一个指令周期内实现一个数的开方,相比那些须要几十甚至成千盈百个 CPU 指令实现的各种算法而言,这速度差别不言而喻。实际上上,古代 CPU 在多外围、流水线、多发射、超标量……等技术的加持下,一般家用 CPU 都能够做到每秒百亿次的浮点运算。

__生存处处有惊喜__,当我关上 python math 模块的源码时,没有发现浮点数的求根(预计也是间接用的 CPU 级指令),但我发现了一个更骚的对 64 位整数求根的操作,所以这里再补充介绍一个 python 的近似求根算法。

python 中的_approximate_isqrt()

上面这段代码能够返回输出值求根后的整数局部,但齐全不晓得是什么原理。

_approximate_isqrt(uint64_t n)
{uint32_t u = 1U + (n >> 62);
    u = (u << 1) + (n >> 59) / u;
    u = (u << 3) + (n >> 53) / u;
    u = (u << 7) + (n >> 41) / u;
    return (u << 15) + (n >> 17) / u;
}

源码中有正文,这段 C 代码能够翻译为以下 python 代码,不过我仍旧看不懂。

    def isqrt(n):
        """Return the integer part of the square root of the input."""
        n = operator.index(n)
        if n < 0:
            raise ValueError("isqrt() argument must be nonnegative")
        if n == 0:
            return 0
        c = (n.bit_length() - 1) // 2
        a = 1
        d = 0
        for s in reversed(range(c.bit_length())):
            # Loop invariant: (a-1)**2 < (n >> 2*(c - d)) < (a+1)**2
            e = d
            d = c >> s
            a = (a << d - e - 1) + (n >> 2*c - e - d + 1) // a
        return a - (a*a > n)

结语

这篇博客从立题到实现经验了好几天的工夫,期间整顿思路、编码、绘图、查阅材料、批改欠缺总累计耗时近 8h。写作不易,如果文章对你有用欢送素质三连(点赞、珍藏加关注)。

援用链接

[1] 牛顿迭代法: https://baike.baidu.com/item/…
[2] 卡马克疾速平方根倒数算法: http://jcf94.com/2016/01/14/2…
[3] Fast Inverse Square Root: http://read.pudn.com/download…
[4] From 百度百科: https://baike.baidu.com/item/…
[5] 这条汇编指令: https://code.woboq.org/usersp…
[6] Rics- v 中文手册: http://crva.ict.ac.cn/documen… Chinese-v2p1.pdf
[7] python math 模块: https://github.com/python/cpy…

退出移动版