[CF]codeforces1106F 矩阵乘法，BSGS+exgcd求K次剩余

一道不错的K次剩余（有原根版）板子题 1106F 题意简述：原本有一个f数列，满足∀i>k,$f_i = \Pi_{j=1}^k f_{i-j} ^ {b_j} (mod \ P)$,P=998244353，现在已知f1~fk−1=1，以及fn(n≥k)=m，求fk的可能值。 我们如果假设fk = a , 那么可以想到有m = a^u (u为某一个数)。对于某一个f[i] (i>k) 也为a^v(v为某一个数)，并且我们设定指数为g,那么g[i] = \sum g[i-k] * b[k]。 那么这是一个常系数线性递推，矩阵快速幂解决。 也就是说我们知道fk^u = m ，已知u,m，求fk。这是经典的k次剩余。我们转化为原根的形式，那么有g^(x * u) = g^r (mod p) , r我们通过g^r = m ，BSGS求出，然后我们就是要解x,那么又有 x * u = r (mod p-1)。这个通过扩展欧基里德求逆元（可能无解），就得到x,之后fk = g^k (g 表示998244353的原根3) code:

#include<stdio.h>
#include<cstdio>
#include<algorithm>
#include<iostream>
#include<cmath>
#include<map>
#define int long long
using namespace std;
const int mod = 998244353;
int add(int x,int y) { x+=y; return x>=mod?x-mod:x; }
int sub(int x,int y) { x-=y; return x<0?x+mod:x; }
int mul(int x,int y) { return 1ll*x*y%mod; }
int ksm(int a,int b) {
    int ans = 1;
    for(;b;b>>=1,a=mul(a,a))
        if(b&1) ans=mul(ans,a);
    return ans;
}
int k;
struct mat {
    int o[101][101];
}g,CS;
mat operator*(mat aa,mat bb) {
    static mat tmp;
    for(int i=1;i<=k;i++) {
        for(int j=1;j<=k;j++) {
            int oo = 0;
            for(int o=1;o<=k;o++) {
                oo = ( oo+1ll*aa.o[i][o]*bb.o[o][j]%(mod-1) )%(mod-1);
            }
            tmp.o[i][j] = oo;
        }
    }
    return tmp;
}
int jzksm(int b) {
    for(;b;b>>=1,g=g*g)
        if(b&1) CS = CS*g;
    return CS.o[1][1];
}
int n,m;
int gcd(int a,int b) {
    return !b?a:gcd(b,a%b);
}
void exgcd(int &x,int &y,int a,int b) {
    if(!b) { x=1;y=0; return; }
    exgcd(x,y,b,a%b);
    int tx = x; x = y; y = tx - (a/b)*y;
}
map<int,int>ma;
int BSGS(int x) {
    int sq = ceil(sqrt(mod));
    int Ai = 1; int Aj = 1;
    ma[x] = 0;
    for(int i=1;i<sq;i++) {
        Ai = mul(Ai,3);
        ma[ mul(Ai,x) ] = i;
    }
    Ai = mul(Ai,3);
    for(int i=1;i<=sq;i++) {
        Aj = mul(Aj,Ai);
        if(ma.count(Aj)) return i*sq-ma[Aj];
    }
}
main() {
    scanf("%I64d",&k);
    for(int i=1;i<=k;i++) {
        int x; scanf("%I64d",&x);
        g.o[i][1] = x;
    }
    for(int j=2;j<=k;j++) {
        g.o[j-1][j] = 1;
    }
    scanf("%I64d%I64d",&n,&m);
    CS.o[1][1] = 1;
    int p = jzksm(n-k);
    m = BSGS(m);
    cerr<<p<<' '<<m<<endl;
    if(m%gcd(p,mod-1)) puts("-1");
    else {
        int dd = gcd(p,mod-1);
        int x,y, a = p , b = mod - 1;
        m /= dd; a/=dd; b/=dd;
        exgcd(x,y,a,b);
        x = (x+mod-1)%(mod-1);
        x = x*m%(mod-1);
        printf("%I64d",ksm(3,x));
    }
}