CF1106F Lunar New Year and a Recursive Sequence - 矩阵快速幂、数论

题解

第一次不看 editorial 做出来 CF 难度 \(2400\) 的题！！！

这个推 \(f_i\) 的式子一看就很像矩阵加速，但它里面是乘号。于是我们考虑维护每个 \(f_i\) 里面含有多少个 \(f_1,f_2,\dots,f_k\)。

记 \(f_{i,j}\) 为 \(f_i\) 里面含有的 \(f_j(1\le j\le k)\) 个数。于是可以构造转移矩阵：

\[\begin{bmatrix} f_{n,1} & f_{n,2} & \dots & f_{n,k} \\ f_{n+1,1} & f_{n+1,2} & \dots & f_{n+1,k} \\ \vdots & \vdots & \vdots & \vdots \\ f_{n+k-1,1} & f_{n+k-1,2} & \dots & f_{n+k-1,k} \end{bmatrix}=\begin{bmatrix} 0 & 1 & 0 & 0 & \dots & 0 \\ 0 & 0 & 1 & 0 & \dots & 0 \\ \vdots & \vdots & \vdots &\vdots &\vdots &\vdots \\ b_k & b_{k-1} & b_{k-2} & b_{k-3} & \dots & b_1 \end{bmatrix} \times \begin{bmatrix} f_{n-1,1} & f_{n-1,2} & \dots & f_{n-1,k} \\ f_{n,1} & f_{n,2} & \dots & f_{n,k} \\ \vdots & \vdots & \vdots & \vdots \\ f_{n+k-2,1} & f_{n+k-2,2} & \dots & f_{n+k-2,k} \end{bmatrix} \]

矩阵快速幂完了，取第 \(1\) 行第 \(k\) 列那个值（记为 \(a\)），那么用原根和 BSGS 解同余方程 \(x^a\equiv m\pmod{998,244,353}\) 即可。

代码出人意料地短。

代码

#include <cstdio>
#include <cstring>
#include <cctype>
#include <map>
#include <cmath>
using namespace std;
#define For(Ti,Ta,Tb) for(int Ti=(Ta);Ti<=(Tb);++Ti)
#define Dec(Ti,Ta,Tb) for(int Ti=(Ta);Ti>=(Tb);--Ti)
template<typename T> void Read(T &x){
	x=0;int _f=1;
	char ch=getchar();
	while(!isdigit(ch)) _f=(ch=='-'?-1:_f),ch=getchar();
	while(isdigit(ch)) x=x*10+(ch^48),ch=getchar();
	x=x*_f;
}
template<typename T,typename... Args> void Read(T &x,Args& ...others){
	Read(x);Read(others...);
}
typedef long long ll;
const int K=105,Mod=998244353,Root=3;
typedef ll Mat[K][K];
int len;
void MatMul(Mat &c,const Mat &a,const Mat &b){
	Mat res;memset(res,0,sizeof res);
	For(k,1,len) For(i,1,len) For(j,1,len){
		res[i][j]=(res[i][j]+a[i][k]*b[k][j])%(Mod-1);
	}memcpy(c,res,sizeof res);
}
void MatPow(Mat &c,const Mat &x,ll b){
	Mat res,a;
	memset(res,0,sizeof res);memcpy(a,x,sizeof x);
	For(i,1,len) res[i][i]=1;
	while(b){
		if(b&1) MatMul(res,res,a);
		b>>=1,MatMul(a,a,a);
	}memcpy(c,res,sizeof res);
}
ll Pow(ll a,ll b,ll p){
	ll res=1;
	while(b){
		if(b&1) res=res*a%p;
		b>>=1,a=a*a%p;
	}return res;
}
ll BSGS(ll a,ll b,ll p){//a^x=b(mod p)
	map<ll,int> mp;
	ll t=ceil(sqrt(p+.5));
	for(ll cur=1,x=0;x<=t;++x,cur=cur*a%p){
		mp[cur*b%p]=x;
	}
	ll temp=Pow(a,t,p);
	for(ll cur=temp,x=1;x<=t;++x,cur=cur*temp%p){
		if(mp.count(cur)) return x*t-mp[cur];
	}return -1;
}
int n;ll b[K],m;Mat mat;
int main(){
	Read(len);For(i,1,len) Read(b[i]);
	Read(n,m);
	For(i,1,len-1) mat[i][i+1]=1;
	For(j,1,len) mat[len][j]=b[len-j+1];
	MatPow(mat,mat,n-1);
	ll a=mat[1][len];
	ll res=BSGS(Pow(Root,a,Mod),m,Mod);
	if(res==-1) return puts("-1"),0;
	printf("%lld\n",Pow(Root,res,Mod));
	return 0;
}