01 Algorithms with Numbers

基本算数

一个整数可以用物体的个数来表示。比如用8个点我们就能表示数字8。但这样的话对于大的数字我们就必须用非常多的点，造成了不方便。因此为了更方便地表示一个数，我们通常要选择一个进制。我们在日常生活中常用的是10进制，在计算机中常用的是2进制，等等。

数$N$在$a$进制下需要$\log_a N$位，在$b$进制下需要$\log_b N$，二者的比值$\log_b^a$是常数。因此我们可以认为，对进制的选择是不影响算法的时间复杂度的。

加法

从我们在小学时做的列竖式计算中就可以看出，逐位计算并进位的复杂度是正比于数的位数的。设数的位数为$n$，则复杂度是$O(n)$。同时，仅仅读入数据就要消耗$O(n)$，因此这就是加法的最优复杂度了。

乘法

按照列竖式的方法，乘法的复杂度是$O(n^2)$。

我们也可以通过这样的递归方法来完成乘法：

$x \cdot y=\left\{\begin{array}{ll} 2(x \cdot\lfloor y / 2\rfloor) & \text { if } y \text { is even } \\ x+2(x \cdot\lfloor y / 2\rfloor) & \text { if } y \text { is odd } \end{array}\right.$

递归次数正比于$y$的位数$n$，每次完成加法的复杂度也是$n$，总的复杂度依然是$O(n^2)$。

乘法本质上是多项式乘法，通过快速傅里叶变换是可以优化到$O(n \log n)$的。我们在讨论“分治”算法时还会来讨论这个问题。

除法

我们想算出$x$除以$y$的商$q$和余数$r$。考虑递归，假设已经算出了$\lfloor x/2\rfloor$除以$y$的余数$q'$和$r'$。则有$\lfloor x/2 \rfloor=q' \cdot y + r'$。两边同时乘以2，得到$2\lfloor x/2 \rfloor=2q' \cdot y + 2r'$。如果$x$是偶数，那么有$x=2q' \cdot y + 2r'$，否则有$x = 2q' \cdot y + 2r'+1$。注意如果$2r'$或$2r'+1$超出了$y$，要减一个$y$。这样就推出了$q,r$。

递归的次数为$n$，减$y$的复杂度为$n$，总复杂度是$O(n^2)$。

模运算

对于$x$除以$N$的商为$q$，余数为$r$，可以写出$x = qN+r$，其中$r$被限制在$\{0,\cdots,N-1\}$内。如果对于$y$也有$y=q'N+r$，那么就称$x$和$y$是同余的。记作

$x \equiv y \pmod N$

由此可见$x-y$一定是$N$的倍数，即$x-y \equiv 0 \pmod N$。

如果$x \equiv x' \pmod N$，同时$y \equiv y' \pmod N$。那么此时有$x-x' = sN$, $y-y' = tN$。那么有$x+y=(s+t)N+x'+y'$，所以$x+y \equiv x'+y' \pmod N$。而$xy = (x'+sN)(y'+tN) = x'y'+(tx'+sy'+stN)N$，所以$xy \equiv x'y' \pmod N$。这告诉我们，在模运算中加法和乘法是可以被同余的项代换的。

模运算中由于每个数都被限制在$\{0,\cdots,N-1\}$内，因此加法不会超过$2N-2$，因此加法以后最多只需要减一次$N$就能保证模运算的范围。两数相乘不会超过$(N-1)^2$，因此需要做一次除法（区域）才能保证范围，但我们注意到数字的位数最多只会扩展一倍。

注意，在模运算下两个都不为0的数相乘可能得到0。比如5*5在模25意义下为0。

快速幂

在模意义下计算$x^y$要用到快速幂算法。核心是：

$x^{y}=\left\{\begin{array}{ll} \left(x^{2}\right)^{\lfloor y / 2\rfloor} & \text { if } y \text { is even } \\ x \cdot\left(x^{2}\right)^{\lfloor y / 2\rfloor} & \text { if } y \text { is odd } \end{array}\right.$

递归的次数是$n$，每次乘法用去$n^2$，复杂度为$O(n^3)$.

欧几里得算法(GCD)

对于$x \geq y$，有$gcd(x,y)=gcd(x\mod y,y)$。一直这样递归下去，直到$x$成为$y$的倍数为止，我们就求出了最大公约数。这称为求解最大公约数的辗转相除法。要证明这个关系，只需证明辗转相减法$gcd(x,y)=gcd(x-y,y)$，因为模运算只是减法运算的一种“捷径”。这是自然的，因为$x,y$都可以看作是由它们的公约数这个基本积木拼成的，它们的差依然是由同样的基本积木组成的。严格地，设$x,y$有一个公约数为$k$，那么有$x = s \cdot k$，$y = t \cdot k$。于是$x-y = (s-t) \cdot k$。因此$x,y$的任意一个公约数一定也是$x-y,y$的一个公约数。所以一定有$gcd(x,y) \leq gcd(x-y,y)$。反之同理，$x-y,y$的任意一个公约数一定也是$x,y$的一个公约数，所以又有$gcd(x-y,y) \leq gcd(x,y)$。联立即得$gcd(x,y)=gcd(x-y,y)$。这样我们就证明了辗转相除法。

为了分析辗转相除法的复杂度，我们来估计一下$x \mod y$的大小。如果$y \leq x/2$，那么$x \mod y < y \leq x/2$；如果$y>x/2$，那么$x \mod y = x-y < x/2$。因此我们每一次至少把一个数缩小了一半，递归的次数是$O(n)$的。而模运算的计算复杂度就是除法的复杂度$O(n^2)$，因此总的复杂度是$O(n^3)$。

扩展欧几里得算法(EXGCD)

对于方程$ax+by=d$，左边可以提取出$a,b$的最大公约数，因此$d$一定是$gcd(a,b)$的倍数。换句话说，$ax+by$构成的集合一定是“量子化”的，$gcd(a,b)$是最小的单位。如果$d$不是$gcd(a,b)$的倍数，那么方程一定是无解的。

而我们要说明，$ax+by=gcd(a,b)$一定是有解的。由欧几里得算法得，$gcd(a,b)=gcd(b,a \mod b)$。我们列一个新的方程$bx'+(a\mod b)y'=gcd(b,a \mod b)$。联立可得$ax+by=bx'+(a\mod b)y'$。注意到$a \mod b = a-b\lfloor a/b \rfloor$。因此有$ax+by=ay'+b(x'-\lfloor a/b \rfloor y')$。假设我们已经递归地解得了$x',y'$，就可知$x=y',y=x'-\lfloor a/b\rfloor y'$是一组可行的解。而由于$gcd(a,b)$递归的尽头$b=0$，此时$x$一定有解，因此返回去最初的$x,y$一定有解。这样求解方程的方法就成为扩展欧几里得算法。

我们通过Exgcd已经求出了$(x_0,y_0)$是$ax+by=gcd(a,b)$的一组解，如何用它去找到方程的所有可行解呢？如果$a(x_0+p)+b(y_0+q)=gcd(a,b)$也是一组解，那么必须满足$ap+bq=0$（很像线性代数中的零空间的解）。反过来，所有满足$ap+bq=0$的$p,q$一定能够保证$(x+p,y+q)$还是方程的解。因此我们只需解出$ap+bq=0$的所有可行解$(p,q)$即可。设$a,b$除掉彼此的最大公约数得到$a',b'$，$ap+bq=0$当且仅当$a'p+b'q=0$。要使得$a'p$是$b'$的倍数，$a'p$必须是$lcm(a',b')$的倍数，而$gcd(a',b')=1$，所以$lcm(a',b')=a'b'$，$a'p$必须是$a'b'$的倍数，即$p$必须是$b'$的倍数，因此必须满足$p=\dfrac{kb}{gcd(a,b)}$。而容易验证所有这样的$p$都是满足的。综上，$ax+by=gcd(a,b)$的全部解就是$(x_0+\dfrac{kb}{gcd(a,b)},y_0-\dfrac{ka}{gcd(a,b)})$。

对于一般的$ax+by=d$，我们可以解出$ax+by=gcd(a,b)$，再把得到的解乘上$\dfrac{d}{gcd(a,b)}$倍，得到一组可行的解$(x_0,y_0)$。而后，完全相同的，我们依然利用$ap+bq=0$把特解拓展到通解，这个过程中只用到了关于$a,b$的性质，与等式右侧是无关的，因为我们所作的修正是完全相同的。得到的通解依然是$(x_0+\dfrac{kb}{gcd(a,b)},y_0-\dfrac{ka}{gcd(a,b)})$。

乘法逆元

我们想要定义模运算下的除法，本质上就是定义模运算下的倒数。如果$ax \equiv 1 \pmod N$，就称$x$是$a$在模$N$意义下的乘法逆元（其中$a,x$都在$0$到$N-1$之间）。

并不是所有$a$都有乘法逆元的。比如$a=0$就找不到这样的$x$；如果$a=2,N=6$，也可以验证没有任何的$x$是满足的。由扩展欧几里得，我们发现乘法逆元存在当且仅当$gcd(a,N)=1$：如果$ax \equiv 1\pmod N$有解，则$ax+Ny=1$有解，这要求$gcd(a,N)=1$。而只要$gcd(a,N)=1$满足，由EXGCD就可以给出$ax+Ny=1$的一组解。这样得到的解可能会使得$x$没有落在$0$到$N-1$内，但根据通解的表达式，$x$只能间隔$\dfrac{N}{gcd(a,N)}=N$地变化。所以我们把$x$加上若干个$N$使它落在$0$到$N-1$就可以了，这一定是能够做到的。

对于给定的$a$，逆元$x$是唯一的。因为如果存在$x' \neq x$使得$x'$也是$a$的乘法逆元，那么就有$ax \equiv ax' \equiv 1 \pmod N$。于是$a(x-x')=tN$。由于$gcd(a,N)=1$，所以$a(x-x')$一定是$a,N$的最小公倍数的倍数，即$a(x-x')=kaN$，于是$x-x'=kN$，这与$0 \leq x,x'<N$矛盾。

素数

费马小定理

费马小定理指出：如果$p$是素数，那么对任意满足$1 \leq a < p$的整数$a$，一定有$a^{p-1} \equiv 1 \pmod p$

我们发现， 对于$i=1,\cdots,p-1$，$a \cdot i$一定是两两不同且不为0的（刚好构成了一个permutation）：首先，如果$a \cdot i \equiv 0\pmod p$，根据$gcd(a,p)=1$，$a$存在乘法逆元，两边同时乘上$a$的乘法逆元，就会得到$i \equiv 0$，矛盾；如果对于$i \neq j$成立$a \cdot i \equiv a \cdot j$，那么同样乘以乘法逆元会得到$i \equiv j$，矛盾。因此，$a \cdot i$刚好取遍了$1,\cdots,p-1$。把他们全都乘起来，就得到$a^{p-1}!(p-1)! \equiv (p-1)! \pmod p$。而$1$到$p-1$的每个数都没有$p$这个质因子，因此$(p-1)!$一定是与$p$互质的，因此也有乘法逆元。于是两边同时乘以这个逆元，得到$a^{p-1} \equiv 1 \pmod p$。这样就证明了费马小定理。

遗憾的是，费马小定理的逆命题并不成立——我们不能由$a^{N-1} \equiv 1 \pmod N$推出$N$是素数。反例：合数$341=11\times 31$，却有$2^{340} \equiv 1 \pmod {341}$。甚至存在这样的合数$N$，使得所有比$N$小的与它互质的数$a$全都满足$a^{N-1} \equiv 1 \pmod N$。（$N=561=3\times 11 \times 17$）这样的数被称为Carmichael数。

$a^{N-1} \equiv 1 \pmod N$虽然无法推出$N$是素数，但可以推出$gcd(a,N)=1$。如果$gcd(a,N)>1$，那么$a^{N-1}-tN$一定是$gcd(a,N)$的倍数，因此永远不可能等于1。

根据费马小定理，我们又得到了另一种求乘法逆元的方式：$p$是素数时，$a$在模$p$下的逆就是$a^{p-2} \mod p$。

素数判断

我们可以证明，如果$N$不是Carmichael数，即存在$a$（与$N$互质）使得$a^{N-1} \not\equiv 1 \pmod N$，那么这样的$a$在$1$到$N-1$间至少有一半。我们先任意固定这样的一个$a$，它是一个常数，满足$N-1$方不为1。假如我们能找到一个$b$满足$b^{N-1} \equiv 1 \pmod{N}$，那么对于$ab$（模意义下）有$(ab)^{N-1} \equiv a^{N-1}b^{N-1} \equiv a^{N-1} \not\equiv 1 \pmod {N}$。因为$a$与$N$互质，所以有乘法逆元，所以对于某个$b'$满足$b \not\equiv b' \pmod{N}$，一定成立$ab \not\equiv ab' \pmod N$。因此对于每个这样的$b$，我们一定可以对应地找到一个$ab$。前者满足$N-1$次方等于1的条件，后者满足不等于$N-1$次方等于1的条件。由于$b^{N-1} \equiv 1 \pmod{N}$成立，$b$一定存在乘法逆元，而显然$a \not \equiv 1 \pmod N$，所以$b \not\equiv ab \pmod N$一定成立。我们一对对选这样的$b,ab$，会形成两个独立的集合，一个全都是$N-1$次方等于1的，一个全都是$N-1$次方不等于1的。也就意味着，每找到一个等于1的就一定能找到一个不等于1的，所以不等于1的至少占了一半。

所以我们可以设计这样一个素数判定算法，对任意给定的$N$，我们在$2$到$N-1$间随机一个$a$，计算$a^{N-1} \mod N$。如果答案为1就返回yes否则返回no。如果$N$是素数，那么程序一定正确；如果$N$不是素数且$N$不是Carmichael数，此时返回yes的概率小于$1/2$。也就是说此时程序出错的概率小于$1/2$。如果我们连着做$k$次，出错的概率就被降到了$1/2^k$。

由于Carmichael数是罕见的，我们直接忽略这种情况。这样我们就有了一个高效的算法，复杂度就是计算快速幂的复杂度$O(n^3)$。考虑到$k$后的复杂度为$O(kn^3)$。

素数生成

Lagrange指出$N$以内的素数个数同阶于$\dfrac{N}{\ln N}$。这说明素数是非常密集的。我们可以不停地随机一个$N$以内的数，然后用上面的算法判定它是否是一个素数，直到它通过判定为止。这样我们就大概率能得到一个素数了。

每个随机到的数是素数的概率为$\dfrac{1}{\ln N}$，因此期望$\ln N$次随机后算法结束，共$O(n)$次。每次判定的复杂度是$O(n^3)$。因此总复杂度为$O(n^4)$。

密码学

一切发送出去的（通信、电波……）信息都是可能被监听的。想要达成加密的目的就必须有从未发送出去过的“密钥”。如果Alice和Bob都有一个二进制串$r$（密钥），那么Alice可以先把要发送的串与$r$异或，Bob收到后再与$r$异或一次就能得到明文了。但这要求$r$很长，而且两人要在不被窃听的情况下交换$r$是不容易的。

RSA

RSA是公钥加密的典范。在这里，Bob手上有两个密钥，一个密钥全世界只有他知道，称为私钥；一个密钥向全世界公开，称为公钥。Alice用公钥加密自己的信息，Bob利用私钥解密，这样就永远都不可能有人破解了。

在这里，Bob随机生成两个质数$p,q$。定义$N=pq$，再任意找一个与$(p-1)(q-1)$互质的数$e$（这是容易办到的，比如取$e=3$）。$N,e$构成了我们的公钥。求出$e$在模$(p-1)(q-1)$意义下的乘法逆元$d$，把$d$作为我们的私钥。

由于$ed \equiv 1 \pmod{(p-1)(q-1)}$，可以写出$ed = k(p-1)(q-1)+1$。而$x^{p-1} \equiv 1\pmod p$，因此$x^{ed}\equiv x^{k(p-1)(q-1)+1} \equiv x \cdot (x^{p-1})^{k(q-1)} \equiv x \pmod p$。同理，也有$x^{ed} \equiv x \pmod q$。因此$x^{ed}-x=t_1p=t_2q$。说明$x^{ed}-x$是$p,q$的公倍数，因此有$x^{ed} \equiv x \pmod{pq}$，即$(x^{e})^d \equiv x \pmod N$。这个式子对于任意的$0$到$N-1$的$x$都成立，如果把$f(x)=x^e$看作一个映射，$g(x)=x^d$看作一个映射，那么$g(f(x))=x$。这说明$g(x)$的象集至少有$N$个元素，那么其定义域也至少要有$N$个元素。而$x$最多只有$N$个元素。因此$f,g$都必须是双射。所以对于任意的$x$，我们都可以用$f$加密（它依赖于$e,N$，因此是公钥加密），接收者都可以用$g$来解密（由于涉及到$d$，因此是私钥解密）。

RSA的核心在于，作为公钥外人只知道$N$而不知道$p,q$。在$N$很大时是难以找出$p,q$的。因此也就难以找出$(p-1)(q-1)$，也就难以求出逆元$d$。我们把要发送的信息以$x^e$的方式加密，不知道$d$就无法解密，只有Bob能够解密。