首先我们不难得出一个矩阵快速幂的做法。不妨设转移矩阵为一个$n*n$的矩阵$A$，那么我们要求的即为$\begin{bmatrix}1 1 … 1\end{bmatrix}*A^T* \begin{bmatrix}1\\1\\…\\1\\\end{bmatrix}$，由$Hamilton-Cayley$定理可知，一个矩阵存在一个$n-1$次特征多项式$A^{n}=c_{1}*A^{n-1}+c_{2}*A^{n-2}+…+c_{n}*A^{0}$。但是直接求矩阵的特征多项式是$N*4$的，在这里并不行。考虑到我们要求一个矩阵乘以$A^T$再乘以一个矩阵，发现这样其实就是在特征多项式两边乘了两个矩阵，这样并不会改变这个$A^{n}$的递推性质，即这递推式依旧成立。现在我们将矩阵转化成了一个数，这样求递推式就可以更快得求了。先通过$N*M$的$DP$预处理出答案前几项，然后我们通过$BM$算法($Berlekamp-Massey$)可以在$N^2$的时间复杂度内求出这个递推式$c_{i}$。求出递推式后，便是一个这样的问题，已知一个数列前$n$项和$n$阶递推式，求次数列的第$t$项。
这个问题可以用倍增+多项式取模完成。即求$x^{T} mod f(x)$,$f(x)$为特征多项式，求出来后乘以数列前$n$项即可（模板题是$BZOJ4162: shlw loves matrix II$）。那么这题就做完了。
时间复杂度$O(N*M+NlogNlogT)$（第一次出题，感觉很有趣。多项式倍增+取模常数较大，原来$T$想放$2^{10000}$，但是跑得超级慢，于是$T$放了$2^{100}$）。
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 |
#include<cstdio> #include<algorithm> #include<cmath> #include<cstring> #include<vector> #define ll long long #define G 3 #define mo 998244353 using namespace std; const int N=30010,M=5010; int n,m,tot,Next[M*2],head[M],tree[M*2],f[2][N],fail[N],delta[N],ans[N],K,MOD[N],Len,T[10010]; vector<int>A,B; char s[10010]; namespace poly { int a[N],b[N],w[N],L,INV; ll mi(ll x,int y) { ll ans=1; while (y) { if (y&1) ans=ans*x%mo; x=x*x%mo; y>>=1; } return ans; } void NTT(int *A,int f) { for (int i=0,j=0;i<L;i++) { if (i>j) swap(A[i],A[j]); for (int k=L>>1;(j^=k)<k;k>>=1); } for (int len=2;len<=L;len<<=1) { int l=len>>1; int W=mi(G,(mo-1)/len); for (int i=1;i<l;i++) w[i]=(ll)w[i-1]*W%mo; for (int i=0;i<L;i+=len) for (int j=0;j<l;j++) { int x=A[i+j],y=(ll)A[i+j+l]*w[j]%mo; A[i+j]=(x+y)%mo; A[i+j+l]=(x-y+mo)%mo; } } if (f==-1) { for (int i=1;i<L/2;i++) swap(A[i],A[L-i]); for (int i=0;i<L;i++) A[i]=(ll)A[i]*INV%mo; } } void pre(int l) { L=1; while (L<l) L<<=1; INV=mi(L,mo-2); w[0]=1; } void Get_INV(int *A,int n,int *B) { if (n==1) { B[0]=mi(A[0],mo-2);return;} Get_INV(A,n>>1,B); pre(n*2); static int C[N]; for (int i=0;i<L;i++) C[i]=0; for (int i=0;i<n;i++) C[i]=A[i]; NTT(C,1);NTT(B,1); for (int i=0;i<L;i++) B[i]=((B[i]+B[i])%mo-(ll)C[i]*B[i]%mo*B[i]%mo+mo)%mo; NTT(B,-1); for (int i=n;i<L;i++) B[i]=0; } void Calc_INV(int *A,int n,int *B) { ll X=1; while (X<=n) X<<=1; Get_INV(A,X,B); } void mul(int *A,int *B,int *C,int l) { pre(l); NTT(A,1);NTT(B,1); for (int i=0;i<L;i++) C[i]=(ll)A[i]*B[i]%mo; NTT(B,-1);NTT(C,-1); } void division(int *A,int *B,int n,int m,int *Ans) { reverse(A,A+n); reverse(B,B+m); static int C[N]; Calc_INV(B,n-m+1,C); mul(C,A,C,n+n); for (int i=n-m+1;i<L;i++) C[i]=0; reverse(C,C+n-m+1); pre(n+n); reverse(A,A+n); reverse(B,B+m); NTT(B,1);NTT(C,1); for (int i=0;i<L;i++) B[i]=(ll)B[i]*C[i]%mo; NTT(B,-1); for (int i=0;i<m-1;i++) Ans[i]=(A[i]-B[i]+mo)%mo; for (int i=0;i<L;i++) A[i]=B[i]=C[i]=0; } void mul(int *A,int *B) { static int a[N],b[N],Mod[N]; for (int i=0;i<K;i++) a[i]=A[i],b[i]=B[i],Mod[i]=MOD[i]; mul(a,b,a,K+K-1); division(a,Mod,K+K-1,K,A); } void solve(int *ans) { static int f[N]; ans[0]=1;f[1]=1; for (int i=1;i<=Len;i++) { if (T[i]) mul(ans,f); mul(f,f); } } } void add(int x,int y) { tot++; Next[tot]=head[x]; head[x]=tot; tree[tot]=y; } ll mi(ll x,int y) { ll ans=1; while (y) { if (y&1) ans=ans*x%mo; x=x*x%mo; y>>=1; } return ans; } void Add(int &x,int y) { x+=y;if (x>=mo) x-=mo;} vector<int> BM(vector<int> &A) { vector<int> B,S; int best=0,cnt=0; for (int i=0;i<(int)A.size();i++) { ll x=A[i]; for (int j=0;j<(int)B.size();j++) x=(x-(ll)B[j]*A[i-j-1]%mo+mo)%mo; delta[i]=x; if (!x) continue; fail[cnt++]=i; if (B.empty()) { B.resize(i+1); best=cnt-1; continue; } ll k=x*mi(delta[fail[best]],mo-2)%mo; vector<int> C; C.resize(i-fail[best]-1); C.push_back(k); for (int j=0;j<(int)S.size();j++) C.push_back((mo-k*S[j]%mo)%mo); if (C.size()<B.size()) C.resize(B.size()); for (int j=0;j<(int)B.size();j++) C[j]=(C[j]+B[j])%mo; if (B.size()<=i-fail[best]+S.size()) best=cnt-1,S=B; B=C; } return B; } void DP() { f[0][1]=1; int sum=1,now=0; for (int t=1;t<=n+n+5;t++) { now^=1; for (int i=1;i<=n;i++) f[now][i]=f[now^1][i]; for (int i=1;i<=n;i++) for (int j=head[i];j;j=Next[j]) { int v=tree[j]; Add(f[now][v],f[now^1][i]); } for (int i=1;i<=n;i++) Add(sum,f[now][i]); A.push_back(sum); } } void dec() { for (int i=1;i<=Len;i++) T[i]=s[i]-'0'; T[1]--; int now=1; while (T[now]<0) { T[now]+=2;T[now+1]--;now++;} } int calc() { int sum=0; for (int i=Len;i>=1;i--) sum=sum<<1|T[i]; return sum; } int main() { scanf("%d%d",&n,&m); for (int i=1;i<=m;i++) { int x,y; scanf("%d%d",&x,&y); add(x,y);add(y,x); } scanf("%s",s+1); Len=strlen(s+1); dec(); DP(); if (Len<=20&&calc()<=n+n+4) { printf("%d\n",A[calc()]);return 0;} B=BM(A); K=B.size(); for (int i=0;i<K;i++) MOD[i]=mo-B[K-1-i]; MOD[K++]=1; poly::solve(ans); int Ans=0; for (int i=0;i<K;i++) Add(Ans,(ll)ans[i]*A[i]%mo); printf("%d\n",Ans); return 0; } |
您好～我是腾讯云开发者社区运营，关注了您分享的技术文章，觉得内容很棒，我们诚挚邀请您加入腾讯云自媒体分享计划。完整福利和申请地址请见：https://cloud.tencent.com/developer/support-plan
作者申请此计划后将作者的文章进行搬迁同步到社区的专栏下，你只需要简单填写一下表单申请即可，我们会给作者提供包括流量、云服务器等，另外还有些周边礼物。