旭穹の陋室

扩展欧几里得算法 费马定理（Fermat’s theorem）使我们能够通过二进制指数法（binary exponentiation）在O(log⁡n)O(\log n)O(logn)的操作内计算模数除法的乘法逆元（modular multiplicative inverses），但这只适用于对质数取模。它的一个推广是欧拉定理（Euler’s theorem），此定理表明，如果 m 和 a 互质（即他们之间除了1，没有其他公约数），那么 aϕ(m)≡1(mod  m)a^{\phi(m)} ≡ 1(\mod m) aϕ(m)≡1(modm) 其中，ϕ(m)\phi(m)ϕ(m)是欧拉函数，定义为小于mmm的正整数xxx与mmm互质的个数。当mmm是一个质数时，所有小于mmm的m−1m-1m−1个正整数与mmm都是互质的，所以ϕ(m)=m−1\phi(m)=m-1ϕ(m)=m−1，这就得出了费马定理。 这允许我们通过ϕ(m)\phi (m)ϕ(m)来计算aaa的模除乘法逆元，记作aϕ(m)−1a^{\phi (m)−1}aϕ(m)−1，但反过来，计算ϕ(m)\phi (m)ϕ(m)并不那 ...

二进制指数法 在模数运算和计算代数中，你经常需要把一个数增大到它的nnn次方。这在进行模数除法，进行素数测试，或计算一些组合值时非常常见。而且你通常希望在Θ(n)次操作内完成计算。 二进制指数计算，也被称为通过平方进行指数运算，是一种用O(log⁡n)O(\log n)O(logn)次乘法计算nnn次幂的方法。这种方法基于： a2k=(ak)2a2k+1=(ak)2⋅a\begin{align} a^{2k} &= (a^k)^2 \\ a^{2k+1} &= (a^k)^2 \cdot a \end{align} a2ka2k+1​=(ak)2=(ak)2⋅a​​ 为了计算ana^nan，我们可以递归地计算a⌊n/2⌋a^{⌊n/2⌋}a⌊n/2⌋，然后对其平方，如果nnn为奇数，则再乘上aaa，对应如下递归计算： an=f(a,n)={1,n=0f(a,n2)2,2 ∣ nf(a,n−1)⋅a,2∤na^n = f(a,n) = \begin{cases} 1,& n = 0\\ f(a,\frac{n}{2})^2, & 2\space |\spa ...

模运算 计算机通常以从1970年1月1日（即“Unix era” 的开始）至今所经过的秒数来存储时间，并在所有和时间有关的计算中使用这些时间戳。 而人类则以历史的某一时刻为参照，记录时间，这个时刻通常具有政治或宗教意义。例如，写这段话的时候大约从公元1年（这是六世纪东罗马僧侣对耶稣基督出生日期的最佳估计）至今已经过了大约63882260594秒。 相比计算机要处理数位可能达到11位的整个时间戳，人类情境中需要的时间信息通常更为简单和局限。例如，人们可能只需要知道现在是下午2点，该去吃饭了；或者今天是周四，所以赛百味的每日特价三明治是意大利BMT。人们通常只从整个时间戳中提取出自己所需要的信息，而不需要了解全部内容。这种处理方式（只处理1-2位数字）比处理整个时间戳（11位数字）来得更简便。 问题是：今天是周四，一年后的今天是周几？ 首先，我们将一周的每一天赋以一个数值，星期一是0，星期二是1，以此类推，直到星期日是6。在这个系统下，星期四是3。然后，我们添加365（也就是一年的天数）到这个数字上，然后再对7取余数。更方便的算法是，因为365除以7的余数是1，所以一年后的这一天是星期五（ ...

数值理论 1940年，英国数学家 G. H. Hardy 发表了一篇名为《一个数学家的辩解》的著名论文。在这篇论文中，Hardy展开了对“应该为了数学本身而追求数学，而不是为了它的应用”这一观念的讨论。 与数学类似，计算机科学的各个领域也形成了一个谱状结构，在这个“光谱”的一端是数学逻辑和计算性理论，而另一端是网页编程和应用开发。 作者推测阅读者更倾向于应用方面：这本书就是为了证明应当有更多的人从事于实践性的算法设计而不是理论计算机科学的研究，而且，既然你读到了第七章，那么你可能也赞同这个观点。 Hardy 作为一个62岁的数学家，在他撰写本文时，他见证了第一次世界大战和正在进行的第二次世界大战所造成的破坏。而科学的武器化对这些破坏起到了极大的放大作用。作为一个研究数论的人，哈代觉得在一个“无用”的领域工作可以让他平静，并不需要面对任何道德困境。他写道： 到目前为止，还没有人发现数论或相对论可以为战争目的服务，而且在未来的很多年里，人们似乎也不太可能找到这样的用途。 讽刺的是，仅仅5年后，这个声明就被彻底推翻了。因为如果没有对相对论的理解，原子弹就不可能造出来。而计算机时代的密码学 ...

整数除法 与其他算术运算相比，除法在x86等普通计算机上的执行效率非常差。在硬件层面上实现整数和浮点数除法是非常困难的。这涉及到算术逻辑单元（ALU）的设计，除法运算电路会占据大量的空间，同时计算过程也需要很多步骤，因此导致完成div除法运算通常需要10-20个计算周期。对于处理小数据类型的除法运算，其延迟稍微少一些。 x86上的除法和取模 没有人愿意为了独立的模数运算复制所有的这些乱七八糟的东西，所以div指令既用于除法运算也用于取模运算。要执行32位的整数除法运算，你需要特定地将被除数（dividend）放在eax寄存器中，并调用div指令，除数（divisor）作为唯一操作数。这样操作后，商数（quotient）将被存储在eax寄存器中，余数（remainder）将被存储在edx寄存器中。 需要注意的是，被除数实际上需要被存储在eax和edx两个寄存器中：这种机制使得64位除以32位或者甚至128位除以64位的除法运算成为可能，这类似于128位乘法的工作方式。当执行常规的32位除以32位有符号的除法运算时，我们需要将eax进行符号扩展到64位，并存储其高位部分在edx中。 div ...

整数 如果你从一开始就按顺序阅读本章，你可能会想：为什么我要在介绍完浮点数之后才介绍整数运算？整数不是应该更容易吗？ 没错：整数表示更简单。但是，与直觉相反的是，这种简单性反而使得可以有更多的操作方式。相比之下，浮点数的表示往往需要硬件的支持，而有效地操作整数则需要更多创新性的指令集使用方式。 二进制格式 无符号整数就是写作二进制的自然数： 510=1012=4+14210=1010102=32+8+225610=1000000002=28\begin{align} 5_{10} &= 101_2 = 4+1 \\ 42_{10} &= 101010_2 = 32 + 8 + 2 \\ 256_{10} &= 100000000_2 = 2^8 \end{align} 510​4210​25610​​=1012​=4+1=1010102​=32+8+2=1000000002​=28​​ 当运算的结果超出当前整数类型尺寸所能表示的范围（如，对32位无符号整型而言，大于或等于2322^{32}232），就会溢出（overflow），只留下结果的低32位。同样的，如果 ...

快速平方根倒数 浮点数的平方根倒数1x\frac{1}{\sqrt{x}}x​1​在计算规范化向量时会被用到，规范化向量的运算在众多模拟场景中应用广泛，比如在计算机图形学中（例如，用来确定入射角和反射角以模拟光照效果）。 v^=v⃗vx2+vy2+vz2\hat{v} = \frac{\vec{v}}{\sqrt{v_x^2+v_y^2+v_z^2}} v^=vx2​+vy2​+vz2​​v​ 直接计算一个数的平方根倒数的操作非常慢，因为这需要先计算出这个数的平方根，然后再用1除以这个平方根。尽管这两步操作都是由硬件来完成的，但因为它们本身的计算过程就非常慢，所以总的计算过程也就慢。 有一个利用了浮点数在内存中存储的方式的近似算法，效果出乎意料的好。实际上，这种算法已经被实现在硬件中，因此它对于软件工程师来说并不再那么重要，但我们仍然准备详细介绍这个算法，因为它本身具有令人欣赏的优雅性和极高的教育价值。 除了方法本身，该方法的发明历史也很有趣。快速平方根倒数的计算方法被广泛的归功于游戏工作室"id Software"，该工作室在其1999年的标志性游戏"Q ...

牛顿法 在实际的算法中，很少需要达到最大可能的精度。在现实世界的数据中，建模和测量误差通常比浮点数舍入导致的误差大了几个数量级，我们通常很愿意选择一种近似方法，以牺牲精度为代价换取速度。 接下来，我们将介绍这种近似的数值算法中最重要的组成部分之一：牛顿法。 Newton’s Method 牛顿方法是一种简单但非常强大的算法，用于寻找实值函数的近似根，即以下通用方程的解： f(x)=0f(x) = 0 f(x)=0 对于函数fff，我们只假设它至少有一个根，并且在搜索区间上f(x)f(x)f(x)是连续且可微的。虽然也有一些令人厌烦的异常情况，但在实际应用中几乎从未出现过这种情况，所以我们可以非正式地说这个函数是"好”的。 牛顿法的主要思想是以某个初始点x0x_0x0​开始，然后不断通过在函数图像上各点处切线与x轴交点的x坐标作为新的近似点来迭代优化。这个算法的直观理解是，如果f是“良好”的（光滑且单调），并且当前的近似点xix_ixi​已经足够接近根了，那么下一个近似点xi+1x_{i+1}xi+1​会更接近这个根。 为了获取交点xnx_nxn​，我们需要另其切线函数等于零 ...

舍入误差 硬件浮点数的四舍五入方式异常简单：只有在操作的结果不能精确表示时，四舍五入才会发生，而且默认情况下，会采取四舍五入到最接近的可表示数字（在平局的情况下，优先选择以零结尾的数字）。 考虑以下代码片段： float x = 0;for (int i = 0; i < (1 << 25); i++) x++;printf("%f\n", x); 结果输出1677216=2241677216 = 2^{24}1677216=224，而不是数学上期望的结果225=335544322^{25} = 33554432225=33554432。为什么呢？ 当重复递增一个浮点数xxx时，我们最终会到达这样一个状态，浮点数变得如此之大以至于(x+1)(x+1)(x+1)会被舍入回xxx。递增得到的第一个这样的数是2242^{24}224（即尾数位数加1），因为2242^{24}224作为浮点数递增（尾数末位加1）的下一个值是： (float)224+1=224⋅1.0...01⏟×23=1677218(float)2^{24} + 1 = 2^{24} ...

IEEE-754 浮点数 当我们设计自定义的浮点类型时，省略了很多重要的小细节： 在整个数据中我们应该分配多少比特用于表示尾数部分和指数部分？ 符号位为0表示的是正数还是负数？ 这些比特是如何在内存中存储的？ 我们怎么表示0？ 舍入操作是怎么进行的？ 如果我们试图除以0会发生什么？ 如果我们试图对一个负数进行开方会发生什么？ 如果对最大可表示的数继续增大会发生什么？ 我们可以检测到上述行为的发生吗？ 早期的计算机并不支持浮点运算，当供应商开始添加浮点协处理器时，它们对这些问题应该如何解答有些不同的看法。由于实现方式多样，使得浮点运算难以可靠和便携地使用，尤其对于那些开发编译器的人来说。 1985年，电气和电子工程师协会（IEEE）发布了一项名为IEEE 754的标准，这项标准为浮点数应该如何工作提供了正式的规范，这一规范很快就被供应商接受并在几乎所有的通用计算机中使用。 浮点格式 与我们自己实现的浮点类型类似，硬件浮点类型也使用一个位来表示符号（sign），同时使用可变数量的位来表示指数（exponent）和尾数（mantissa）部分。例如，在标准的32位浮点数编码中，它使用第一 ...