显然,我们只需考虑$0 \leq a < n $ 的情形。这个问题应该很早就被人考虑好了,不过无所谓吧(反正只是个博文而已)。以下内容都是独立完成的。并可看作二次剩余和Gauss互反律 这篇博文的延续。最后再给出方程的一个解(如果存在的话)。

记$A_n = \{0<x<n \;|\; \gcd(x,n) = 1\}$,则$A_n$ 关于$\mod n$ 的乘法构成一个群(根据Euler定理$|A_n| = \psi(n)$)

所以若$\gcd(a,n)= 1$,则任意$x \in A_n$,存在唯一$y \in A_n$ 使得$xy \equiv a \mod n$

若 $x^2 \equiv a \mod n$ 无解,则 $\prod_{x \in A_n} x = a^{\frac{\psi(n)}{2}} \mod n$

若 $x^2 \equiv a \mod n$ 有解,设$0<a_1<a_2<\cdots<a_r<n$ 是全部解。则 $\prod_{x \in A_n} x = a^{\frac{\psi(n)-r}{2}} a_1 \cdots a_r \mod n$,注意到$r$只与$n$有关,与$a$ 无关。

用$[\frac{a}{n}] = 1$表示有解,用$[\frac{a}{n}] = 0$ 表示无解。

为了考虑这个问题,我们分情况讨论

$x^2 \equiv a \mod 2^m, \gcd(a,p^m) = 1$ 何时有解

当$m=1,2$时,有解当且仅当$a=1$

当$m \geq 3$ 时,$2^m | (a_2 -a_1)(a_2+a_1)$, 若$(a_2-a_1)$和$a_2+a_1$都是$4$的倍数,则$a_1,a_2$ 也是$2$的倍数矛盾于$\gcd(a,2^m) = 1$。所以$2^{m-1}|(a_2+a_1)$或者 $2^{m-1}|(a_2-a_1)$。所以$0 < a_1 < 2^{m-2}$,$a_2 = 2^{m-1}-a_1,a_3 = 2^{m-1}+a_1,a_4 = 2^m-a_1$。

但是,若$x$是奇数,则$x^2 \equiv 1\mod 8$ ,令一方面$\gcd(a,2^m) = 1$ 当且仅当$a \equiv \mod 2$,而且每一个解对应4个$x$。

即 $1^2,2^2,\cdots,(2^m-1)^2 \mod 2^m$,只有 $2^{m-3}$ 个数,但是 小于$2^m$且模8为1的数也只有,$2^{m-3}$个。因此,

有解当且仅当 $a \equiv 1 \mod 8$

$x^2 \equiv a \mod p^m,\;p>2, \gcd(a,p^m) = 1$ 何时有解

因为 $p^m|(a_2-a_1)(a_2+a_1)$,若$(a_2-a_1)$和$a_2+a_1$都是$p$的倍数,则$a_1,a_2$ 也是$p$的倍数矛盾于$\gcd(a,p^m) = 1$,所以 $a_2 = p^m - a_1$。从而$r=2,\; a_1a_2 \equiv -a \mod p^m$,所以

  • 若 $x^2 \equiv a \mod p^m$ 无解,则 $\prod_{x \in A_n} x = a^{\frac{p^{m-1}(p-1)}{2}} \mod p^m$

  • 若 $x^2 \equiv a \mod p^m$ 有解,则 $\prod_{x \in A_n} x = -a^{\frac{p^{m-1}(p-1)}{2}} \mod p^m$

取$a=1$,则$\prod_{x \in A_n} x = -1 \mod p^m$,即

$x^2 \equiv a p^i \mod p^m, \; i<m, \; \gcd(a,p^m) = 1$ 何时有解

若$i=2k-1$ 为奇数,则$p^k|x$,从而 $p^{i+1} | ap^i$,矛盾于$ \gcd(a,p^m) = 1$

所以$i$ 为偶数,显然此时上式有解当且仅当 $x^2 \equiv a \mod p^{m-i},\gcd(a,p^{m-i}) = 1$ 有解

综上

$x^2 \equiv a \mod m_1m_2,\gcd(m_1,m_2)=1$何时有解

考虑

则,存在 $t_1,t_2$ 使得,$t_1m_1 + t_2m_2 = 1$,所以 $t_1^2m_1^2 + t_2^2 m_2^2 + 2t_1m_1t_2m_2=1$, 所以 $x \equiv (t_2m_2 a_1)^2$ 且$x \equiv (t_1m_1a_2)^2$,所以

令一方面若 $x \equiv a^2 \mod m_1m_2$,则 $x \equiv (a \mod m_i)^2 \mod m_i$

即 $x^2 = a \mod m_1m_2$ 有解当且仅当 $[\frac{a}{m_1}] = [\frac{a}{m_2}]=1$

最终方案:$n = p_1^{m_1} \cdots p_r^{m_r},\;p_1< \cdots <p_r$, $[\frac{a}{m}] = \prod [\frac{a}{p_i^{m_i}}]$

判断是否有解的sagemath代码:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
# modSqrtSymbol.ipynb
def modSqrtSymbol(an,n):
for p,m in factor(n):
a = an%p^m
if a == 0: continue
flag = True
while a%p == 0:
a//=p
m-=1
flag = not flag
if not flag: return False
if p == 2:
if a%8 != 1: return False
else:
if power_mod(a,p^(m-1)*(p-1)//2,p^m) != 1:
return False
return True

因为factor复杂度是$O(\sqrt{n})$, 所以整体复杂度 $O(\sqrt{n})$

求出 $x^2 \equiv a \mod n$ 的一个解

本来是想求全部解的,但是考虑到如果 $\gcd(a,n) \neq 1$ 则解的个数实在太多,就不想求全部解了,当然如果真的有这个需要,求全部解也不难。

另外,无论是是否有解还是求全部解都可以直接暴力做,只是复杂度为$O(n)$ ,效率太低。

求$x^2 \equiv a \mod 2^m$ 的全部解,不妨设 $a \equiv 1 \mod 8 $

若 $m=1$,则$x=1$,若$m=2$,则$x=1,3$,若$m=3$,则$x=1,3,5,7$,若$m \geq 4$,不妨设$x_{m-1}$ 是满足方程$x^2 \equiv a \mod 2^{m-1}$的最小正整数解,即 $x_{m-1}^2 = a + 2^{m-1}k$

  • 若$k$ 是偶数,则必然是$x_{m-1}^2 \equiv a \mod 2^m$ 的解,最小的解就是 $\min(x_{m-1},2^{m-2}-x_{m-1})$,但由于$x_{m-1} \leq 2^{m-3}$,所以 $x_m=x_{m-1}$

  • 若$k$为奇数,则$(2^{m-2}-x_{m-1})^2 = 2^{2m-4}-2^{m-1}+x_{m-1}^2 \equiv a \mod 2^m$,且 $0<2^{m-2}-x_{m-1}<2^{m-2}$,所以$x_m =2^{m-2}-x_{m-1}$

所以经过 $m$ 步就可找到所有解就是 $x_m,2^{m-1}-x_m,2^{m-1}+x_m,2^m-x_m$

求$x^2 \equiv a \mod p^m,\;p>2$ 的全部解,其中 $\gcd(a,p^m) = 1$

令$q = \frac{\psi(p^m)}{2} = p^{m-1}\frac{p-1}{2}$,若$q$ 是奇数,则$(a^{\frac{q+1}{2}})^2 = a^q a \equiv a\mod p^m$,从而解就是$a^{\frac{q+1}{2}},p^m-a^{\frac{q+1}{2}}$

否则,记 $q = 2^{i}q’$,(想模仿上面q是奇数的情况),我们找一个非二次剩余$b$,即$b^q \equiv -1 \mod p^m$。我们记 $x = a^{\frac{q’+1}{2}},\;y = b^{q’}, \; t= a^{q’}$。 则

$x^2 = ta, t^{2^i} = 1, y^{2^i}=-1$。我们想保持$x^2 = ta$,让$t$ 变成1。就打到了我们的目的。

若$t^{2^{s}} = -1, t^{2^{s+1}} = 1$, 则$s<i$,那么$(y^{2^{i-s}}t)^{2^{s}} = y^{2^i}t^{2^{s}} = 1$,并且 $(y^{2^{i-s-1}}x)^2 = (y^{2^{i-s}}t)a$,我们更新$t = y^{2^{i-s}}t$,更新 $x=y^{2^{i-s-1}}x$。重复上述过程直到 $t=1$。

我们也可以类似$p=2$的情况,经过$m$步得到全部解,只是$m=1$的时候是困难的,并且和上述分析一致,所以就没用$p=2$时的方法。若$x_{m-1}$满足 $x_{m-1}^2 = a + p^{m-1}k$

  • 若$k$ 是偶数,$(x_{m-1} - \frac{k}{2}p^{m-1})^2 \equiv a \mod p^m$

  • 若$k$ 是奇数,$(x_{m-1} + \frac{p-k}{2}p^{m-1})^2 \equiv a \mod p^m$

对一般的 $n = p_1^{m_1} \cdots p_r ^{m_r}$,我们对每个$(a, p_i ^m)$ 给出一个解,再用中国剩余定理 拼出最终解。

$x^2 \equiv a \mod n$ 的一个解的Sagemath 代码

需要用到上面的 modSqrtSymbol 函数

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
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
# modSqrt.ipynb

# find one x s.t. x^2 = a mod 2^m, gcd(a,2) = 1 and 0 < a < 2^m and must have an answer
def modSqrt2Core(a,m):
x = 1
for i in range(4,m+1):
if (x*x-a)%(1<<i) != 0:
x = (1<<i-2)-x
return x

# find one x s.t. x^2 = a mod p^m, gcd(a,p)=1, p>2 and 0 < a < p^m and must have an answer
def modSqrtPCore(a,p,m):
n = p^m
q = p^(m-1)*(p-1)>>1
if q%2 == 1: return power_mod(a,(q+1)>>1,n)
b=randint(0,n-1)
while(power_mod(b,q,n)==1): b=randint(0,n-1)
ni = q&(-q)
q //= ni
x = power_mod(a,(q+1)>>1,n)
y = power_mod(b,q,n)
t = power_mod(a,q,n)
# x^2 = ta mod n, y^ni = -1 mod n
while t != 1:
ns = 1;tt=t*t;
while tt!=1:
tt=(tt*tt)%n
ns<<=1
# t^ns = -1 mod n
t = power_mod(y,ni//ns,n)*t%n
x = power_mod(y,ni//(ns*2),n)*x%n
return x

# find one x s.t. x^2 = a mod p^m and must have an answer
def modPowerSqrt(a,p,m):
n = p^m
a %= n
if a == 0: return 0
time = 0
while a%p==0:
a//=p
m-=1
time+=1
if a == 0: return 0
return 2**(time//2)*(modSqrt2Core(a,m) if p==2 else modSqrtPCore(a,p,m))

def modSqrt(a,n):
if not modSqrtSymbol(a,n): return []
fact = list(factor(n))
x = [modPowerSqrt(a,p,m) for p,m in fact]
m = [p^m for p,m in fact]
return crt(x,m)


# 测试几个样例
a,m = 17,2**8
x = modSqrt(a,m)
print(x,x*x%m,a,m)

a,m = 15,7**8
x = modSqrt(a,m)
print(x,x*x%m,a,m)

a,m = 16,3*5*5*7*7*7
x = modSqrt(a,m)
print(x,x*x%m,a,m)