CF1106F Lunar New Year and a Recursive Sequence——矩阵快速幂&&bsgs

Posted 2022-10-17 lfri

题意

设 $$f_i = \\left\\\\beginmatrix
1 , \\ \\ \\ \\ \\ \\ \\ \\ \\ \\ \\ \\ \\ \\ \\ \\ \\ \\ \\ \\ \\ \\ \\ \\ \\  i < k\\\\
\\prod_j=1^k f_i-j^b_j \\ mod \\ p, \\ \\ \\ \\ \\ i > k
\\endmatrix\\right.$$

求 $f_k$（$1 \\leq f_k < p$），使得 $f_m = n$.（$1 \\leq k\\leq 100$）

分析

$f_n$ 可以表示成 $f_k^x$ 的形式，也就是指数的线性递推式，用矩阵快速幂求出最终 $f_n$ 中的次数就行了。

$$\\beginbmatrix f_k\\\\  f_k-1\\\\   \\vdots \\\\  f_1 \\endbmatrix =
\\beginbmatrix b_1 & b_2 & \\cdots  & b_k\\\\  1 & 0 & 0 & 0\\\\  \\vdots  & \\ddots & \\vdots  & \\vdots \\\\  0 & 0 & 1 & 0 \\endbmatrix \\cdot
\\beginbmatrix f_k-1\\\\  f_k-2\\\\   \\vdots \\\\  f_0 \\endbmatrix$$

即 $F_n = B\\cdot F_n-1 = B^n-kF_k$

那么就是 $f_k^x \\equiv f_n \\ (mod p) $ 形式了，其中 $x$ 是已经用矩阵快速幂算出来的。

于是就是关于形如 $x^a\\equiv b\\pmodp$ 方程的求解，直接用模板。

其中998244353的原根为3，算常识了吧。

注意算矩阵快速幂时，模并不是 $p$，由欧拉定理，模是 $p-1$.

#include<bits/stdc++.h>
using namespace std;

typedef long long ll;
struct matrix

    int r, c;
    ll mat[101][101];
    matrix()
        memset(mat, 0, sizeof(mat));
    
;
const ll p = 998244353;
int k, b[110], n, m;

matrix mul(matrix A, matrix B, ll p)   //矩阵相乘

    matrix ret;
    ret.r = A.r; ret.c = B.c;
    for(int i = 0;i < A.r;i++)
        for(int k = 0;k < A.c;k++)
            for(int j = 0;j < B.c;j++)
            
                ret.mat[i][j] = (ret.mat[i][j] + A.mat[i][k] * B.mat[k][j]) % p;
            
    return ret;


matrix mpow(matrix A, int n, int p)

    matrix ret;
    ret.r = A.r; ret.c = A.c;
    for(int i = 0;i < ret.r;i++)  ret.mat[i][i] = 1;
    while(n)
    
        if(n & 1)  ret = mul(ret, A, p);
        A = mul(A, A, p);
        n >>= 1;
    
    return ret;


ll gcd(ll a, ll b)

    return b ? gcd(b, a%b) : a;


ll qpow(ll a, ll b, ll p)

    a = a % p;
    ll ret = 1;
    while(b)
    
        if(b&1)  ret = ret * a % p;
        a = a * a %p;
        b >>= 1;
    
    return ret % p;




map<int,int>mp;
int bsgs(int a, int b, int p)    //a^x = b (mod P),(a,p)=1，返回x,x>=1
    int m=sqrt(p)+1;mp.clear();
    for(register int i=0,res=b;i<m;++i,res=1ll*res*a%p)mp[res]=i;
    for(register int i=1,tmp=qpow(a,m,p),res=tmp;i<=m+1;++i,res=1ll*res*tmp%p)
        if(mp.count(res))return i*m-mp[res];
    return -1;


int  main()

    scanf("%d", &k);
    for(int i = 1;i <= k;i++)  scanf("%d", &b[i]);
    scanf("%d%d", &n, &m);
    matrix B;
    B.r = B.c = k;
    for(int i = 0;i < k;i++)  B.mat[0][i] = b[i+1];
    for(int i = 1;i < k;i++)  B.mat[i][i-1] = 1;

    B = mpow(B, n-k, p-1);     
    int a = B.mat[0][0] % (p-1);    //注意，是模p-1 而非p

    int c = bsgs(qpow(3, a, p), m, p);
    if(c == -1)  printf("-1\\n");
    else
    
        int fk = qpow(3, c, p);
        printf("%d\\n", fk);
    

 

 

参考链接：https://www.cnblogs.com/bztMinamoto/p/10348641.html