上一篇:
一、对数函数\(\ln(A)\)
求一个多项式\(B(x)\),满足\(B(x)=\ln(A(x))\)。
这里需要一些最基本的微积分知识(不会?戳我(暂时戳不动):【知识总结】微积分初步挖坑待填)。
另外,\(n\)次多项式\(A(x)\)可以看成关于\(x\)的\(n\)次函数,可以对其求导。显然,\(A(x)=\sum\limits_{i=0}^{n-1}a_ix^i\)的导数是\(A'(x)=\sum\limits_{i=0}^{n-2}a_{i+1}x^i(i+1)\),积分是\(\int A(x)\mathrm{d} x=\sum\limits_{i=1}^{n}\frac{a_{i-1}}{i}x^i\)。可以写出如下代码(非常简单):
void derivative(const int *a, int *b, const int n){ for (int i = 1; i < n; i++) b[i - 1] = (ll)a[i] * i % p; b[n - 1] = 0;}void integral(const int *a, int *b, const int n){ for (int i = n - 1; i >= 0; i--) b[i + 1] = (ll)a[i] * inv(i + 1) % p; b[0] = 0;}
\(f(x)=\ln(x)\)的导数是\(f'(x)=\frac{1}{x}\)。回到原问题,对两边同时求导,得到(要用一下链式法则\(g(f(x))\)的导数是\(g'(f(x))f'(x)\)):
\[B'(x)=A'(x)\frac{1}{A(x)}=\frac{A'(x)}{A(x)}\]
求个\(A(x)\)的逆元(多项式求逆)然后乘上\(A'(x)\),最后把\(B(x)\)积分回去就好了。
至于代码……下面算多项式指数函数的时候要算对数函数,所以暂时省略。
二、指数函数\(\exp(x)\)
求多项式\(B(x)\)满足\(B(x)=e^{A(x)}\)。
首先,这个式子相当于求\(\ln B(x)=A(x)\)即\(\ln B(x)-A(x)=0\)
设关于多项式的函数\(F(B(x))=\ln B(x)-A(x)\),那么问题就是求这个函数的零点(\(A(x)\)是给定的,视作常数)。
求函数零点的方法之一是牛顿迭代,公式如下(\(i\)是迭代次数,\(x\)是自变量,\(F(x)\)是要求零点的函数,\(F'(x_0)\)是\(F(x)\)在\(x_0\)处的导数):
\[x_{i+1}=x_i-\frac{F(x_i)}{F'(x_i)}\]
把\(F(B(x))=\ln B(x)-A(x)\)求导,得到\(F'(B(x))=\frac{1}{B(x)}\)(注意自变量是\(B(x)\)不是\(x\)。这不是一个\(F(x)\)和\(B(x)\)的复合函数)。然后代入上面的公式:
\[ \begin{aligned} B_{i+1}(x)&=B_i(x)-\frac{\ln B_i(x)-A(x)}{\frac{1}{B_i(x)}}\\ &=B_i(x)-B_i(\ln B_i(x)-A(x))\\ &=B_i(x)(1-\ln B_i(x)-A(x)) \end{aligned} \]
由于多项式乘法的存在,每迭代一次\(B\)的有效长度会增加一倍。
下一篇:
代码():
#include#include #include #include #undef i#undef j#undef k#undef true#undef false#undef min#undef max#undef swap#undef sort#undef if#undef for#undef while#undef printf#undef scanf#undef putchar#undef getchar#define _ 0using namespace std;namespace zyt{ template inline bool read(T &x) { char c; bool f = false; x = 0; do c = getchar(); while (c != EOF && c != '-' && !isdigit(c)); if (c == EOF) return false; if (c == '-') f = true, c = getchar(); do x = x * 10 + c - '0', c = getchar(); while (isdigit(c)); if (f) x = -x; return true; } template inline void write(T x) { static char buf[20]; char *pos = buf; if (x < 0) putchar('-'), x = -x; do *pos++ = x % 10 + '0'; while (x /= 10); while (pos > buf) putchar(*--pos); } typedef long long ll; const int N = 1e5 + 10, LEN = (N << 2), p = 998244353, g = 3; namespace Polynomial { inline int power(int a, int b) { int ans = 1; while (b) { if (b & 1) ans = (ll)ans * a % p; a = (ll)a * a % p; b >>= 1; } return ans; } inline int inv(const int a) { return power(a, p - 2); } int omega[LEN], winv[LEN], rev[LEN]; void init(const int n, const int lg2) { int w = power(g, (p - 1) / n), wi = inv(w); omega[0] = winv[0] = 1; for (int i = 1; i < n; i++) { omega[i] = (ll)omega[i - 1] * w % p; winv[i] = (ll)winv[i - 1] * wi % p; } for (int i = 0; i < n; i++) rev[i] = ((rev[i >> 1] >> 1) | ((i & 1) << (lg2 - 1))); } void ntt(int *a, const int *w, const int n) { for (int i = 0; i < n; i++) if (i < rev[i]) swap(a[i], a[rev[i]]); for (int l = 1; l < n; l <<= 1) for (int i = 0; i < n; i += (l << 1)) for (int k = 0; k < l; k++) { int tmp = (a[i + k] - (ll)w[n / (l << 1) * k] * a[i + l + k] % p + p) % p; a[i + k] = (a[i + k] + (ll)w[n / (l << 1) * k] * a[i + l + k] % p) % p; a[i + l + k] = tmp; } } void mul(const int *a, const int *b, int *c, const int n) { static int x[LEN], y[LEN]; int m = 1, lg2 = 0; while (m < (n << 1) - 1) m <<= 1, ++lg2; init(m, lg2); memcpy(x, a, sizeof(int[n])); memset(x + n, 0, sizeof(int[m - n])); memcpy(y, b, sizeof(int[n])); memset(y + n, 0, sizeof(int[m - n])); ntt(x, omega, m), ntt(y, omega, m); for (int i = 0; i < m; i++) x[i] = (ll)x[i] * y[i] % p; ntt(x, winv, m); int invm = inv(m); for (int i = 0; i < m; i++) x[i] = (ll)x[i] * invm % p; memcpy(c, x, sizeof(int[n])); } void _inv(const int *a, int *b, const int n) { if (n == 1) b[0] = inv(a[0]); else { static int tmp[LEN]; _inv(a, b, (n + 1) >> 1); int m = 1, lg2 = 0; while (m < (n << 1) + 1) m <<= 1, ++lg2; init(m, lg2); memcpy(tmp, a, sizeof(int[n])); memset(tmp + n, 0, sizeof(int[m - n])); memset(b + ((n + 1) >> 1), 0, sizeof(int[m - ((n + 1) >> 1)])); ntt(tmp, omega, m); ntt(b, omega, m); for (int i = 0; i < m; i++) b[i] = (b[i] * 2LL % p - (ll)tmp[i] * b[i] % p * b[i] % p + p) % p; ntt(b, winv, m); int invm = inv(m); for (int i = 0; i < m; i++) b[i] = (ll)b[i] * invm % p; memset(b + n, 0, sizeof(int[m - n])); } } void inv(const int *a, int *b, const int n) { static int tmp[LEN]; memcpy(tmp, a, sizeof(int[n])); _inv(tmp, b, n); } void derivative(const int *a, int *b, const int n) { for (int i = 1; i < n; i++) b[i - 1] = (ll)a[i] * i % p; b[n - 1] = 0; } void integral(const int *a, int *b, const int n) { for (int i = n - 1; i >= 0; i--) b[i + 1] = (ll)a[i] * inv(i + 1) % p; b[0] = 0; } void ln(const int *a, int *b, const int n) { static int tmp[LEN], inva[LEN]; derivative(a, tmp, n); inv(a, inva, n - 1); mul(inva, tmp, b, n - 1); integral(b, b, n - 1); } void _exp(const int *a, int *b, const int n) { if (n == 1) b[0] = 1; else { static int tmp[LEN]; _exp(a, b, (n + 1) >> 1); ln(b, tmp, n); for (int i = 0; i < n; i++) tmp[i] = (-tmp[i] + a[i] + p) % p; tmp[0] = (tmp[0] + 1) % p; mul(b, tmp, b, n); } } void exp(const int *a, int *b, const int n) { static int tmp[LEN]; memcpy(tmp, a, sizeof(int[n])); _exp(tmp, b, n); } } int work() { static int a[LEN]; int n; read(n); for (int i = 0; i < n; i++) read(a[i]); Polynomial::exp(a, a, n); for (int i = 0; i < n; i++) write(a[i]), putchar(' '); return (0^_^0); }}int main(){ return zyt::work();}