二次剩余

对于素数 $p$ 和数 $a$，满足 $(a,p)=1$。（注意 $a$ 不一定小于 $p$）

若 ${\exists }_x{x}^2\equiv a(\mathrm{m}\mathrm{o}\mathrm{d}p)$，则称 $a$ 是模 $p$ 意义下的二次剩余，$x$ 称为该二次剩余方程式的解。

否则称 $a$ 是模 $p$ 意义下的非二次剩余。

特判掉 $p=2$ 的情况，接下来默认 $p$ 均为非 $2$ 质数，即奇素数。

勒让德符号

定义：
$$\left(\frac{a}{p}\right)=\{\begin{array}{l}0p|a\\ 1a在(\mathrm{m}\mathrm{o}\mathrm{d}p)意义下是二次剩余\\ -1a在(\mathrm{m}\mathrm{o}\mathrm{d}p)意义下不是二次剩余\end{array}$$

欧拉准则

$$\left(\frac{a}{p}\right)\equiv a^{\frac{p-1}{2}}(\mathrm{m}\mathrm{o}\mathrm{d}p)$$

证明：

• $p|a$。显然成立。

• 当 $a$ 是在模 $p$ 意义下的二次剩余。

  即 ${\exists }_x{x}^2\equiv a(\mathrm{m}\mathrm{o}\mathrm{d}p)$。

  有 $(x^2{)}^{\frac{p-1}{2}}\equiv a^{\frac{p-1}{2}}(\mathrm{m}\mathrm{o}\mathrm{d}p)\iff x^{p-1}\equiv a^{\frac{p-1}{2}}(\mathrm{m}\mathrm{o}\mathrm{d}p)$。

  $\because (a,p)=1\therefore (x,p)=1$。

  $\Rightarrow x^{p-1}\equiv 1\equiv a^{\frac{p-1}{2}}(\mathrm{m}\mathrm{o}\mathrm{d}p)$。

• 当 $a$ 在模 $p$ 意义下不是二次剩余。

  则有 ${\forall }_{i\in [1,p)}ij\equiv a(\mathrm{m}\mathrm{o}\mathrm{d}p)$ 的 $j$ 是唯一的，且 $j\ne i$。

  考虑用证明逆元唯一的思路来证明这里的唯一性。反证法。

  假设 $0<j_1<j_2<p$ 满足 $i{j}_1\equiv i{j}_2\equiv a(\mathrm{m}\mathrm{o}\mathrm{d}p)$。

  则 $i(j_2-j_1)\equiv 0(\mathrm{m}\mathrm{o}\mathrm{d}p)\Rightarrow p|i(j_2-j_1)$。

  $(i,p)=1\Rightarrow p|(j_2-j_1)$，与 $0<j_2-j_1<p\Rightarrow (j_2-j_1,p)=1$ 矛盾。

  $p-1$ 是偶数，那么我们考虑可以把这些数两两配对，配出 $\frac{p-1}{2}$ 对。

  即 $a^{\frac{p-1}{2}}\equiv (p-1)!\equiv -1(\mathrm{m}\mathrm{o}\mathrm{d}p)$（威尔逊定理）。

---

换言之：

$x^2\equiv a(\mathrm{m}\mathrm{o}\mathrm{d}p),(a,p)=1$。

$a^{\frac{p-1}{2}}\equiv \pm 1(\mathrm{m}\mathrm{o}\mathrm{d}p)$。

方程有解当且仅当 $a^{\frac{p-1}{2}}\equiv 1(\mathrm{m}\mathrm{o}\mathrm{d}p)$。

Cipolla

lemma ：设 $\omega$ 不是模 $p$ 的二次剩余，即 $/{\exists }_x{x}^2\equiv \omega (\mathrm{m}\mathrm{o}\mathrm{d}p)$，且 $\omega =t^2-a$。则 $x=(t+\sqrt{\omega }{)}^{\frac{p+1}{2}}$ 是二次剩余方程 $x^2\equiv a(\mathrm{m}\mathrm{o}\mathrm{d}p)$ 的解。

证明：

$(t+\sqrt{\omega }{)}^p={\sum }_{i=0}^p(\genfrac{}{}{0ex}{}{p}{i})t^i(\sqrt{\omega }{)}^{p-i}=t^p+{\omega }^{\frac{p}{2}}+{\sum }_{i=1}^{p-1}\left(\genfrac{}{}{0ex}{}{p}{i}\right)t^i(\sqrt{\omega }{)}^{p-i}$。

$i\in (0,p),\left(\genfrac{}{}{0ex}{}{p}{i}\right)=\frac{p!}{i!(p-i)!}$。

组合数是整数，$i!,(p-i)!$ 里的每个数都 $<p$，又 $p$ 是奇素数，所以 $i!,(p-i)!$ 里面均不会筛到 $p$。

所以 $p\left|\right(\genfrac{}{}{0ex}{}{p}{i})$。也就是说后面的求和在模 $p$ 意义下是 $0$。

$\Rightarrow (t+\sqrt{\omega }{)}^p\equiv t^p+(\sqrt{\omega }{)}^p(\mathrm{m}\mathrm{o}\mathrm{d}p)$。

$t^p+(\sqrt{\omega }{)}^p=t^p+{\omega }^{\frac{p-1}{2}}\sqrt{\omega }$，因为 $\omega$ 不是模 $p$ 意义下的二次剩余，所以 $\left(\frac{\omega }{p}\right)=-1={\omega }^{\frac{p-1}{2}}$，$t^p=t^{p-1}t\equiv t(\mathrm{m}\mathrm{o}\mathrm{d}p)$。

$\Rightarrow (t+\sqrt{\omega }{)}^p\equiv t-\sqrt{\omega }(\mathrm{m}\mathrm{o}\mathrm{d}p)$

$x^2\equiv (t+\sqrt{\omega }{)}^{p+1}\equiv (t+\sqrt{\omega }{)}^p(t+\sqrt{\omega })\equiv (t-\sqrt{\omega })(t+\sqrt{\omega })\equiv t^2-\omega \equiv a(\mathrm{m}\mathrm{o}\mathrm{d}p)$。

在算法实现中，随机 $t$ 后求 $\omega$，大约有一半的数都是非二次剩余。期望两次便可求出 $x$。

---

考虑另一种二次剩余方程：$x^2\equiv a(\mathrm{m}\mathrm{o}\mathrm{d}{p}^n)$，其中 $(a,p)=1$，$p$ 为奇质数。

先求出 $x^2\equiv a(\mathrm{m}\mathrm{o}\mathrm{d}p)$ 的一个解 $\tau$，有 ${\tau }^2-a=kp\Rightarrow ({\tau }^2-a{)}^n\equiv k^n{p}^n\equiv 0(\mathrm{m}\mathrm{o}\mathrm{d}{p}^n)$。

$\{\begin{array}{l}(r+\sqrt{a}{)}^n=f+g\sqrt{a}\\ (r-\sqrt{a}{)}^n=f-g\sqrt{a}\end{array}$

两个式子的结果答案长相一定可以这么写。

证明：二项式定理展开 $\left(\genfrac{}{}{0ex}{}{n}{i}\right){\tau }^{n-i}(\pm \sqrt{a}{)}^i$，$(\pm \sqrt{a}{)}^i$ 如果 $i$ 是偶数，第 $i$ 项就在 $f$ 里面，否则就在 $g$ 里面。

也就是说，$({\tau }^2-a{)}^n\equiv [(\tau +a)(\tau -a){]}^n\equiv (\tau +a{)}^n(\tau -a{)}^n\equiv (f+g\sqrt{a})(f-g\sqrt{a})\equiv f^2-g^{2}a\equiv 0(\mathrm{m}\mathrm{o}\mathrm{d}{p}^n)$。

$(f,p)=1,(g,p)=1\Rightarrow f^2\text{inv}(g^2)\equiv a(\mathrm{m}\mathrm{o}\mathrm{d}n)$。

由于 $p^n$ 不是质数，所以用 $\text{exgcd}$ 求逆元即可。

代码

#include <cstdio>
#include <iostream>
#include <ctime>
#include <random>
using namespace std;
#define int long long
int T, n, p, omega;
mt19937 wwl(time(0));struct complex {int x, y;complex(){}complex( int X, int Y ) { x = X, y = Y; }complex operator * ( complex v ) {int ansx = ( x * v.x % p + y * v.y % p * omega % p ) % p;int ansy = ( x * v.y % p + y * v.x % p ) % p;return complex( ansx, ansy );}
};int qkpow( complex x, int y ) {complex ans( 1, 0 );while( y ) {if( y & 1 ) ans = ans * x;x = x * x;y >>= 1;}return ans.x;
}int qkpow( int x, int y, int mod ) {int ans = 1;while( y ) {if( y & 1 ) ans = ans * x % mod;x = x * x % mod;y >>= 1;}return ans;
}int cipolla( int n, int p ) {n %= p;if( p == 2 ) return n;if( qkpow( n, p - 1 >> 1, p ) == p - 1 ) return -1;int t;while( 1 ) {t = wwl() % p;omega = ( ( t * t % p - n ) % p + p ) % p;if( qkpow( omega, p - 1 >> 1, p ) == p - 1 ) break;}complex ans( t, 1 );return qkpow( ans, p + 1 >> 1 );
}signed main() {scanf( "%lld", &T );while( T -- ) {scanf( "%lld %lld", &n, &p );if( n == 0 ) { puts("0"); continue; }int ans = cipolla( n, p );if( ! ~ ans ) { puts("Hola!"); continue; }if( ans == p - ans ) printf( "%lld\n", ans );else printf( "%lld %lld\n", min( ans, p - ans ), max( ans, p - ans ) );}return 0;
}