El pueblo unido jamas serà vencido!

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

前言

昨晚与神 \(\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)\) 了.

然后是普通做法 :

\[\large\begin{aligned} \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\).

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

但是不够普适对吧.

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

\[\large\begin{aligned} \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\) 次剩余乘回去即可.

\[\large x^k \equiv b \pmod p \]

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

式子化为 :

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

然后代入值 :

\[\large 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;
}
posted @ 2022-01-30 19:56  AstatineAi  阅读(48)  评论(0编辑  收藏  举报