P17226 [Math×Girl²] 终末电台题解
2021hych
·
·
题解
::::info[题意]{open}
记 \bmod p 意义下的非零二次剩余组成的集合为 S。对于每组询问 (k,n)，求不定方程：
x_1+x_2+\cdots+x_k \equiv n (\bmod \, p)\\
x_i \in S,i=1,2,\cdots,k
\end{cases}
在 \bmod p 意义下解的个数，记为 f_{n,k}。
::::
解答
你选择由浅入深地考虑问题。k=0 平凡无趣，考虑 k=1，这变成了一个判断 n 是否为模 p 的二次剩余。学过 \text{Legendre} 符号的你将其表示为：
f_{n,1}=\frac{1}{2}\left(1+\left( \frac{n}{p}\right)-\delta_{n,0}\right)
其中 \delta_{n,0} 在 n=0 时取值为 1，其余情况均为 0。
$$f_{n,2}=\sum_{a\in S}f_{n-a,1}=\frac{1}{2}\left[\sum_{a\in S}1+\sum_{a\in S}\left( \frac{n-a}{p}\right)-\sum_{a\in S}\delta_{n-a,0}\right]$$
依次观察括号中的三项，第一项显然等于 $|S|$，熟知缩系中模 $p$ 的二次剩余恰好占一半（后续的和式化简会反复用到这个性质，本文中，如果一个和式突然变成了 $0$ 就意味着该性质被使用），从而第一项的结果为 $\dfrac{p-1}{2}$。
先看第三项，显然这又变成了判断 $n$ 是否为模 $p$ 的二次剩余，其值为 $\dfrac{1}{2}\left(1+\left( \dfrac{n}{p}\right)-\delta_{n,0}\right)$。
最后看第二项，幸好，你见过如何处理这种求和式，给出了套路化的化简如下：(用到了 $\left( \dfrac{n}{p}\right)$ 的完全积性)
$$\sum_{a\in S}\left( \frac{n-a}{p}\right)=\sum_{a=1}^{p-1}\left(1+\left( \frac{a}{p}\right)\right)\left( \frac{n-a}{p}\right)=\dfrac{1}{2}\left[\sum_{a=0}^{p-1}\left( \frac{n-a}{p}\right)-\left( \frac{n}{p}\right)+\sum_{a=1}^{p-1}\left( \frac{na-a^2}{p}\right)\right]$$
$$=\frac{1}{2}\left[\sum_{a=1}^{p-1}\left( \frac{(a^{-1})^2}{p}\right)\left( \frac{na-a^2}{p}\right)-\left( \frac{n}{p}\right)\right]=\frac{1}{2}\left[\sum_{a=1}^{p-1}\left( \frac{na^{-1}-1}{p}\right)-\left( \frac{n}{p}\right)\right]$$
$$=\frac{1}{2}\left[\sum_{b=0}^{p-1}\left( \frac{nb-1}{p}\right)-\left( \frac{-1}{p}\right)-\left( \frac{n}{p}\right)\right]$$
如果 $n=0$，最后一个求和式就是 $p$ 个 $\left( \dfrac{-1}{p}\right)$ 的和，否则 $nb-1$ 会通过模 $p$ 的一个缩系和 $0$，求和式为 $0$。于是第二项的值为：
$$\frac{1}{2}\left((p\delta_{n,0}-1)\left( \frac{-1}{p}\right)-\left( \frac{n}{p}\right)\right)$$
你不着急代入，选择观察一下式子有什么特点。你发现，如果把 $\left( \dfrac{-1}{p}\right)$ 视作常数的话，似乎 $f_{n,1}$ 和 $f_{n,2}$ 都可以写成 $a\left( \dfrac{n}{p}\right)+b\delta_{n,0}+c$ 的形式，其中 $a,b,c$ 为常数。你又观察到 $k \leqslant 10^{18}$ 的数据范围。一个大胆的想法呼之欲出：你想把 $f_{n_,k}$ 写成 $a_k\left( \dfrac{n}{p}\right)+b_k\delta_{n,0}+c_k$ 的形式尝试，找出 $a_k,b_k,c_k$ 与 $a_{k+1},b_{k+1},c_{k+1}$ 之间的线性递推关系，使用矩阵快速幂加速递推！
你的观察很敏锐，事实上:
$$f_{n,k+1}=\sum_{a\in S}f_{n-a,k}=a_k\sum_{a\in S}\left( \frac{n-a}{p}\right)+b_k\sum_{a\in S}\delta_{n-a,0}+c_k\sum_{a\in S}1$$
$$=a_k \cdot \frac{1}{2}\left((p\delta_{n,0}-1)\left( \frac{-1}{p}\right)-\left( \frac{n}{p}\right)\right)+b_k \cdot \dfrac{1}{2}\left(1+\left( \dfrac{n}{p}\right)-\delta_{n,0}\right) +c_k \cdot \frac{p-1}{2}$$
$$=\left(-\frac{1}{2}a_k+\frac{1}{2}b_k\right)\left( \dfrac{n}{p}\right)+\left(\frac{1}{2}\left( \frac{-1}{p}\right)pa_k-\frac{1}{2}b_k\right)\delta_{n,0}+\left[-\frac{1}{2}\left( \frac{-1}{p}\right)a_k+\frac{1}{2}b_k+\frac{p-1}{2}c_k\right]$$
故：
$$\begin{cases}
a_{k+1}=-\dfrac{1}{2}a_k+\dfrac{1}{2}b_k\\
b_{k+1}= \dfrac{1}{2}\left( \dfrac{-1}{p}\right)pa_k-\dfrac{1}{2}b_k\\
c_{k+1}=-\dfrac{1}{2}\left( \dfrac{-1}{p}\right)a_k+\dfrac{1}{2}b_k+\dfrac{p-1}{2}c_k
\end{cases}$$
写成矩阵形式如下：
$$\begin{pmatrix}
a_{k+1} & b_{k+1} & c_{k+1}
\end{pmatrix}=
\begin{pmatrix}
a_k & b_k & c_k
\end{pmatrix}
\begin{pmatrix}
-\dfrac{1}{2} & \dfrac{1}{2}\left( \dfrac{-1}{p}\right)p & -\dfrac{1}{2}\left( \dfrac{-1}{p}\right)\\
\dfrac{1}{2} & -\dfrac{1}{2} & \dfrac{1}{2} \\
0 & 0 & \dfrac{p-1}{2}\\
\end{pmatrix}
$$
使用矩阵快速幂即可。
---
至此，你的猜测得到了验证，总时间复杂度是 $O(Q(27\log k+\log p))$，里面的两个 $\log$ 分别对于矩阵快速幂和 $\text{Legendre}$ 符号计算。
# 参考代码
```cpp
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int mod=998244353;
const int inv2=(mod+1)/2;
int Q,p,n,s,ans; ll k;
struct mat {
int a[4][4];
int r,c;
mat() {
memset(a,0,sizeof(a));
}
};
mat operator * (mat x,mat y) {
mat z;
z.r=x.r,z.c=y.c;
for(int i=1;i<=x.r;i++)
for(int j=1;j<=y.c;j++)
for(int k=1;k<=x.c;k++)
z.a[i][j]=(z.a[i][j]+(ll)x.a[i][k]*y.a[k][j])%mod;
return z;
}
mat ksm(mat x,ll y) {
mat ret=x;
y--;
while(y) {
if(y&1) ret=ret*x;
x=x*x;
y>>=1;
}
return ret;
}
int power(int a,int b,int p) {
int ans=1;
for(;b;b>>=1) {
if(b&1) ans=(ll)ans*a%p;
a=(ll)a*a%p;
}
return ans;
}
mat sta,Ans,G;
int main() {
scanf("%d%d",&Q,&p);
if(p%4==1) s=1;
else s=-1;
G.r=G.c=3;
G.a[1][1]=G.a[2][2]=(-inv2+mod)%mod;
G.a[2][1]=G.a[2][3]=inv2;
G.a[1][2]=((ll)(inv2*s+mod)%mod*p)%mod;
G.a[1][3]=(-inv2*s+mod)%mod;
G.a[3][3]=(ll)(p-1)*inv2%mod;
sta.r=1,sta.c=3;
sta.a[1][1]=sta.a[1][3]=inv2;
sta.a[1][2]=(-inv2+mod)%mod;
while(Q--) {
scanf("%lld%d",&k,&n);
if(!k) {
if(!n) ans^=1;
continue;
}
if(k==1) Ans=sta;
else Ans=sta*ksm(G,k-1);
s=power(n,(p-1)/2,p);
if(s==p-1) s=-1;
ans^=((ll)((ll)(Ans.a[1][1]*s+mod)%mod+Ans.a[1][2]*(n==0))%mod+Ans.a[1][3])%mod;
}
printf("%d",ans);
return 0;
}
```