[可能有用科技]多项式任意次幂&任意次根

望月Asta

2022-01-11 17:33:21

Personal

前言

昨晚与神 \mathrm{{\color{black}w}{\color{red}ind\_whisper}} 讨论到我做丢失的题面 #6 用的是暴力多项式三次根而不是分析数据特殊性然后使用组合数.

于是就聊到了这种方式是否可以应用到开任意 k 次根,遂有此博客.

多项式 k 次幂

模板 : \texttt{Link.}

对于一个多项式 F(x), 求 G(x) \equiv (F(x))^k \bmod x^n.

首先看一个 naive 的做法 : 参考普通的快速幂,即可做到 \mathcal{O}(n\log n \log k) 了.

然后是普通做法 :

\ln (F(x))^k &= kF(x)\\ (F(x))^k &\equiv \exp (k \ln F(x)) \bmod x^n \end{aligned}

但是这个做法有一个问题,那就是因为自然对数的底 e 没有模意义下的值,所以多项式 \ln 必须保证 F(0) = 1 才行.

那么既然已经有模数了,就让多项式的每一项都乘以 (F(0))^{-1} \bmod {998244353} 即可.

但是众所周知 0 没有逆元, F(0) = 0 就寄了.

那么就把前面的 0 全部吞掉,最后前面把吞了的 0 输出一下即可.

写完发现其实把模数作为一个 char 型指针传进去比较文明,遂参考了一下学长的实现.

inline void PolyPow(const int *a,int *b,int n,const char *k,int m = -1) {
    int m1 = 0,m2 = 0,m3 = 0;
    ll t = 0;while(!a[t] && t < n) ++t;
    if(~m) m1 = m2 = m3 = m;
    else {
        int l = strlen(k);
        repl(i,0,l) {
            m1 = ((ll)m1 * 10 + (k[i] ^ 48)) % MOD;
            m2 = ((ll)m2 * 10 + (k[i] ^ 48)) % (MOD - 1);
            if(((ll)m3 * 10 + (k[i] ^ 48)) <= MOD)
                m3 = (ll)m3 * 10 + (k[i] ^ 48);
        }
        if((ll)m3 * t >= n) {
            repl(i,0,n) b[i] = 0;
            return ;
        }
    }
    static int x[N],c[N];
    repl(i,0,n) x[i] = a[i + t];
    t = t * m1;
    const int res = qpow(x[0],m2),invx = qpow(x[0]);
    repl(i,0,n) x[i] = (ll)x[i] * invx % MOD;
    PolyLn(x,b,n);repl(i,0,n) b[i] = (ll)b[i] * m1 % MOD;
    PolyExp(b,c,n);repl(i,0,n) c[i] = (ll)c[i] * res % MOD;
    repl(i,0,t) b[i] = 0;
    repl(i,t,n) b[i] = c[(ll)i - t];
}

多项式 k 次根

模板 : \texttt{Link.}

对于一个多项式 F(x), 求 (G(x))^k \equiv F(x) \bmod x^n.

众所周知这是可以推式子然后直接牛顿迭代的.

但是不够普适对吧.

那么现在给出一个新(旧)方法.

\ln \sqrt[k]{F(x)} &= \frac{F(x)}{k}\\ \sqrt[k]{F(x)} &\equiv \exp (k^{-1} \ln F(x)) \bmod x^n \end{aligned}

无敌的对数函数

发现第一项不是 1, \ln\exp 是错的.

于是先全部项乘上 F(0) 的逆元,\exp 完后直接求 F(0)k 次剩余乘回去即可.

x^k \equiv b \pmod p

但是现在模数是 998244353,原根为 3 ,这代表着可以使用 3^r 表示一个数.

式子化为 :

(3^{r})^k \equiv b \pmod p

然后代入值 :

3^{rk} \equiv F(0) \pmod{998244353}

BSGS 求解即可.

那么问题来了,如何判断一个数有没有 k 次剩余?

在模 998244353 意义下,先定义 3^x \equiv F(0) \pmod{998244353}

那么仅当 \frac{x}{k} 在模 998244352 有意义才有 k 次剩余.

啥?有多解?

确实可能有多解.

求最小的就暴力枚举呗.

constexpr int SIZ = 100003;
struct HashTable {
    int tot;
    int head[SIZ];
    struct Node {
        int nxt,val;
        ll key;
    }p[N];

    HashTable() : tot(0) {mems(head,0);}

    inline int& operator [] (const int x) {
        int h = x % SIZ;
        for(int i = head[h];i;i = p[i].nxt)
            if(p[i].key == x) return p[i].val;
        p[++tot] = (Node) {head[h],0,x},head[h] = tot;
        return p[tot].val;
    }
}mp;

inline int BSGS(int a,int b) {
    const int mx = ceil(sqrtf(MOD));
    int base = b % MOD;
    rep(i,1,mx) {
        base = (ll)base * a % MOD;
        mp[base] = i;
    }
    base = qpow(a,mx);
    if(a == 0) return b == 0 ? 1 : -1;
    int prod = 1;
    rep(i,1,mx) {
        prod = (ll)prod * base % MOD;
        int px = mp[prod];
        if(px) return (((ll)i * mx - px) % MOD + MOD) % MOD;
    }
    return -1;
}

void exgcd(int a,int b,int& d,int& x,int& y) {
    if(!b) {
        d = a; x = 1; y = 0;
        return ;
    }
    exgcd(b,a % b,d,y,x);
    y -= x * (a / b);
}

int GetInv(int a,int p) {
    int d, x, y;
    exgcd(a, p, d, x, y);
    return d == 1 ? (x + p) % p : -1;
}

inline int gcd(int a,int b) {
    int k = __builtin_ctz(a | b);
    a >>= __builtin_ctz(a);
    b >>= __builtin_ctz(b);
    int p = a % b;
    while(p)
        a = b,b = p,p = a % b;
    return b << k;
}

inline void PolyResidue(const int *a,int *b,int n,int k) {
    static int x[N],c[N];int invk = qpow(k),invc = qpow(a[0]);
    repl(i,0,n) c[i] = (ll)a[i] * invc % MOD;
    PolyLn(c,x,n);
    repl(i,0,n) x[i] = (ll)x[i] * invk % MOD;
    PolyExp(x,b,n);
    int pw = BSGS(G,a[0]),w = INF;
    int M = MOD - 1,g = gcd(M,gcd(k,pw));
    k /= g,pw /= g,M /= g;
    int t = (ll)pw * GetInv(k,M) % M;
    for(;t < MOD - 1;t += M)
        w = std::min(w,qpow(3,t));
    repl(i,0,n) b[i] = (ll)b[i] * w % MOD;
}