CF1106F Lunar New Year and a Recursive Sequence 矩阵加速+BSGS+原根+扩展欧几里得

文章目录

• 一.题目
• 二.Solution
• 三.Code
• Thanks!

一.题目

有一串 n ( n ⩽ 1 0 9 ) n(n\leqslant10^9) n(n⩽109)个数的数列，给你 b 1 ∼ b k ( k ⩽ 100 ) b_1\sim b_k(k\leqslant100) b1​∼bk​(k⩽100)当 i > k i > k i>ki>k i>ki>k时： f i = ( ∏ j = 1 k f i − j b j )   m o d   998244353 f_i=(\prod\limits_{j=1}^kf_{i-j}^{b_j})\bmod{998244353} fi​=(j=1∏k​fi−jbj​​)mod998244353已知 f 1 = f 2 = ⋯ = f k − 1 = 1 , f n = m f_1=f_2=\cdots=f_{k-1}=1,f_n=m f1​=f2​=⋯=fk−1​=1,fn​=m，问最小的正整数 f k f_k fk​可能是多少?

二.Solution

这是一道反正我考场上推不出来的的数学题。

首先因为 f 1 ∼ f k − 1 = 1 f_1\sim f_{k-1}=1 f1​∼fk−1​=1并且递推式是累乘，所以 f n f_n fn​必定可以表示成 f n = f k p f_n=f_k^p fn​=fkp​的形式，然后 p p p又满足递推式： g [ i ] = ∏ j = 1 k b [ j ] ∗ g [ i − j ] g[i]=\prod^{k}_{j=1}b[j]*g[i-j] g[i]=∏j=1k​b[j]∗g[i−j]所以可以用矩阵加速来计算 p p p的大小，构造这样一个递推矩阵 B B B：

0 0 . . . 0 b [ k ] 1 0 . . . 0 b [ k − 1 ] 0 1 . . . 0 b [ k − 2 ] . . . . . . . . . . . . 0 0 . . . 1 b [ 1 ] \begin{matrix} 0 & 0 & ... & 0 & b[k] \\ 1 & 0 & ... & 0 & b[k - 1] \\ 0 & 1 & ... & 0 & b[k - 2] \\ ... &... &... &...\\ 0 & 0 & ... & 1 & b[1] \end{matrix} 010...0​001...0​...............​000...1​b[k]b[k−1]b[k−2]b[1]​

这样一个答案矩阵 A A A：

0 0 . . . 1 \begin{matrix} 0 & 0 & ... & 1 \end{matrix} 0​0​...​1​

最后快速幂完成后 p p p就是 A . j z [ 1 ] [ k ] A.jz[1][k] A.jz[1][k]

那么问题就来到 f k p ≡ m ( m o d   998244353 ) f_{k}^{p}\equiv m (mod\,998244353) fkp​≡m(mod998244353)求 f k f_k fk​。

这里可以用原根做，我们知道 998244353 998244353 998244353的原根是 3 3 3，那么设 f k ≡ 3 s ( m o d   998244353 ) f_k\equiv 3^s(mod\,998244353) fk​≡3s(mod998244353)，设 m ≡ 3 t ( m o d   998244353 ) m\equiv 3^t(mod\,998244353) m≡3t(mod998244353)，因为 m m m已知，所以可以用BSGS求出 t t t，最后又化成： 3 s p ≡ 3 t ( m o d   998244353 ) 3^{sp}\equiv 3^t (mod\,998244353) 3sp≡3t(mod998244353)又变成 s p ≡ t ( m o d   998244352 ) sp\equiv t(mod\,998244352) sp≡t(mod998244352)

这里用扩展欧几里得定理求出 s s s即可，答案就是 3 s p 3^{sp} 3sp

三.Code

#include <cstdio>
#include <cstring>
#include <iostream>
#include <map>
#include <cmath>
using namespace std;

#define M 105
#define LL long long
const LL mod = 998244353;

int k, n;
LL b[M], m;
map <LL, LL> mp;
struct matrix{
	LL jz[M][M];
	matrix clear (int flag){
		if (! flag){
			memset (jz, 0, sizeof jz);
		}
		else{
			for (int i = 1; i <= k; i ++)
				for (int j = 1; j <= k; j ++)
					jz[i][j] = (i == j) ? 1 : 0;
		}
	}
	matrix operator * (const matrix& rhs) const{
		matrix fina;
		fina.clear(0);
		for (int i = 1; i <= k; i ++)
			for (int j = 1; j <= k; j ++)
				for (int z = 1; z <= k; z ++)
					fina.jz[i][j] = (fina.jz[i][j] + jz[i][z] * rhs.jz[z][j] % (mod - 1)) % (mod - 1);
		return fina;
	}
}A, B;

matrix qkpow1 (matrix x, LL y){
	matrix fina;
	fina.clear(1);
	while (y){
		if (y & 1)
			fina = fina * x;
		x = x * x;
		y >>= 1;
	}
	return fina;
}
LL qkpow2 (LL x, LL y){
	LL res = 1;
	while (y){
		if (y & 1)
			res = res * x % mod;
		x = x * x % mod;
		y >>= 1;
	}
	return res;
}
LL bsgs (LL a, LL b){
	mp.clear();
	mp[1] = 0;
	LL sqt = ceil (sqrt (mod - 1)), t1 = 1, t2 = 1;
	for (int i = 1; i < sqt; i ++){
		t1 = t1 * a % mod;
		mp[t1 * b % mod] = 1ll * i;
	}
	t1 = t1 * a % mod;
	for (int i = 1; i <= sqt; i ++){
		t2 = t1 * t2 % mod;
		if (mp[t2])
			return (1ll * i * sqt - mp[t2]);
	}
}
LL get_gcd (LL x, LL y){
	if (! y)
		return x;
	return get_gcd (y, x % y);
}
void exgcd (LL &x, LL &y, LL a, LL b){
	if (! b){
		x = 1, y = 0;
		return ;
	}
	exgcd (y, x, b, a % b);
	y -= a / b * x;
}
int main (){
	scanf ("%d", &k);
	for (int i = 1; i <= k; i ++)
		scanf ("%lld", &b[i]);
	scanf ("%d %lld", &n, &m);
	A.jz[1][k] = 1;
	for (int i = 1; i <= k; i ++)
		B.jz[i][k] = b[k - i + 1];
	for (int i = 1; i < k; i ++)
		B.jz[i + 1][i] = 1;
	A = A * qkpow1 (B, n - k);
	LL p = A.jz[1][k], t = bsgs (3, m);
	if (t % get_gcd (p, mod - 1))
		printf ("-1\n");
	else{
		LL gcd = get_gcd (p, mod - 1);
		LL x, y, b = mod - 1;
		p /= gcd, t /= gcd, b /= gcd;
		exgcd (x, y, p, b);
		x = (x % (mod - 1) + mod - 1) % (mod - 1);
		x = x * t % (mod - 1);
		printf ("%lld\n", qkpow2 (3, x));
	}
	return 0;
}

           

Thanks!