LOJ#6713. 「EC Final 2019」狄利克雷 k 次根 加强版

题目描述

定义两个函数 (f, g: {1, 2, dots, n} ightarrow mathbb Z) 的狄利克雷卷积 (f * g) 为：

[(f * g)(n) = sum_{d | n} f(d)g(frac nd) ]

我们定义 (g = f^k) 即 (k) 次幂为：

[f^{k}=underbrace {f * dots * f} _{k~{ extrm {个}}} ]

在本题中，我们想要解决这个问题的逆问题：给你 (g) 和 (k)，你需要找到一个函数 (f) 使得 (g = f^k)。

另外，保证 (g(1)=1)，你需要保证 (f(1)=1)。所有的运算在 (mathbb F_p) 上进行，其中 (p = 998244353)，这意味着狄利克雷卷积为 ((f*g)(n) = left(sum_{d|n} f(d)g(frac nd) ight) mod p)。

(nle 10^6)

Sol

定义一个数论函数 (f(n)) 的狄利克雷生成函数是 (F(x)=sum_{n=1}^infty f(n)n^x)。

注意这个东西不是关于 (x) 的多项式（一开始我当成了多项式），然而它可以做多项式乘法，由于这个 (a^xcdot b^x ightarrow (ab)^x)，所以对于生成函数系数来说，相当于完成了两个数论函数的 狄利克雷卷积。

既然有类似多项式的性质，那么我们考虑 ln-exp 的方法求解 k 次根。

根据多项式 ln 的推导，有 (B=ln A ightarrow B=int {A'over A})

除法可以用 狄利克雷除法 实现，那么我们只需要支持 导数和积分 运算即可。

考虑这个指数函数的求导法则 (displaystyle { ext{d} a^xover ext{d}x}=a^xln a)。

然而在正整数意义下没有 (ln) 这种运算，然后做法是找一个运算性质跟 (ln) 一样的东西来替代它，比如质因子个数函数 ( ext{cnt}_x)（相同的算作多个），因为它满足 (f(x)+f(y)=f(xy)) 这样的性质。

于是可以实现一个对一个狄利克雷级数求对数函数的代码如下：

void ln(int n, int *a)
{
	b[1] = 0;
	for(int i = 2; i <= n; ++i)
		b[i] = (ll)a[i] * cnt[i] % P;
	for(int i = 2; i <= n; ++i)
	{
		for(int j = 2, k = i * j; k <= n; ++j, k += i)
			b[k] = (b[k] - 1LL * b[i] * a[j]) % P;
		b[i] = (1LL * b[i] * inv[cnt[i]] % P + P) % P;
	}
	for(int i = 1; i <= n; ++i)
		a[i] = b[i];
}

由于 (B=A^{1/k}=exp {{ln A} over k})，那么只要先求 (ln) 再 (exp) 即可。
求指数函数也是类似的，整个代码如下，和暴力 (O(n^2)) 的普通多项式对数/指数函数非常相似。

#include<stdio.h>

typedef long long ll;
const int N = 1000005, P = 998244353;

int fexp(int a, int b)
{
	ll x = 1, o = a;
	for(; b; b >>= 1, o = o * o % P)
		if(b & 1) x = x * o % P;
	return x;
}

int n, K, a[N], b[N], inv[128], cnt[N];
bool np[N];

void init()
{
	cnt[1] = 1;
	for(int i = 2; i <= n; ++i)
		if(!np[i])
		{
			cnt[i] = 1;
			for(int j = 2, k = i * j; k <= n; j++, k += i)
			{
				np[k] = true;
				if(!cnt[k] && cnt[j]) cnt[k] = cnt[j] + 1;
			}
		}
	inv[1] = 1;
	for(int i = 2; i < 128; ++i)
		inv[i] = (ll)inv[P % i] * (P - P / i) % P;
}

void ln(int n, int *a)
{
	b[1] = 0;
	for(int i = 2; i <= n; ++i)
		b[i] = (ll)a[i] * cnt[i] % P;
	for(int i = 2; i <= n; ++i)
	{
		for(int j = 2, k = i * j; k <= n; ++j, k += i)
			b[k] = (b[k] - 1LL * b[i] * a[j]) % P;
		b[i] = (1LL * b[i] * inv[cnt[i]] % P + P) % P;
	}
	for(int i = 1; i <= n; ++i)
		a[i] = b[i];
}

void exp(int n, int *a)
{
	for(int i = 2; i <= n; ++i)
		b[i] = (ll)a[i] * cnt[i] % P, a[i] = 0;
	a[1] = 1;
	for(int i = 1; i <= n; ++i)
	{
		a[i] = 1ll * a[i] * inv[cnt[i]] % P;
		for(int j = 2, k = i * j; k <= n; ++j, k += i)
			a[k] = (a[k] + 1LL * a[i] * b[j]) % P;
	}
}

int main()
{
	scanf("%d %d", &n, &K);
	for(int i = 1; i <= n; ++i)
		scanf("%d", a + i);
	init();
	K = fexp(K, P - 2);
	ln(n, a);
	for(int i = 1; i <= n; ++i) a[i] = 1ll * a[i] * K % P;
	exp(n, a);
	for(int i = 1; i <= n; ++i) printf("%d%c", a[i], " 
"[i == n]);
	return 0;
}