数论知识笔记

发表于 更新于 本文字数： 8k 阅读时长 ≈ 7 分钟
有关数论知识的简要笔记

数论知识笔记

Pell 方程

标准型

Pell 方程是指形如 \[ x^2-Dy^2=1 \] 的二元二次方程，其中 \(D\) 为正整数，主要研究的是满足这个方程的整数解 \((x_i,y_i)\) 和最小整数解 \((x_0,y_0)\)。

通常我们使用连分数来求解其整数解，将 \(\sqrt{D}\) 写成连分数的形式，即 \[ \sqrt{D}=[a_0;a_1,a_2,\cdots,a_{n-1},a_n] \] 如果 Pell 方程有解，总会找到 \(a_n\) 满足 \(a_i=2a_0\)，那么连分数 \[ [a_0;a_1,a_2,\cdots,a_{n-1}] \] 可能是 Pell 方程的一组解 \[ \frac{x}{y}=[a_0;a_1,a_2,\cdots,a_{n-1}] \] > 实际上对于方程 \(x^2-Dy^2=C\) 可以证明解 \((x,y)\) 一定是 \(\sqrt{D}\) 的一个有理近似。

Sagemath 求解代码如下

1
2
3
4
5
6
7
8
9
10
D = 1117
interval = 200
cf = continued_fraction(sqrt(D))
a0 = cf[0]
for i in range(interval):
    if a0 * 2 == cf[i]:
        c = cf.convergent(i - 1)
        x, y = c.as_integer_ratio()
        if x**2 - D * y**2 == 1:
            print((x, y))

一个更特殊的求解方式是，不妨设佩尔方程的解为 \[ (x_1,y_1),(x_2,y_2),\cdots,(x_d,y_d) \] 那么求解第 \(i+j\) 个解的公式为 \[ x_{i+j}=x_ix_j+Dy_iy_j\\ y_{i+j}=x_iy_j+x_jy_i \] 例如当 \(D=1117\) 时的解已知 \((x_1,y_1)\) 和 \((x_2,y_2)\)，那么求解得到 \((x_3,y_3)\) 为

1
2
3
4
5
6
7
8
D = 1117
x1, y1 = (87897747594260774254246835664214545379849, 2629972211566463612149241455626172190220)
x2, y2 = (15452028064288751455987343781568665778822611449292997550027848699978973103390525601, 462337267264377629956081365965549206173112084990468349927211645593553260365753560)
def addition(p1, p2):
    x1,y1 = p1
    x2,y2 = p2
    return ((x1*x2+D*y1*y2), (x1*y2+x2*y1))
addition((x1, y1), (x2, y2))

理论上我们仅需要知道一个解，我们就可以得到无穷多 Pell 方程的解。

推导过程为，考虑两组解满足 \[ x_1^2-Dy_1^2=1\newline x_2^2-Dy_2^2=1 \] 两式相乘结果为 \[ (x_1^2-Dy_1^2)(x_2^2-Dy_2^2)=1 \] 化简可得 \[ (x_1x_2+Dy_1y_2)^2-D(x_1y_2+x_2y_1)^2=1 \] 故显然 \((x_1x_2+Dy_1y_2,x_1y_2+x_2y_1)\) 也是佩尔方程的整数解。

非标准型

Pell 方程是指形如 \[ x^2-Dy^2=C \] 的二元二次方程，其中 \(D\) 为正整数，\(C\) 为任意整数，主要研究的是满足这个方程的整数解 \((x_i,y_i)\)。

考虑标准型方程 \[ r^2-Ds^2=1 \] 我们可以轻松求出一个整数解 \((r,s)\) 满足方程。

我们同样通过连分数求解 \((x_0,y_0)\) 满足 \[ x_0^2-Dy_0^2=C \] 那么两式相乘，得到 \[ (x_0^2-Dy_0^2)(r^2-Ds^2)=C \] 可以化简得到 \[ (x_0r+Dy_0s)^2-D(x_0s+y_0r)=C \] 故显然 \((x_0r+Dy_0s,x_0s+y_0r)\) 也是该非标准型佩尔方程的整数解。

就算是对于方程 \[ Ax^2+Bx+Cy^2+Dy+E=0 \] 我们也可以配方化为 \[ x^2-Dy^2=C \] 来进行求解整数解。

代码

需要用到 Sagemath 的 continued_fraction。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
from sage.all import *

def pell_roots(D: int, C: int = 1):
    """
    求佩尔方程$x^2-Dy^2=C$的整数解(x,y)

    参数：
        D: 佩尔方程参数
        C: 佩尔方程常数，默认为1
    
    返回值：
        迭代器，整数解(x,y)
    """
    intervals = 2**10
    def get_a_root(D: int, C: int):
        """
        获得佩尔方程$x^2-Dy^2=C$的一个整数解(x,y)

        参数：
        D: 佩尔方程参数
        C: 佩尔方程常数，默认为1
    
        返回值：
            整数解(x,y)，如果在intervals限制内求解失败，返回(0,0)
        """
        cf = continued_fraction(sqrt(D))
        for i in range(1, intervals):
            c = cf.convergent(i - 1)
            x, y = c.as_integer_ratio()
            if x**2 - D * y**2 == C:
                return x, y
        for y in range(1, intervals):
            x2 = C + D * y**2
            if x2 <= 0:
                continue
            x = isqrt(x2)
            if x ** 2 == x2:
                return x, y
        return 0, 0
    
    r, s = get_a_root(D, 1)
    x, y = get_a_root(D, C)
    while True:
        yield x, y
        x, y = x * r + D * y * s, r * y + s * x

D = 7
C = 2
x, y = next(pell_roots(D, C))
print(x**2 - D * y**2 == C)

毕达哥拉斯三元组

毕达哥拉斯三元组指的是形如 \[ x^2+y^2=z^2 \] 方程的整数解，即整勾股数。

当 \(x=2mn\) 时，\(y=m^2-n^2,z=m^2+n^2\)。

最简单的情况为，\(x\) 是偶数，那么 \(m=1,n=x/2\) 即可。

奇数情况时，\(y=(a^2-1)/2,z=(a^2+1)/2\)。

卡迈尔数

卡迈尔数是指一系列合数（相对于质数）\(n\) 有得性质：使所有与 \(n\) 互素的整数 \(b\)，满足 \[ b^{n-1}\equiv1\pmod{n} \] 费马小定理说明所有的质数都有这个性质，所以卡迈尔数也被称作伪质数。

同时由于这些数的存在，使得费马质性检验变得不可靠。

一个简单的卡迈尔数为，

1
2
p = 29674495668685510550154174642905332730771991799853043350995075531276838753171770199594238596428121188033664754218345562493168782883
N = p * (313 * (p - 1) + 1) * (353 * (p - 1) + 1)

拓展欧几里得算法

1
2
3
4
5
6
7
8
9
def extendedGCD(a, b):
    # a*xi + b*yi = ri
    if b == 0:
        return 1, 0, a
    else:
        x, y, q = extendedGCD(b, a % b)
        # q = gcd(a, b) = gcd(b, a%b)
        x, y = y, (x - (a // b) * y)
        return x, y, q

可以使用 gmpy2.gcdext(a, b) 替代，返回三元元组 \((g,s,t)\)

其中\(g=\gcd(a, b)=sa+tb\)

存在相关性的定理：裴蜀定理。

裴蜀定理说明了对任意整数 \(a,b\) 和它们的最大公约数 \(\gcd{(a,b)}=d\)，那么对任意整数 \(x,y\) 都有 \(ax+by\) 为 \(d\) 的倍数，即 \[ ax+by=0\pmod{d} \] 特别地，一定存在整数使得 \(ax+by=d\) 成立。

一个重要推论为，\(a,b\) 互质的充分必要条件为存在整数 \(x,y\) 使得 \(ax+by=1\)。

欧拉函数

  • 如果\(p\)是素数且\(k\geq 1\)，则 \[ \phi (p^k)=p^k-p^{k-1} \] 特别地，当\(k=1\)时，表达式即为 \[ \phi (p)=p-1 \]

  • 如果\(gcd(m, n)=1\)，即\(m,n\)互质，则 \[ \phi (mn)=\phi (m)\phi (n) \]

其中与之相关的有欧拉公式：

对任意两个互素整数 \(a,n\)，有 \[ a^{\varphi(n)}=1\pmod{n} \]

中国剩余定理

设\(m\)与\(n\)是整数，\(gcd(m,n)=1\)，\(b\)与\(c\)是任意整数，则同余式组 \[ x\equiv b\ (mod\ m)\ 与\ x\equiv c\ (mod\ n) \] 恰有一个解 \(0\leq x\leq mn\)

详细内容可以在这里了解

威尔逊定理

威尔逊定理给出了判定一个自然数是否为素数的充分必要条件：

当且仅当 \(p\) 为素数时，有 \[ (p-1)!\equiv-1\pmod{p} \] 或推论 \[ (p-1)!+1=0\pmod{p} \] 即 \((p-1)!+1\) 为 \(p\) 的倍数。

一般用于辅助推导，相关的有伽马函数 \[ \Gamma(n)=(n-1)! \]

费马小定理

如果 \(p\) 是一个素数，取任意的非 \(p\) 整数 \(a\)，有 \[ a^{p-1}\equiv1\pmod{p} \] 推论有： \[ na^{p-1}\equiv n\pmod{p},n\in\Z\\ a^{p-2}\equiv a^{-1}\pmod{p}\\ a^b\equiv a^{b\pmod{p-1}}\pmod{p},b\in\Z\\ \]

高次剩余

对于同余式 \(X^n\equiv d\pmod{p}\)

当 \(n|p-1\) 时，满足 \[ d^{\frac{p-1}{n}}\equiv 1\pmod{p} \] 才有解，此时解数为 \(n\)；

当 \(n\) 无法整除 \(p-1\)，满足 \[ d^{\frac{p-1}{k}}\equiv 1\pmod{p} \] 其中 \(k=\gcd(n,p-1)\)，此时解数为 \(k\)。

雅可比符号

首先明确勒让德符号，整数 \(a\) 对整数 \(p\) 的勒让德符号定义如下

  • 如果 \(a\equiv 0\pmod{p}\)，那么 \(\left(\frac{a}{p}\right)=0\)
  • 如果 \(a\ne 0\pmod{p}\)，且存在某个整数 \(x\) 满足 \(x^2\equiv a\pmod{p}\)，那么 \(\left(\frac{a}{p}\right)=1\)
  • 如果 \(a\ne 0\pmod{p}\)，但不存在某个整数 \(x\) 满足 \(x^2\equiv a\pmod{p}\)，那么 \(\left(\frac{a}{p}\right)=-1\)

可以写出公式 \[ \left(\frac{a}{p}\right)\equiv a^{\frac{p-1}{2}}\pmod{p} \] 雅可比符号是勒让德符号的推广，给出扩充定义：

  • 如果 \(m=p_1\cdots p_r\)，那么 \[ \left(\frac{a}{m}\right)=\left(\frac{a}{p_1}\right)\cdots\left(\frac{a}{p_r}\right) \]

  • 如果 \(a=a_1\cdots a_t\)，那么 \[ \left(\frac{a}{m}\right)=\left(\frac{a_1}{m}\right)\cdots\left(\frac{a_t}{m}\right) \]

  • 如果存在 \(\gcd{(b,m)=1}\)，那么有 \[ \left(\frac{ab^2}{m}\right)=\left(\frac{a}{m}\right) \]

显而易见 \[ \left(\frac{b^2}{m}\right)=\left(\frac{b}{m}\right)\left(\frac{b}{m}\right)=(\pm1)^2=1 \]

Hensel Lifting

亨塞尔提升（Hensel Lifting），即亨塞尔引理，规定了素数模幂多项式的根被提升到模更高幂的根的条件。

例如说已知 \(x_1\) 满足如下方程 \[ x_1^e\equiv a\pmod{p} \] 即 \(x_1\) 是 \(a\) 在模 \(p\) 下的 \(e\) 次方根，求解 \[ x_k^e\equiv a\pmod{p^k} \] 其中，\(x_k\) 是 \(a\) 在模 \(p^k\) 下的 \(e\) 次方根。

这里仅考虑首一单项式。

亨塞尔引理告诉我们有如下公式用于求解 \(x_k\) \[ x_{k+1}=x_k-\frac{f(x_k)}{f'(x_1)}\pmod{p^{k+1}} \] 推导如下，

我们考虑在模 \(p^{k+1}\) 下的方程 \[ x_{k+1}^e\equiv a\pmod{p^{k+1}} \] 那么必然在降幂下有同余方程 \[ x_{k+1}^e\equiv a\pmod{p^k} \] 同时存在 \(a\) 在模 \(p\) 下的 \(e\) 次方根为 \(x_1\) 满足 \[ x_k^e\equiv a\pmod{p^k} \] 那么显然，\(x_k\) 与 \(x_{k+1}\) 存在关系为 \[ x_{k+1}=x_{k}+p^kv \] 问题的关键就变为了求解 \(v\)。

运用二项式定理，可以得到 \[ \begin{align} x_{k+1}^e&=(x_k+p^kv)^e\\ &=x_1^e+ex_1^{e-1}p^kv+\cdots+(p^kv)^e\\ &=x_1^e+ex_1^{e-1}p^kv\pmod{p^{k+1}} \end{align} \] 由 \(x_{k+1}^e\equiv a\pmod{p^{k+1}}\) 可得 \[ a-x_k^e\equiv ex_k^{e-1}p^kv\pmod{p^{k+1}} \] 我们知道，\(ex_k^{e-1}\) 不可能整除 \(p\)，必然有 \(\gcd{(ex_1^{e-1}p^k,p^{k+1})}=p^k\)；

同时由于 \(x_k^e\equiv a\pmod{p^k}\)，必然有 \(a-x_k^e\) 整除 \(p^k\)。

所以上式的式子左右两边同时整除 \(p^k\)，那么化简得到 \[ \frac{a-x_k^e}{p^k}\equiv ex_k^{e-1}v\pmod{p} \] 即 \[ v\equiv \frac{a-x_k^e}{ex_k^{e-1}p^k}\pmod{p} \] 所以我们最终得到了一个迭代公式 \[ x_{k+1}=x_k+\frac{a-x_k^e}{ex_k^{e-1}}\pmod{p^{k+1}} \] 代码如下

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
def hensel_lifting(m, c, e, p, k):
    """
    满足pow(m, e, p) == c % p
    求解m % (p ** k)
    """
    assert pow(m, e, p) == c % p
    mk = int(m % p)
    t = int(m * pow(c, -1, p) % p)
    pk = int(p)
    for _ in range(1, k):
        pk *= p
        x = int((c - pow(mk, e, pk)) % pk)
        y = int((x * t * int(pow(e, -1, p))) % pk)
        mk = mk + y
    return mk

可以试着校验

1
2
3
4
5
6
7
8
9
10
from Crypto.Util.number import getPrime

p = getPrime(128)
m = getPrime(512)
e = 0x10001
k = 5
n = p ** k
c = pow(m, e, n)
mk = hensel_lifting(m, c, e, p, k)
print(mk == m)

数学黑洞

123 黑洞

算法

  1. 取一个任意数字串，数出这个数中的偶数个数、奇数个数和总数字个数
  2. 将结果按 偶数个数-奇数个数-总数字个数 排列得到新数
  3. 对新数重复步骤 1、2，最终总会得到数字 123

实例

  1. 取数字串 1234567890，其中偶数个数为 5 (2、4、6、8、0)，奇数个数为 5 (1、3、5、7、9)，总数字个数为 10
  2. 得到新数 5510
  3. 重复步骤，得到新数 134
  4. 重复步骤，得到新数 123

其他

123 数学黑洞（又被称为西西弗斯串），由中国回族学者秋屏先生严格证明，并且推广到六个类似的数字：123、213、312、321、132、231

Kaprekar's routine (Kaprekar 常数)

算法

  1. 取一个任意数字串，将数字串中的数字从大到小排列得到最大值，将数字串中的数字从小到大排列得到最小值
  2. 令最大值减最小值，得到新数
  3. 新数重复步骤 1、2，最终会得到 Kaprekar 不动点（Kaprekar 常数）

实例

当取的任意数字串基于十进制，且为四位数时：

  1. 取数字串 3524，得到最大值 5432，得到最小值 2345
  2. 最大值减最小值，得到新值 3087
  3. 重复步骤，最大值 8730 减最小值 0378，得到新值 8352
  4. 重复步骤，最大值 8532 减最小值 2358，得到新值 6174
  5. 重复步骤，最大值 7641 减最小值 1467，得到新值 6174
  6. 故当任意数串基于十进制，且为四位数时，Kaprekar 不动点（Kaprekar 常数）为 6174

当取的任意数串基于十进制，且为三位数时，Kaprekar 不动点（Kaprekar 常数）为 495

需要注意的是，当位数不足的时候需要补足前导零

其他

实际上，Kaprekar 常数在不同进制下的不同数位中是不同的，我们统一将其称为 Kaprekar 常数。

更多详情可以参考 https://en.wikipedia.org/wiki/Kaprekar%27s_routine

自幂数

算法

取一个共 n 位的数，将各个位上的数的 n 次方相加，得到的累加和结果跟原数相等，我们就称这个数为自幂数。

实例

当取 3 位数时，这样的自幂数称为水仙花数。

例如，取数字 135，各个位上的数的 3 次方分别为 1、27、125，累加和为 135，即跟原数相等。

这样的数只有 135、370、371、407。

除此之外，4 位数的自幂数称为玫瑰花数，这样数有 1634、8208、9474；5 位数的自幂数称为五角星数，这样数有 54748、92727、93084

冰雹猜想

算法

任取一个正整数 N，让它按照以下规则不断递推：

  • 如果这个数是奇数，则下一步变为 3N+1
  • 如果这个数是偶数，则下一步变为 N/2

无论是什么正整数，最终都会使得结果变为 1

实例

冰雹的最大魅力在于不可预知性。

英国剑桥大学教授 John Conway 找到了一个自然数 27。虽然 27 是一个貌不惊人的自然数，但是如果按照上述方法进行运算，则它的上浮下沉异常剧烈：首先，27 要经过 77 步骤的变换到达顶峰值 9232，然后又经过 32 步骤到达谷底值 1。全部的变换过程（称作“雹程”）需要 111 步，其顶峰值 9232，达到了原有数字 27 的 342 倍多，如果以瀑布般的直线下落（2 的 N 次方）来比较，则具有同样雹程的数字 N 要达到 2 的 111 次方。

这个结果是极其惊人的！因为在1到100的范围内，像 27 这样的剧烈波动是没有的（54等27的2的次方倍数的数除外）。

其他

冰雹猜想又叫 Collatz conjecture，之所以翻译为冰雹猜想，是因为这样一个序列像冰雹一样经历多次上升和下降。

同时也被称作角谷猜想，这是因为是一个叫角谷的日本人将这个猜想传入中国的。