多项式初等函数

本页面包含多项式常见的初等函数操作.具体而言,本页面包含如下内容:

  1. 多项式求逆
  2. 多项式开方
  3. 多项式除法
  4. 多项式取模
  5. 多项式指数函数
  6. 多项式对数函数
  7. 多项式三角函数
  8. 多项式反三角函数
初等函数与非初等函数

初等函数的定义如下1:

若域 FF 中存在映射 u→∂uu\to \partial u 满足:

  1. ∂(u+v)=∂u+∂v\partial(u+v)=\partial u+\partial v
  2. ∂(uv)=u∂v+v∂u\partial(uv)=u\partial v+v\partial u

则称这个域为 微分域.

若微分域 FF 上的函数 uu 满足以下的任意一条条件,则称该函数 uu 为初等函数:

  1. uu 是 FF 上的代数函数.
  2. uu 是 FF 上的指数性函数,即存在 a∈Fa\in F 使得 ∂u=u∂a\partial u=u\partial a.
  3. uu 是 FF 上的对数性函数,即存在 a∈Fa\in F 使得 ∂u=∂aa\partial u=\frac{\partial a}{a}.

以下是常见的初等函数:

  1. 代数函数:存在有限次多项式 PP 使得 P(f(x))=0P(f(x))=0 的函数 f(x)f(x),如 2x+12x+1,x\sqrt{x},(1+x2)−1(1+x^2)^{-1},∣x∣|x|.

  2. 指数函数

  3. 对数函数

  4. 三角函数

  5. 反三角函数

  6. 双曲函数

  7. 反双曲函数

  8. 以上函数的复合,如:

    etan⁡x1+x2sin⁡(1+ln⁡2x)\frac{\mathrm{e}^{\tan x}}{1+x^2}\sin\left(\sqrt{1+\ln^2 x}\right) −iln⁡(x+i1−x2)-\mathrm{i} \ln\left(x+\mathrm{i}\sqrt{1-x^2}\right)

以下是常见的非初等函数:

  1. 误差函数:

    erf⁡(x):=2π∫0xexp⁡(−t2)dt\operatorname{erf}(x):=\frac{2}{\sqrt{\pi}}\int_{0}^{x}\exp\left(-t^2\right)\mathrm{d}t

多项式求逆

给定多项式 f(x)f\left(x\right),求 f−1(x)f^{-1}\left(x\right).

解法

倍增法

首先,易知

[x0]f−1(x)=([x0]f(x))−1\left[x^{0}\right]f^{-1}\left(x\right)=\left(\left[x^{0}\right]f\left(x\right)\right)^{-1}

假设现在已经求出了 f(x)f\left(x\right) 在模 x⌈n2⌉x^{\left\lceil\frac{n}{2}\right\rceil} 意义下的逆元 f0−1(x)f^{-1}_{0}\left(x\right). 有:

f(x)f0−1(x)≡1(modx⌈n2⌉)f(x)f−1(x)≡1(modx⌈n2⌉)f−1(x)−f0−1(x)≡0(modx⌈n2⌉)\begin{aligned} f\left(x\right)f^{-1}_{0}\left(x\right)&\equiv 1 &\pmod{x^{\left\lceil\frac{n}{2}\right\rceil}}\\ f\left(x\right)f^{-1}\left(x\right)&\equiv 1 &\pmod{x^{\left\lceil\frac{n}{2}\right\rceil}}\\ f^{-1}\left(x\right)-f^{-1}_{0}\left(x\right)&\equiv 0 &\pmod{x^{\left\lceil\frac{n}{2}\right\rceil}} \end{aligned}

两边平方可得:

f−2(x)−2f−1(x)f0−1(x)+f0−2(x)≡0(modxn)f^{-2}\left(x\right)-2f^{-1}\left(x\right)f^{-1}_{0}\left(x\right)+f^{-2}_{0}\left(x\right)\equiv 0 \pmod{x^{n}}

两边同乘 f(x)f\left(x\right) 并移项可得:

f−1(x)≡f0−1(x)(2−f(x)f0−1(x))(modxn)f^{-1}\left(x\right)\equiv f^{-1}_{0}\left(x\right)\left(2-f\left(x\right)f^{-1}_{0}\left(x\right)\right) \pmod{x^{n}}

递归计算即可.

时间复杂度

T(n)=T(n2)+O(nlog⁡n)=O(nlog⁡n)T\left(n\right)=T\left(\frac{n}{2}\right)+O\left(n\log{n}\right)=O\left(n\log{n}\right)

Newton's Method

参见 Newton's Method.

Graeffe 法

欲求 f−1(x) mod x2nf^{-1}(x)\bmod x^{2n} 考虑

f−1(x) mod x2n=f(−x)(f(x)f(−x))−1 mod x2n=f(−x)g−1(x2) mod x2n\begin{aligned} f^{-1}(x)\bmod x^{2n}&= f(-x)(f(x)f(-x))^{-1}\bmod x^{2n}\\ &=f(-x)g^{-1}(x^2)\bmod x^{2n} \end{aligned}

只需求出 g−1(x) mod xng^{-1}(x)\bmod x^n 即可还原出 g−1(x2) mod x2ng^{-1}(x^2)\bmod x^{2n} 因为 f(x)f(−x)f(x)f(-x) 是偶函数,时间复杂度同上.

代码

多项式求逆
constexpr int MAXN = 262144;
constexpr int mod = 998244353;

using i64 = long long;
using poly_t = int[MAXN];
using poly = int *const;

void polyinv(const poly &h, const int n, poly &f) {
  /* f = 1 / h = f_0 (2 - f_0 h) */
  static poly_t inv_t;
  std::fill(f, f + n + n, 0);
  f[0] = fpow(h[0], mod - 2);
  for (int t = 2; t <= n; t <<= 1) {
    const int t2 = t << 1;
    std::copy(h, h + t, inv_t);
    std::fill(inv_t + t, inv_t + t2, 0);

    DFT(f, t2);
    DFT(inv_t, t2);
    for (int i = 0; i != t2; ++i)
      f[i] = (i64)f[i] * sub(2, (i64)f[i] * inv_t[i] % mod) % mod;
    IDFT(f, t2);

    std::fill(f + t, f + t2, 0);
  }
}

例题

  1. 有标号简单无向连通图计数:「POJ 1737」Connected Graph

多项式开方

给定多项式 g(x)g\left(x\right),求 f(x)f\left(x\right),满足:

f2(x)≡g(x)(modxn)f^{2}\left(x\right)\equiv g\left(x\right) \pmod{x^{n}}

解法

倍增法

首先讨论 [x0]g(x)\left[x^0\right]g(x) 不为 00 的情况.

易知:

[x0]f(x)=[x0]g(x)\left[x^0\right]f(x) = \sqrt{\left[x^0\right]g(x)}

若 [x0]g(x)\left[x^0\right]g(x) 没有平方根,则多项式 g(x)g(x) 没有平方根.

[x0]g(x)\left[x^0\right]g(x) 可能有多个平方根,选取不同的根会求出不同的 f(x)f(x).

假设现在已经求出了 g(x)g\left(x\right) 在模 x⌈n2⌉x^{\left\lceil\frac{n}{2}\right\rceil} 意义下的平方根 f0(x)f_{0}\left(x\right),则有:

f02(x)≡g(x)(modx⌈n2⌉)f02(x)−g(x)≡0(modx⌈n2⌉)(f02(x)−g(x))2≡0(modxn)(f02(x)+g(x))2≡4f02(x)g(x)(modxn)(f02(x)+g(x)2f0(x))2≡g(x)(modxn)f02(x)+g(x)2f0(x)≡f(x)(modxn)2−1f0(x)+2−1f0−1(x)g(x)≡f(x)(modxn)\begin{aligned} f_{0}^{2}\left(x\right)&\equiv g\left(x\right) &\pmod{x^{\left\lceil\frac{n}{2}\right\rceil}}\\ f_{0}^{2}\left(x\right)-g\left(x\right)&\equiv 0 &\pmod{x^{\left\lceil\frac{n}{2}\right\rceil}}\\ \left(f_{0}^{2}\left(x\right)-g\left(x\right)\right)^{2}&\equiv 0 &\pmod{x^{n}}\\ \left(f_{0}^{2}\left(x\right)+g\left(x\right)\right)^{2}&\equiv 4f_{0}^{2}\left(x\right)g\left(x\right) &\pmod{x^{n}}\\ \left(\frac{f_{0}^{2}\left(x\right)+g\left(x\right)}{2f_{0}\left(x\right)}\right)^{2}&\equiv g\left(x\right) &\pmod{x^{n}}\\ \frac{f_{0}^{2}\left(x\right)+g\left(x\right)}{2f_{0}\left(x\right)}&\equiv f\left(x\right) &\pmod{x^{n}}\\ 2^{-1}f_{0}\left(x\right)+2^{-1}f_{0}^{-1}\left(x\right)g\left(x\right)&\equiv f\left(x\right) &\pmod{x^{n}} \end{aligned}

倍增计算即可.

时间复杂度

T(n)=T(n2)+O(nlog⁡n)=O(nlog⁡n)T\left(n\right)=T\left(\frac{n}{2}\right)+O\left(n\log{n}\right)=O\left(n\log{n}\right)

还有一种常数较小的写法就是在倍增维护 f(x)f\left(x\right) 的时候同时维护 f−1(x)f^{-1}\left(x\right) 而不是每次都求逆.

当 [x0]g(x)≠1\left[x^{0}\right]g\left(x\right)\neq 1 时,可能需要使用二次剩余来计算 [x0]f(x)\left[x^{0}\right]f\left(x\right).

上述方法需要知道 f0(x)f_{0}(x) 的逆,所以常数项不能为 00.

若 [x0]g(x)=0\left[x^0\right]g(x) = 0,则将 g(x)g(x) 分解成 xkh(x)x^{k}h(x),其中 [x0]h(x)≠0\left[x^0\right]h(x) \not = 0.

  • 若 kk 是奇数,则 g(x)g(x) 没有平方根.

  • 若 kk 是偶数,则求出 h(x)h(x) 的平方根 h(x)\sqrt{h(x)},然后得到 f(x)≡xk/2h(x)(modxn)f(x) \equiv x^{k/2} \sqrt{h(x)} \pmod{x^{n}}.

洛谷模板题 P5205【模板】多项式开根 参考代码
#include <algorithm>
#include <iostream>
using namespace std;

constexpr int MAXN = 1 << 20, mod = 998244353;

int a[MAXN], b[MAXN], g[MAXN], gg[MAXN];

int qpow(int x, int y) {  // 快速幂
  int ans = 1;

  while (y) {
    if (y & 1) {
      ans = (long long)1 * ans * x % mod;
    }
    x = (long long)1 * x * x % mod;
    y >>= 1;
  }
  return ans;
}

int inv2 = qpow(2, mod - 2);  // 逆元

void change(int *f, int len) {
  for (int i = 1, j = len / 2; i < len - 1; i++) {
    if (i < j) {
      swap(f[i], f[j]);
    }

    int k = len / 2;
    while (j >= k) {
      j -= k;
      k /= 2;
    }
    if (j < k) {
      j += k;
    }
  }
}

void NTT(int *f, int len, int type) {  // NTT
  change(f, len);

  for (int q = 2; q <= len; q <<= 1) {
    int nxt = qpow(3, (mod - 1) / q);
    for (int i = 0; i < len; i += q) {
      int w = 1;

      for (int k = i; k < i + (q >> 1); k++) {
        int x = f[k];
        int y = (long long)1 * w * f[k + (q / 2)] % mod;

        f[k] = (x + y) % mod;
        f[k + (q / 2)] = (x - y + mod) % mod;
        w = (long long)1 * w * nxt % mod;
      }
    }
  }

  if (type == -1) {
    reverse(f + 1, f + len);
    int iv = qpow(len, mod - 2);

    for (int i = 0; i < len; i++) {
      f[i] = (long long)1 * f[i] * iv % mod;
    }
  }
}

void inv(int deg, int *f, int *h) {  // 求逆元
  if (deg == 1) {
    h[0] = qpow(f[0], mod - 2);
    return;
  }

  inv((deg + 1) >> 1, f, h);

  int len = 1;
  while (len < deg * 2) {  // 倍增
    len *= 2;
  }

  copy(f, f + deg, gg);
  fill(gg + deg, gg + len, 0);

  NTT(gg, len, 1);
  NTT(h, len, 1);
  for (int i = 0; i < len; i++) {
    h[i] = (long long)1 * (2 - (long long)1 * gg[i] * h[i] % mod + mod) % mod *
           h[i] % mod;
  }

  NTT(h, len, -1);
  fill(h + deg, h + len, 0);
}

int n, t[MAXN];

// deg:次数
// f:被开根数组
// h:答案数组
void sqrt(int deg, int *f, int *h) {
  if (deg == 1) {
    h[0] = 1;
    return;
  }

  sqrt((deg + 1) >> 1, f, h);

  int len = 1;
  while (len < deg * 2) {  // 倍增
    len *= 2;
  }
  fill(g, g + len, 0);
  inv(deg, h, g);
  copy(f, f + deg, t);
  fill(t + deg, t + len, 0);
  NTT(t, len, 1);
  NTT(g, len, 1);
  NTT(h, len, 1);

  for (int i = 0; i < len; i++) {
    h[i] = (long long)1 * inv2 *
           ((long long)1 * h[i] % mod + (long long)1 * g[i] * t[i] % mod) % mod;
  }
  NTT(h, len, -1);
  fill(h + deg, h + len, 0);
}

int main() {
  cin >> n;

  for (int i = 0; i < n; i++) {
    cin >> a[i];
  }
  sqrt(n, a, b);

  for (int i = 0; i < n; i++) {
    cout << b[i] << ' ';
  }

  return 0;
}

Newton's Method

参见 Newton's Method.

例题

  1. 「Codeforces Round #250」E. The Child and Binary Tree

多项式除法 & 取模

给定多项式 f(x),g(x)f\left(x\right),g\left(x\right),求 g(x)g\left(x\right) 除 f(x)f\left(x\right) 的商 Q(x)Q\left(x\right) 和余数 R(x)R\left(x\right).

解法

发现若能消除 R(x)R\left(x\right) 的影响则可直接 多项式求逆 解决.

考虑构造变换

fR(x)=xdeg⁡ff(1x)f^{R}\left(x\right)=x^{\operatorname{deg}{f}}f\left(\frac{1}{x}\right)

观察可知其实质为反转 f(x)f\left(x\right) 的系数.

设 n=deg⁡f,m=deg⁡gn=\operatorname{deg}{f},m=\operatorname{deg}{g}.

将 f(x)=Q(x)g(x)+R(x)f\left(x\right)=Q\left(x\right)g\left(x\right)+R\left(x\right) 中的 xx 替换成 1x\frac{1}{x} 并将其两边都乘上 xnx^{n},得到:

xnf(1x)=xn−mQ(1x)xmg(1x)+xn−m+1xm−1R(1x)fR(x)=QR(x)gR(x)+xn−m+1RR(x)\begin{aligned} x^{n}f\left(\frac{1}{x}\right)&=x^{n-m}Q\left(\frac{1}{x}\right)x^{m}g\left(\frac{1}{x}\right)+x^{n-m+1}x^{m-1}R\left(\frac{1}{x}\right)\\ f^{R}\left(x\right)&=Q^{R}\left(x\right)g^{R}\left(x\right)+x^{n-m+1}R^{R}\left(x\right) \end{aligned}

注意到上式中 RR(x)R^{R}\left(x\right) 的系数为 xn−m+1x^{n-m+1},则将其放到模 xn−m+1x^{n-m+1} 意义下即可消除 RR(x)R^{R}\left(x\right) 带来的影响.

又因 QR(x)Q^{R}\left(x\right) 的次数为 (n−m)<(n−m+1)\left(n-m\right)<\left(n-m+1\right),故 QR(x)Q^{R}\left(x\right) 不会受到影响.

则:

fR(x)≡QR(x)gR(x)(modxn−m+1)f^{R}\left(x\right)\equiv Q^{R}\left(x\right)g^{R}\left(x\right)\pmod{x^{n-m+1}}

使用多项式求逆即可求出 Q(x)Q\left(x\right),将其反代即可得到 R(x)R\left(x\right).

时间复杂度 O(nlog⁡n)O\left(n\log{n}\right).

多项式对数函数 & 指数函数

给定多项式 f(x)f(x),求模 xnx^{n} 意义下的 ln⁡f(x)\ln{f(x)} 与 exp⁡f(x)\exp{f(x)}.

解法

普通方法

多项式对数函数

首先,对于多项式 f(x)f(x),若 ln⁡f(x)\ln{f(x)} 存在,则由其 定义,其必须满足:

[x0]f(x)=1[x^{0}]f(x)=1

对 ln⁡f(x)\ln{f(x)} 求导再积分,可得:

dln⁡f(x)dx≡f′(x)f(x)(modxn)ln⁡f(x)≡∫dln⁡f(x)≡∫f′(x)f(x)dx(modxn)\begin{aligned} \frac{\mathrm{d} \ln{f(x)}}{\mathrm{d} x} & \equiv \frac{f'(x)}{f(x)} & \pmod{x^{n}} \\ \ln{f(x)} & \equiv \int \mathrm{d} \ln{f(x)} \equiv \int\frac{f'(x)}{f(x)} \mathrm{d} x & \pmod{x^{n}} \end{aligned}

多项式的求导,积分时间复杂度为 O(n)O(n),求逆时间复杂度为 O(nlog⁡n)O(n\log{n}),故多项式求 ln⁡\ln 时间复杂度 O(nlog⁡n)O(n\log{n}).

多项式指数函数

首先,对于多项式 f(x)f(x),若 exp⁡f(x)\exp{f(x)} 存在,则其必须满足:

[x0]f(x)=0[x^{0}]f(x)=0

否则 exp⁡f(x)\exp{f(x)} 的常数项不收敛.

对 exp⁡f(x)\exp{f(x)} 求导,可得:

dexp⁡f(x)dx≡exp⁡f(x)f′(x)(modxn)\frac{\mathrm{d} \exp{f(x)}}{\mathrm{d} x} \equiv \exp{f(x)}f'(x)\pmod{x^{n}}

比较两边系数可得:

[xn−1]dexp⁡f(x)dx=∑i=0n−1([xi]exp⁡f(x))([xn−i−1]f′(x))[x^{n-1}]\frac{\mathrm{d} \exp{f(x)}}{\mathrm{d} x} = \sum_{i = 0}^{n - 1} \left([x^{i}]\exp{f(x)}\right) \left([x^{n-i-1}]f'(x)\right) n[xn]exp⁡f(x)=∑i=0n−1([xi]exp⁡f(x))((n−i)[xn−i]f(x))n[x^{n}]\exp{f(x)} = \sum_{i = 0}^{n - 1} \left([x^{i}]\exp{f(x)}\right) \left((n - i)[x^{n - i}]f(x)\right)

使用分治 FFT 即可解决.

时间复杂度 O(nlog⁡2n)O(n\log^{2}{n}).

Newton's Method

使用 Newton's Method 即可在 O(nlog⁡n)O(n\log{n}) 的时间复杂度内解决多项式 exp⁡\exp.

代码

多项式 ln/exp
constexpr int MAXN = 262144;
constexpr int mod = 998244353;

using i64 = long long;
using poly_t = int[MAXN];
using poly = int *const;

void derivative(const poly &h, const int n, poly &f) {
  for (int i = 1; i != n; ++i) f[i - 1] = (i64)h[i] * i % mod;
  f[n - 1] = 0;
}

void integrate(const poly &h, const int n, poly &f) {
  for (int i = n - 1; i; --i) f[i] = (i64)h[i - 1] * inv[i] % mod;
  f[0] = 0; /* C */
}

void polyln(const poly &h, const int n, poly &f) {
  /* f = ln h = ∫ h' / h dx */
  assert(h[0] == 1);
  static poly_t ln_t;
  const int t = n << 1;

  derivative(h, n, ln_t);
  std::fill(ln_t + n, ln_t + t, 0);
  polyinv(h, n, f);

  DFT(ln_t, t);
  DFT(f, t);
  for (int i = 0; i != t; ++i) ln_t[i] = (i64)ln_t[i] * f[i] % mod;
  IDFT(ln_t, t);

  integrate(ln_t, n, f);
}

void polyexp(const poly &h, const int n, poly &f) {
  /* f = exp(h) = f_0 (1 - ln f_0 + h) */
  assert(h[0] == 0);
  static poly_t exp_t;
  std::fill(f, f + n + n, 0);
  f[0] = 1;
  for (int t = 2; t <= n; t <<= 1) {
    const int t2 = t << 1;

    polyln(f, t, exp_t);
    exp_t[0] = sub(pls(h[0], 1), exp_t[0]);
    for (int i = 1; i != t; ++i) exp_t[i] = sub(h[i], exp_t[i]);
    std::fill(exp_t + t, exp_t + t2, 0);

    DFT(f, t2);
    DFT(exp_t, t2);
    for (int i = 0; i != t2; ++i) f[i] = (i64)f[i] * exp_t[i] % mod;
    IDFT(f, t2);

    std::fill(f + t, f + t2, 0);
  }
}

例题

  1. 计算 fk(x)f^{k}(x)

    普通做法为多项式快速幂,时间复杂度 O(nlog⁡nlog⁡k)O(n\log{n}\log{k}).

    当 [x0]f(x)=1[x^{0}]f(x)=1 时,有:

    fk(x)=exp⁡(kln⁡f(x))f^{k}(x)=\exp{\left(k\ln{f(x)}\right)}

    当 [x0]f(x)≠1[x^{0}]f(x)\neq 1 时,设 f(x)f(x) 的最低次项为 fixif_{i}x^{i},则:

    fk(x)=fikxikexp⁡(kln⁡f(x)fixi)f^{k}(x)=f_{i}^{k}x^{ik}\exp{\left(k\ln{\frac{f(x)}{f_{i}x^{i}}}\right)}

    时间复杂度 O(nlog⁡n)O(n\log{n}).

多项式三角函数

给定多项式 f(x)f\left(x\right),求模 xnx^{n} 意义下的 sin⁡f(x),cos⁡f(x)\sin{f\left(x\right)}, \cos{f\left(x\right)} 与 tan⁡f(x)\tan{f\left(x\right)}.

解法

首先由 Euler's formula (eix=cos⁡x+isin⁡x)\left(\mathrm{e}^{\mathrm{i}x} = \cos{x} + \mathrm{i}\sin{x}\right) 可以得到 三角函数的另一个表达式:

sin⁡x=eix−e−ix2icos⁡x=eix+e−ix2\begin{aligned} \sin{x} &= \frac{\mathrm{e}^{\mathrm{i}x} - \mathrm{e}^{-\mathrm{i}x}}{2\mathrm{i}} \\ \cos{x} &= \frac{\mathrm{e}^{\mathrm{i}x} + \mathrm{e}^{-\mathrm{i}x}}{2} \end{aligned}

那么代入 f(x)f\left(x\right) 就有:

sin⁡f(x)=exp⁡(if(x))−exp⁡(−if(x))2icos⁡f(x)=exp⁡(if(x))+exp⁡(−if(x))2\begin{aligned} \sin{f\left(x\right)} &= \frac{\exp{\left(\mathrm{i}f\left(x\right)\right)} - \exp{\left(-\mathrm{i}f\left(x\right)\right)}}{2\mathrm{i}} \\ \cos{f\left(x\right)} &= \frac{\exp{\left(\mathrm{i}f\left(x\right)\right)} + \exp{\left(-\mathrm{i}f\left(x\right)\right)}}{2} \end{aligned}

直接按上述表达式编写程序即可得到模 xnx^{n} 意义下的 sin⁡f(x)\sin{f\left(x\right)} 与 cos⁡f(x)\cos{f\left(x\right)}.再由 tan⁡f(x)=sin⁡f(x)cos⁡f(x)\tan{f\left(x\right)} = \frac{\sin{f\left(x\right)}}{\cos{f\left(x\right)}} 可求得 tan⁡f(x)\tan{f\left(x\right)}.

代码

多项式三角函数

注意到我们是在 Z998244353\mathbb{Z}_{998244353} 上做 NTT,那么相应地,虚数单位 i\mathrm{i} 应该被换成 8658371886583718 或 911660635911660635:

i=−1≡998244352(mod998244353)  ⟹  ori≡86583718(mod998244353)ori≡911660635(mod998244353)\begin{aligned} & \mathrm{i} = \sqrt{-1} \equiv \sqrt{998244352} \pmod{998244353} \\ \implies & \phantom{\text{or}} \quad \mathrm{i} \equiv 86583718 \pmod{998244353} \\ & \text{or} \quad \mathrm{i} \equiv 911660635 \pmod{998244353} \end{aligned}
constexpr int MAXN = 262144;
constexpr int mod = 998244353;
constexpr int imgunit = 86583718; /* sqrt(-1) = sqrt(998233452) */

using i64 = long long;
using poly_t = int[MAXN];
using poly = int *const;

void polytri(const poly &h, const int n, poly &sin_t, poly &cos_t) {
  /* sin(f) = (exp(i * f) - exp(- i * f)) / 2i */
  /* cos(f) = (exp(i * f) + exp(- i * f)) / 2 */
  /* tan(f) = sin(f) / cos(f) */
  assert(h[0] == 0);
  static poly_t tri1_t, tri2_t;

  for (int i = 0; i != n; ++i) tri2_t[i] = (i64)h[i] * imgunit % mod;
  polyexp(tri2_t, n, tri1_t);
  polyinv(tri1_t, n, tri2_t);

  if (sin_t != nullptr) {
    const int invi = fpow(pls(imgunit, imgunit), mod - 2);
    for (int i = 0; i != n; ++i)
      sin_t[i] = (i64)(tri1_t[i] - tri2_t[i] + mod) * invi % mod;
  }
  if (cos_t != nullptr) {
    for (int i = 0; i != n; ++i) cos_t[i] = div2(pls(tri1_t[i], tri2_t[i]));
  }
}

多项式反三角函数

给定多项式 f(x)f\left(x\right),求模 xnx^{n} 意义下的 arcsin⁡f(x),arccos⁡f(x)\arcsin{f\left(x\right)}, \arccos{f\left(x\right)} 与 arctan⁡f(x)\arctan{f\left(x\right)}.

解法

仿照求多项式 ln⁡\ln 的方法,对反三角函数求导再积分可得:

ddxarcsin⁡x=11−x2arcsin⁡x=∫11−x2dxddxarccos⁡x=−11−x2arccos⁡x=−∫11−x2dxddxarctan⁡x=11+x2arctan⁡x=∫11+x2dx\begin{aligned} \frac{\mathrm{d}}{\mathrm{d} x} \arcsin{x} &= \frac{1}{\sqrt{1 - x^{2}}} \\ \arcsin{x} &= \int \frac{1}{\sqrt{1 - x^{2}}} \mathrm{d} x \\ \frac{\mathrm{d}}{\mathrm{d} x} \arccos{x} &= - \frac{1}{\sqrt{1 - x^{2}}} \\ \arccos{x} &= - \int \frac{1}{\sqrt{1 - x^{2}}} \mathrm{d} x \\ \frac{\mathrm{d}}{\mathrm{d} x} \arctan{x} &= \frac{1}{1 + x^{2}} \\ \arctan{x} &= \int \frac{1}{1 + x^{2}} \mathrm{d} x \end{aligned}

那么代入 f(x)f\left(x\right) 就有:

ddxarcsin⁡f(x)=f′(x)1−f2(x)arcsin⁡f(x)=∫f′(x)1−f2(x)dxddxarccos⁡f(x)=−f′(x)1−f2(x)arccos⁡f(x)=−∫f′(x)1−f2(x)dxddxarctan⁡f(x)=f′(x)1+f2(x)arctan⁡f(x)=∫f′(x)1+f2(x)dx\begin{aligned} \frac{\mathrm{d}}{\mathrm{d} x} \arcsin{f\left(x\right)} &= \frac{f'\left(x\right)}{\sqrt{1 - f^{2}\left(x\right)}} \\ \arcsin{f\left(x\right)} &= \int \frac{f'\left(x\right)}{\sqrt{1 - f^{2}\left(x\right)}} \mathrm{d} x \\ \frac{\mathrm{d}}{\mathrm{d} x} \arccos{f\left(x\right)} &= - \frac{f'\left(x\right)}{\sqrt{1 - f^{2}\left(x\right)}} \\ \arccos{f\left(x\right)} &= - \int \frac{f'\left(x\right)}{\sqrt{1 - f^{2}\left(x\right)}} \mathrm{d} x \\ \frac{\mathrm{d}}{\mathrm{d} x} \arctan{f\left(x\right)} &= \frac{f'\left(x\right)}{1 + f^{2}\left(x\right)} \\ \arctan{f\left(x\right)} &= \int \frac{f'\left(x\right)}{1 + f^{2}\left(x\right)} \mathrm{d} x \end{aligned}

直接按式子求就可以了.

代码

多项式反三角函数
constexpr int MAXN = 262144;
constexpr int mod = 998244353;

using i64 = long long;
using poly_t = int[MAXN];
using poly = int *const;

void derivative(const poly &h, const int n, poly &f) {
  for (int i = 1; i != n; ++i) f[i - 1] = (i64)h[i] * i % mod;
  f[n - 1] = 0;
}

void integrate(const poly &h, const int n, poly &f) {
  for (int i = n - 1; i; --i) f[i] = (i64)h[i - 1] * inv[i] % mod;
  f[0] = 0; /* C */
}

void polyarcsin(const poly &h, const int n, poly &f) {
  /* arcsin(f) = ∫ f' / sqrt(1 - f^2) dx  */
  static poly_t arcsin_t;
  const int t = n << 1;
  std::copy(h, h + n, arcsin_t);
  std::fill(arcsin_t + n, arcsin_t + t, 0);

  DFT(arcsin_t, t);
  for (int i = 0; i != t; ++i) arcsin_t[i] = sqr(arcsin_t[i]);
  IDFT(arcsin_t, t);

  arcsin_t[0] = sub(1, arcsin_t[0]);
  for (int i = 1; i != n; ++i)
    arcsin_t[i] = arcsin_t[i] ? mod - arcsin_t[i] : 0;

  polysqrt(arcsin_t, n, f);
  polyinv(f, n, arcsin_t);
  derivative(h, n, f);

  DFT(f, t);
  DFT(arcsin_t, t);
  for (int i = 0; i != t; ++i) arcsin_t[i] = (i64)f[i] * arcsin_t[i] % mod;
  IDFT(arcsin_t, t);

  integrate(arcsin_t, n, f);
}

void polyarccos(const poly &h, const int n, poly &f) {
  /* arccos(f) = - ∫ f' / sqrt(1 - f^2) dx  */
  polyarcsin(h, n, f);
  for (int i = 0; i != n; ++i) f[i] = f[i] ? mod - f[i] : 0;
}

void polyarctan(const poly &h, const int n, poly &f) {
  /* arctan(f) = ∫ f' / (1 + f^2) dx  */
  static poly_t arctan_t;
  const int t = n << 1;
  std::copy(h, h + n, arctan_t);
  std::fill(arctan_t + n, arctan_t + t, 0);

  DFT(arctan_t, t);
  for (int i = 0; i != t; ++i) arctan_t[i] = sqr(arctan_t[i]);
  IDFT(arctan_t, t);

  inc(arctan_t[0], 1);
  std::fill(arctan_t + n, arctan_t + t, 0);

  polyinv(arctan_t, n, f);
  derivative(h, n, arctan_t);

  DFT(f, t);
  DFT(arctan_t, t);
  for (int i = 0; i != t; ++i) arctan_t[i] = (i64)f[i] * arctan_t[i] % mod;
  IDFT(arctan_t, t);

  integrate(arctan_t, n, f);
}

参考资料与链接

Footnotes

  1. Elementary function——Wikipedia ↩

本文内容来自 OI Wiki《多项式初等函数》;作者与修订历史见原文。AIOJ 对排版、站内链接和资源路径做了格式适配。

内容依 CC BY-SA 4.0 与 SATA 许可共享;附加条款以原文声明为准。感谢 OI Wiki 作者与贡献者;请按 SATA 约定前往 OI Wiki 仓库 Star 支持。