二次剩余（Cipolla）

二次剩余 / 模意义下开根

求解 \(x ^ 2 \equiv n \pmod {p}\)

BSGS

注意到 \(\log x = \frac 1 2 \log n\)，先找出 \(p\) 的原根，再利用 BSGS 求解出 \(\log n\) 即可，若 \(\log n\) 为奇则无解。

复杂度 \(O(\sqrt n)\)。

Cipolla

注意到 BSGS 要找原根麻烦且很慢，所以不好用。

首先对于判断无解，我们知道对于 \(\forall x\) 有 \(x ^ {p - 1} \equiv 1\)。对于 \(n \ne 0\)，设其二次剩余为 \(x\)，那么 \(x \equiv n ^ {\frac 1 2}\)，那么显然 \(x ^ {p - 1} \equiv \pm n ^ \frac {p - 1} 2 \equiv \pm 1\)，显然当 \(x ^ {p - 1} \equiv 1\) 才有解，所以若 \(n ^ \frac {p - 1} 2 \equiv -1\) 则 \(x\) 无解。

现在我们来算 \(x\)。首先我们可以找到一个 \(a\)，使得 \(a ^ 2 - n\) 无二次剩余，记 \(w ^ 2 \equiv a ^ 2 - n\)，由于 \(a ^ 2 - n\) 无二次剩余所以 \(w\) 可以看作一个虚数，计算法则类似复数。我们设 \(x\) 两个解为 \(x_1\) 和 \(x_2\)，显然 \(x_1 \equiv -x_2\)，那么我们只需求解出一个即可。我们注意到设 \(x\) 为其中一个解，那么：

\[x \equiv (a + w) ^ \frac {p + 1} 2 \]

接下来是证明。

将原式化为：

\[x ^ 2 \equiv (a + w) (a + w) ^ p \]

二项式定理展开：

\[x ^ 2 \equiv (a + w) \sum_{k = 0} ^ p \binom{p}{k} a ^ k w ^ {p - k} \]

由于我们的计算是模 \(p\) 意义下的，所以对于 \(\binom{p}{k}\) 只有在 \(k\) 为 \(0\) 或 \(p\) 两项时不为 \(0\) 而为 \(1\)，带入得：

\[x ^ 2 \equiv (a + w) (a ^ 0 w ^ p + a ^ p w ^ 0) \]

对于 \(w ^ p\) 有：

\[w ^ p \equiv w ^ {p - 1} \times w \equiv -w \]

带入有：

\[x ^ 2 \equiv (a + w) (a - w) \equiv a ^ 2 - w ^ 2 \]

带入 \(w ^ 2 \equiv a ^ 2 - n\) 得：

\[x ^ 2 \equiv a ^ 2 - a ^ 2 + n \equiv n \]

得证。

code

int mod;

il int ksm(int a, int k = mod - 2) { int res = 1; for(; k; k >>= 1, a = a * a % mod) if(k & 1) res = res * a % mod; return res; }

int w2;
struct num {
    int x, y;
    num(int X = 0, int Y = 0) { x = X, y = Y; }
    il friend num operator +(const num &a, const num &b) { return num((a.x + b.x) % mod, (a.y + b.y) % mod); }
    il friend num operator -(const num &a, const num &b) { return num((a.x + mod - b.x) % mod, (a.y + mod - b.y) % mod); }
    il friend num operator *(const num &a, const num &b) { return num((a.x * b.x % mod + a.y * b.y % mod * w2 % mod) % mod, (a.x * b.y % mod + a.y * b.x % mod) % mod); }
    il friend num operator ^(num a, int k) { num res = num(1, 0); for(; k; k >>= 1, a = a * a) if(k & 1) res = res * a; return res; }
} ;

il pii cipolla(int n) {
    if(!n) return {0, 0};
    if(ksm(n, (mod - 1) / 2) ^ 1) return {-1, -1};
    int a;
    do { a = 1ll * rand() * rand() % mod; } while(ksm(Mod(a * a % mod - n), (mod - 1) / 2) == 1);
    w2 = Mod(a * a % mod - n);
    num x = (num(a, 1) ^ ((mod + 1) / 2));
    int x1 = x.x, x2 = mod - x1;
    if(x2 < x1) swap(x1, x2);
    return {x1, x2};
}