题解 P4726 【模板】多项式指数函数(多项式 exp)
Gauss0320
2019-05-23 17:41:29
先从牛顿迭代讲起。
已知多项式函数$G(z)$,求多项式函数$F(x)$满足
$$G(F(x))\equiv0 \pmod{x^n}$$
考虑用迭代求解,假设我们已经求得$F_0(x)$满足
$$G(F_0(x))\equiv0\pmod{x^{\left\lceil\frac{n}{2}\right\rceil}}$$
将函数$G$在$z=F_0(x)$处进行泰勒展开
$$G(F(x))=\sum_{i=1}^{\infty}\frac{G^i(F_0(x))}{i!}(F(x)-F_0(x))^i$$
其中$G^i$为$G$的$i$阶导函数.
取前两项
$$G(F(x))\equiv G(F_0(x))+G'(F_0(x))(F(x)-F_0(x))\pmod{x^n}$$
考虑到$G(F(x))\equiv 0\pmod{x^n}$
$$F(x)\equiv F_0(x)-\frac{G(F_0(x))}{G'(F_0(x))}\pmod{x^n}$$
边界条件即$f[0]=1$,向上迭代即可.
回到本题,考虑到
$$B(x)\equiv e^{A(x)}\pmod{x^n}$$
即
$$\ln B(x)-A(x)\equiv0\pmod{x^n}$$
于是令
$$G(B(x))\equiv\ln B(x)-A(x)\pmod{x^n}$$
由于$A(x)$为常数
$$G'(B(x))=B^{-1}(x)$$
套牛顿迭代
$$B(x)\equiv B_0(x)(1-\ln B_0(x)+A(x))\pmod{x^n}$$
这就很好求了,代码如下
```cpp
#include <bits/stdc++.h>
#define re register
using namespace std;
typedef long long ll;
const int N = 1<<20;
int read() {
int x = 0, f = 1; char c = getchar();
while(c < '0' || c > '9') {if(c == '-') f = -1; c = getchar();}
while(c >= '0' && c <= '9') x = (x<<1) + (x<<3) + c - '0', c = getchar();
return x * f;
}
inline void write(int a,char ed='\n')
{
static short s[13],tp;
if(!a){putchar('0'),putchar(ed);return;}
for(tp=0;a;a/=10)s[++tp]=a%10;
for(;tp;putchar(s[tp--]^48));
putchar(ed);
}
namespace Polynomial
{
const ll P = 998244353, g = 3, gi = 332748118;
static int rev[N];
int lim, bit;
int add(int a, int b)
{
return (a += b) >= P ? a - P : a;
}
int qpow(int a, int b)
{
int prod = 1;
while(b)
{
if(b & 1) prod = (ll)prod * a % P;
a = (ll)a * a % P;
b >>= 1;
}
return (prod + P) % P;
}
void calrev() {
for(re int i = 1; i < lim; ++i)
rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << (bit - 1));
}
void NTT(int *A, int inv)
{
for(re int i = 0; i < lim; ++i)
if(i < rev[i]) swap(A[i], A[rev[i]]);
for(re int mid = 1; mid < lim; mid <<= 1)
{
int tmp = qpow(inv == 1 ? g : gi, (P - 1) / (mid << 1));
for(re int j = 0; j < lim; j += (mid << 1))
{
int omega = 1;
for(re int k = 0; k < mid; ++k, omega = (ll)omega * tmp % P)
{
int x = A[j + k], y = (ll)omega * A[j + k + mid] % P;
A[j + k] = (x + y) % P;
A[j + k + mid] = (ll)(x - y + P) % P;
}
}
}
if(inv == 1) return;
int invn = qpow(lim, P - 2);
for(re int i = 0; i < lim; ++i)
A[i] = (ll)A[i] * invn % P;
}
static int x[N], y[N];
void mul(int *a, int *b)
{
memset(x, 0, sizeof x);
memset(y, 0, sizeof y);
for(re int i = 0; i < (lim >> 1); ++i)
x[i] = a[i], y[i] = b[i];
NTT(x, 1), NTT(y, 1);
for(re int i = 0; i < lim; ++i)
x[i] = (ll)x[i] * y[i] % P;
NTT(x, -1);
for(re int i = 0; i < lim; ++i)
a[i] = x[i];
}
static int c[2][N];
void Inv(int *a, int n)
{
int p = 0;
memset(c, 0, sizeof c);
c[0][0] = qpow(a[0], P - 2);
lim = 2, bit = 1;
while(lim <= (n << 1))
{
lim <<= 1, bit++;
calrev();
p ^= 1;
memset(c[p], 0, sizeof c[p]);
for(re int i = 0; i <= lim; ++i)
c[p][i] = add(c[p^1][i], c[p^1][i]);
mul(c[p^1], c[p^1]);
mul(c[p^1], a);
for(re int i = 0; i <= lim; ++i)
c[p][i] = add(c[p][i], P - c[p^1][i]);
}
for(re int i = 0; i < lim; ++i)
a[i] = c[p][i];
}
void derivative(int *a, int n)
{
for(re int i = 1; i <= n; ++i)
a[i - 1] = (ll)a[i] * i % P;
a[n] = 0;
}
void inter(int *a, int n)
{
for(re int i = n; i >= 1; --i)
a[i] = (ll)a[i - 1] * qpow(i, P - 2) % P;
a[0] = 0;
}
static int b[N], T[N], K[N];
void ln(int *a, int n)
{
memcpy(b, a, sizeof b);
Inv(b, n), derivative(a, n);
while(lim <= (n << 2)) lim <<= 1, bit++;
calrev();
mul(a, b);
inter(a, n);
for(re int i = n + 1; i <= lim; ++i)
a[i] = 0;
}
void exp(int *a, int n)
{
int z, d;
z = lim = 2, d = bit = 1;
memset(K, 0, sizeof K);
K[0] = 1;
while(z <= (n<<1))
{
z <<= 1, d++;
for(re int i = 0; i < (z>>1); ++i)
T[i] = K[i];
ln(T, (z>>1) - 1);
for(re int i = 0; i < (z>>1); ++i)
T[i] = add(a[i] + (i == 0), P - T[i]);
lim = z, bit = d;
calrev();
mul(K, T);
for(re int i = z; i <= (z<<1); ++i)
K[i] = T[i] = 0;
}
for(re int i = 0; i <= n; ++i)
a[i] = K[i];
}
}
using namespace Polynomial;
int n;
int F[N];
int main()
{
n = read();
for(re int i = 0; i < n; ++i)
F[i] = read();
exp(F, n);
for(re int i = 0; i < n; ++i)
write(F[i], ' ');
return 0;
}
```