[萌熊]染色（题解）

题意

给定正整数 \(n\)，素数 \(p\) 和一个矩阵 \(A\)，求 \(\sum\limits_{u\in [0,p)^n} \left[p\mid (u^TAu)\right]\)，对 \(10^9+7\) 取模。

其中 \(n\le 500,p<10^9+7,A_{i,j}\in [0,p)\)。

题解

类似题目 ABC220H。

线性代数部分

由于这是个二次型，有：\(u^TAu=\sum\limits_{i,j} A_{i,j}u_iu_j\)。

所以我们可以始终钦定 \(A_{j,i}\gets A_{i,j}+A_{j,i}(i<j)\)，这不影响结果。

考虑可以对 \(u\) 乘任意可逆矩阵 \(Q\)，此时 \(v=Qu\) 依然遍历 \([0,p)^n\)

此时 \(v^TAv=u^T(Q^TAQ)u=u^TA'u\)，其中 \(A'=Q^TAQ\)。

我们有结论：

• 当 \(p=2\) 时，存在 \(Q\) 使得 \(A'\) 为下双对角矩阵。

• 当 \(p>2\) 时，存在 \(Q\) 使得 \(A'\) 为对角矩阵。

具体请看这里哦！

---

当 \(p=2\) 时，由于只有 \(A_{i+1,i},A_{i,i}\) 处有值，直接 dp 即可。

\(f_{i,j,k}\) 表示考虑前 \(i\) 位，之前和 \(\bmod 2\) 是 \(j\)，上一个数是 \(k\) 的方案数，dp 是 \(O(n)\) 的，总复杂度 \(O(n^3)\)。

---

当 \(p\) 是奇素数时，设 \(A\) 矩阵已经消成对角阵 \(\text{diag}(d_1,d_2,\cdots ,d_n)\)。

欲求 \(\sum\limits_{u\in [0,p)^n} \left[p\mid \sum\limits_{i=1}^n d_iu_i^2\right]\)，对 \(10^9+7\) 取模。

数论部分

下面默认多项式 \(\bmod p\) 意义下循环：即指数可任意 \(\bmod p\)​。乘法也是循环卷积。

考虑多项式：\(f_d(x)=\sum\limits_{k=0}^{p-1} x^{dk^2}\)，我们转换为求 \([x^0]\prod\limits_{i=1}^n f_{d_i}(x)\)。

考察二次剩余的性质，显然：

• 当 \(d\) 为二次剩余时，\(dk^2\) 遍历 \(0\) 和两遍二次剩余集。

• 当 \(d\) 为非二次剩余时，\(dk^2\) 遍历 \(0\) 和两遍非二次剩余集。

任取一个非二次剩余 \(r\)，此时记 \(f(x)=f_1(x),g(x)=f_r(x)\)。

设 \(d\) 中有 \(u\) 个二次剩余，\(v\) 个非二次剩余，\(w\) 个 \(0\)，此时答案为：\(p^w[x^0]f^u(x)g^v(x)\)。

此时直接多项式乘法做，能做到 \(O(p\log n\log p)\)。

但是我们大胆猜测 \(f^n(x)\) 和 \(g^n(x)\) 有特殊性质。

下面暂时不加证明的给出结论，过程太冗长后面再写，这些结论确实可以打表发现规律：

\[\boxed{P=f^2(x)=g^2(x)=A+B\sum\limits_{i=1}^{p-1}x^i\\Q=f(x)g(x)=A^{*}+B^{*}\sum\limits_{i=1}^{p-1}x^i} \]

• 当 \(p=4k+1\) 时，\(A=2p-1,B=p-1,A^*=1,B^*=p+1\)

• 当 \(p=4k+3\) 时，\(A^*=2p-1,B^*=p-1,A=1,B=p+1\)

于是 \(f^ug^v=P^U\times f^{o_1}g^{o_2},U=u/2+v/2,o_1=u\bmod 2,o_2=v\bmod 2\)。

我们再做个记号：\(P(A,B)=A+B\sum\limits_{i=1}^{p-1} x^i,T(A,B)=A\sum\limits_{i=0}^{p-1} x^{i^2}+B\sum\limits_{i=0}^{p-1} x^{ri^2}\)。

---

怎么算 \(P^n\) 呢？容易打表观察到 \(P^n\) 也是符合 \(P(A,B)\) 形式，于是求解如下问题：

给定 \(p,A,B\)，令 \(f(x)=P(A,B)\)，已知 \(f^n(x)=P(A_n,B_n)\)，求 \(A_n,B_n\)。

记 \(S=A+(p-1)B,D=A-B\)

\[\boxed{ \begin{aligned} A_n \;=\; \frac{\,S^n \;+\;(p-1)\,D^n\,}{p}, \qquad B_n \;=\; \frac{\,S^n \;-\; D^n\,}{p}. \end{aligned} } \]

类似的，我们有：\(P(A,B)P(C,D)=P(AC+(p-1)BD,AD+BC+(p-2)BD)\)

---

于是若 \(o_1=o_2=0\)，则 \(p^uq^v=P^{U}\)，直接计算即可。

若 \(o_1=o_2=1\)，则 \(p^uq^v=P^{U}\times Q\)，直接计算即可。

若 \(o_1,o_2\) 恰有一个 \(1\)，则 \(p^uq^v=P^U\times T(o_1,o_2)\)。此时我们依然有结论：

\[\boxed{ P(A,B)\times T(C,D) \;=\;\sum_{k=0}^{p-1} h_k\,x^k, \quad \text{其中} \quad h_k = \begin{cases} (A-B + B\,p)\,(C+D), &k\equiv0\pmod p,\\[6pt] 2C\,(A-B)\;+\;B\,p\,(C+D), &k\neq0\text{ 且 }k\text{ 为二次剩余},\\[6pt] 2D\,(A-B)\;+\;B\,p\,(C+D), &k\text{ 为二次非剩余}. \end{cases} } \]

于是答案就是 \(A+(p-1)B\)。

至此我们总复杂度 \(O(n^3)\) 的完成了这题，可能有一些 \(\text{polylog}(n,p)\) 加在后面，当然这是小量直接忽略。

注意：\(p^w\) 不要忘记乘了！！！

$\bf{code}$

#include<bits/stdc++.h>
#define LL long long
#define fr(x) freopen(#x".in","r",stdin);freopen(#x".out","w",stdout);
using namespace std;
const int N=505,mod=1e9+7;
int n,p,a[N][N],f[N][79][79],g[N][2][2],A[N],B[N],s;
inline void AD(int &x,int y){x+=y;(x>=p)&&(x-=p);}
inline void ad(int &x,int y){x+=y;(x>=mod)&&(x-=mod);}
inline int ksm(int x,int p,const int mod){int s=1;for(;p;(p&1)&&(s=1ll*s*x%mod),x=1ll*x*x%mod,p>>=1);return s;}
inline int inv(int x,const int mod){return ksm(x,mod-2,mod);}
int main()
{
    fr(graph)
    ios::sync_with_stdio(0);cin.tie(0);cout.tie(0);cin>>n>>p;
    for(int i=1;i<=n;i++) for(int j=i;j<=n;j++) cin>>a[j][i];
    for(int i=1;i<n;i++)
    {
        for(int j=i+1;j<=n;j++) AD(a[j][i],a[i][j]),a[i][j]=0;
        for(int j=i+2;j<=n;j++) if(a[j][i])
        {
            swap(a[i+1],a[j]);
            for(int k=1;k<=n;k++) swap(a[k][i+1],a[k][j]);
        }
        if(!a[i+1][i]) continue;
        int I=p-inv(a[i+1][i],p);
        for(int j=i+2;j<=n;j++)
        {
            int t=a[j][i]*I%p;
            for(int k=1;k<=n;k++) a[j][k]=(a[j][k]+a[i+1][k]*t)%p;
            for(int k=1;k<=n;k++) a[k][j]=(a[k][j]+a[k][i+1]*t)%p;
        }
    }
    if(p^2)
    {
        const int I2=(p+1)>>1;
        for(int i=1;i<=n;i++) A[i]=a[i][i],B[i]=a[i+1][i]*I2%p,a[i][i]=a[i+1][i]=0;
        for(int i=1,x,t;i<n;i++) if(x=B[i])
        {
            if(!A[i]&&!A[i+1])
            {
                A[i]=x=2*x%p,A[i+1]=p-x,t=B[i+1]*inv(x,p)%p,B[i+1]=p-B[i+1],A[i+2]=(A[i+2]+B[i+1]*t)%p;
            }
            else
            {
                if(!A[i]) swap(A[i],A[i+1]),B[i+1]=0;
                A[i+1]=(A[i+1]+(p-x)*x*inv(A[i],p))%p;
            }
        }
		int u=0,v=0,w=0,t,c;
		for(int i=1;i<=n;i++) A[i]?(ksm(A[i],(p-1)/2,p)==1?u++:v++):w++;
		t=u/2+v/2;u&=1,v&=1;
		int _A=1,_B=p+1;
		if(p%4==1) _A=2*p-1,_B=p-1;
		int S=(_A+(p-1)*_B)%mod,D=(_A-_B+mod)%mod,ip=inv(p,mod);
		_A=(ksm(S,t,mod)+1ll*(p-1)*ksm(D,t,mod))%mod*ip%mod;
		_B=1ll*(ksm(S,t,mod)-ksm(D,t,mod)+mod)*ip%mod;
		if(u&v)
		{
			int _a=1,_b=p+1;
			if(p%4==3) _a=2*p-1,_b=p-1;
			s=(1ll*_a*_A+1ll*_b*_B%mod*(p-1))%mod;
		}
		else s=(u|v)?(_A+1ll*(p-1)*_B)%mod:_A;
		cout<<1ll*s*ksm(p,w,mod)%mod;
    }
	else
	{
		for(int i=0;i<2;i++) g[1][i][a[1][1]&i]=1;
		for(int i=1;i<n;i++) for(int j:{0,1}) for(int s:{0,1})
			if(g[i][j][s]) for(int k:{0,1})
				ad(g[i+1][k][s^(a[i+1][i]&j&k)^(a[i+1][i+1]&k)],g[i][j][s]);
		for(int i:{0,1}) ad(s,g[n][i][0]);cout<<s;
	}
    return 0;
}

Bonus

为啥我们要算整个多项式呢？

因为做这个题，事实上我们只需要知道答案多项式的 \([x^0]\) 即可。

下面我们解决一般性问题：给定 \(n\times n\) 矩阵 \(A\)，素数 \(p\)，一个函数 \(F\)，求 \(\sum\limits_{u\in [0,p)^n} F\left((u^TAu)\bmod p\right)\) 的问题。

\(p=2\) 无所谓，反正都是 \(O(n^3)\) 的，考虑 \(p>2\)。

由于答案多项式部分，\([x^0],[x \text{ 的二次剩余幂次}],[x \text{ 的非二次剩余幂次}]\) 分别有有三个相同的系数 \(A,B,C\)。

于是我们若能知道 \(F(0),\sum F(\text{所有二次剩余}),\sum F(\text{所有非二次剩余})\)，直接乘上 \(A,B,C\) 三个系数即可。

构造特殊的 \(F\) 函数即可加强此题。一般的 \(F\) 也能做到 \(O(n^3+p)\) 的复杂度，十分优秀。

关于那些二次剩余相关多项式的式子证明

请看这里。这里证明了一些东西的乘法封闭，并给出乘完的结果！可以通过那里最后的式子稍加推到得到本文的结论。