模板 - 二次剩余

求解方程 $ x^2 = a \space mod \space p $ ，其中p是质数，无解返回-1，有解则返回两个解之中较小的那一个。事实上假如p是奇质数的话，两个解相等仅当x=0的时候，完全可以特判掉。

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;

ll p, w;
struct Complex {
    ll x, y;
    Complex() {};
    Complex(ll x, ll y): x(x), y(y) {}
    Complex operator*(const Complex&c) {
        return Complex((x * c.x % p + y * c.y % p * w % p) % p, (x * c.y % p + y * c.x % p) % p);
    }
};

Complex qpow(Complex x, ll n) {
    Complex res(1, 0);
    while(n) {
        if(n & 1)
            res = res * x;
        x = x * x;
        n >>= 1;
    }
    return res;
}

ll qpow(ll x, ll n, ll m) {
    ll res = 1;
    while(n) {
        if(n & 1)
            res = res * x % m;
        x = x * x % m;
        n >>= 1;
    }
    return res;
}

//求解 x^2 = a mod p 的x的两个p范围内的解中的其中一个解，其中p是质数（2被特判掉了），无解返回-1
ll solve(ll a, ll p) {
    a %= p; //可能可以省略
    if(a == 0)
        return 0;
    if(p == 2)
        return a;
    if(qpow(a, (p - 1) / 2, p) == p - 1)
        return -1;
    ll b, t;
    while(1) {
        b = rand() % p;
        t = b * b - a;
        w = (t % p + p) % p;
        if(qpow(w, (p - 1) / 2, p) == p - 1)
            break;
    }
    Complex res(b, 1);
    res = qpow(res, (p + 1) / 2);
    ll x = res.x, y = (p - x) % p;
    return x <= y ? x : y;
}

int main() {
#ifdef Inko
    freopen("Inko.in", "r", stdin);
#endif // Inko
    int T;
    scanf("%d", &T);
    while(T--) {
        ll a;
        scanf("%lld%lld", &a, &p);
        ll x = solve(a, p);
        printf("%lld\n", x);
    }
}