多项式学习笔记

符号与约定

若无特殊说明,多项式的大小均默认为 \(n=2^k\)

我们定义 \([x^i]F(x)\) 表示 \(F(x)\) 的第 \(i\) 项系数,其中 \(F(x)\) 为多项式。那么有 \(F(x)=\displaystyle\sum_{i=0}^{n-1}[x^i]F(x)\cdot x^i\)

我们使用 \(F(x)\bmod x^n\) 来限定 \(F(x)\) 的次数上界,例如 \((x^2+2x+1)\bmod x^2=2x+1\)。显然模 \(x^n\) 的多项式域对加减乘封闭(因为只会从低位转移到高位)

基础模板

const int kM = 998244353;

struct Mi {
  int v;

  Mi(int v = 0) : v(v) {}
  Mi operator=(int v) { return *this = Mi(v); }
  Mi operator+(const Mi &o) const { return Mi(v + o.v - (v + o.v >= kM ? kM : 0)); }
  Mi operator-(const Mi &o) const { return Mi(v - o.v + (v - o.v < 0 ? kM : 0)); }
  Mi operator*(const Mi &o) const { return Mi(1LL * v * o.v % kM); }
  Mi operator~() const;
  Mi operator/(const Mi &o) const { return *this * ~o; }
  Mi operator+=(const Mi &o) { return *this = *this + o; }
  Mi operator-=(const Mi &o) { return *this = *this - o; }
  Mi operator*=(const Mi &o) { return *this = *this * o; }
  Mi operator/=(const Mi &o) { return *this = *this / o; }
};
istream &operator>>(istream &in, Mi &o) {
  LL v;
  in >> v;
  o = v;
  return in;
}
ostream &operator<<(ostream &out, const Mi &o) { return out << o.v; }
Mi P(Mi b, int e) {
  Mi s = 1;
  for (; e; e >>= 1, b *= b) {
    if (e & 1) {
      s *= b;
    }
  }
  return s;
}
Mi Mi::operator~() const { return P(this->v, kM - 2); }
const Mi kG = 3, kiG = ~kG;
using Poly = vector<Mi>;

int Scale(int n) {
  int l = 1;
  for (; l < n; l <<= 1) {
  }
  return l;
}
const int kL = 22;
Mi kO[kL], kiO[kL];
int _Init_omega = []() {
  for (int i = 0; i < kL; ++i) {
    kO[i] = P(kG, (kM - 1) / (1 << i));
    kiO[i] = P(kiG, (kM - 1) / (1 << i));
  }
  return 0;
}();
int kBr[1 << kL];
void InitBr(int n) {
  for (int i = 0; i < n; ++i) {
    kBr[i] = kBr[i >> 1] >> 1 | (i & 1) * (n / 2);
  }
}
int mxn;
void Init(int mx) {
  mxn = Scale(mx);
  InitBr(mxn);
}

加、减、多项式乘单项式

\[[x^i](F\pm G)(x)=[x^i]F(x)\pm [x^i]G(x)\\ [x^i](vF)(x)=v[x^i]F(x) \]

Poly operator+(const Poly &a, const Poly &b) {
  Poly s(max(a.size(), b.size()));
  for (int i = 0; i < s.size(); ++i) {
    s[i] = (i < a.size() ? a[i] : 0) + (i < b.size() ? b[i] : 0);
  }
  return s;
}
Poly operator-(const Poly &a, const Poly &b) {
  Poly s(max(a.size(), b.size()));
  for (int i = 0; i < s.size(); ++i) {
    s[i] = (i < a.size() ? a[i] : 0) - (i < b.size() ? b[i] : 0);
  }
  return s;
}
Poly operator*(const Poly &a, const Mi &b) {
  Poly s = a;
  for (Mi &v : s) {
    v *= b;
  }
  return s;
}

多项式乘法/卷积

详见我的 深入理解 FFT

知周所众,浮点运算 \(\approx\) 精度误差+大常熟,所以拉了/ruo。

考虑在模意义下能够替代单位根的东西。

可以发现,原根的性质和单位根很像,设 \(g\) 是某个模数的原根,那么可以发现 \(g^k\) 的阶数即为 \(\dfrac{p-1}{\gcd(k,p-1)}\),即 \(g^k=\omega_{(p-1)/\gcd(k,p-1)}\),代入 \(k=(p-1)/n\) 可得 \(g^{(p-1)/n}=\omega_{(p-1)/\gcd((p-1)/n,p-1)}=\omega_{(p-1)/((p-1)/n)}=\omega_n\)。由于我们默认 \(n=2^k\),我们只需要让模数 \(p\) 形如 \(2^k\cdot r+1\),其中 \(k\) 是个不小的数即可。经典模数 \(998244353=2^{23}\times 7\times 17+1\) 就符合此条件。

void NTT(const Poly &a, Poly &b, bool iv) {
  int n = a.size();
  b.resize(n);
  for (int i = 0; i < n; ++i) {
    b[i] = a[kBr[i]];
  }
  for (int l = 1, c = 0; l < n; l <<= 1, ++c) {
    Mi bo = (iv ? kiO : kO)[c + 1], so = 1;
    for (int i = 0; i < l; ++i) {
      for (int j = 0; j < n; j += l << 1) {
        Mi v = so * b[j + l + i];
        b[j + l + i] = b[j + i] - v;
        b[j + i] += v;
      }
      so *= bo;
    }
  }
  if (iv) {
    Mi _v = ~Mi(n);
    for (Mi &v : b) {
      v *= _v;
    }
  }
}
void iNTT(Poly &a, bool iv) {
  Poly b;
  NTT(a, b, iv);
  a = b;
}
Poly operator|(const Poly &x, const Poly &y) {
  Poly s = x;
  for (int i = 0; i < x.size(); ++i) {
    s[i] *= y[i];
  }
  return s;
}
Poly operator*(const Poly &x, const Poly &y) {
  Poly _x, _y;
  NTT(x, _x, 0), NTT(y, _y, 0);
  NTT(_x | _y, _x, 1);
  return _x;
}

多项式求逆

我们有这样一个问题:给定 \(F(x)\),求一个多项式 \(G(x)\) 使得 \(F(x)\times G(x)\equiv 1\pmod {x^n}\)

考虑对界进行倍增,边界情况为 \(F(x)\times G(x)\equiv 1\pmod{x}\),此时对 \(F(x)\) 的常数项取逆即可。

假设我们现在已经知道了 \(G'(x)\equiv F(x)^{-1}\pmod{x^{n/2}}\),设 \(G(x)\equiv F(x)^{-1}\pmod{x^n}\),显然有 \(G(x)\equiv G'(x)\pmod{x^{n/2}}\),于是有:

\[G(x)\equiv G'(x)\pmod{x^{n/2}}\\ G(x)-G'(x)\equiv 0\pmod{x^{n/2}}\\ G(x)^2+G'(x)^2-2G(x)G'(x)\equiv 0\pmod{x^n}\\ G(x)+F(x)G'(x)^2-2G'(x)\equiv 0\pmod{x^n}\\ G(x)\equiv 2G'(x)-F(x)G'(x)^2\pmod{x^n}\\ G(x)\equiv G'(x)(2-F(x)G'(x))\pmod{x^n}\\ \]

Poly operator~(const Poly &a) {
  int n = a.size();
  Poly s(1, ~a[0]);
  for (int l = 2; l <= n; l <<= 1) {
    InitBr(l << 1);
    Poly b = a;
    b.resize(l), b.resize(l << 1);
    s.resize(l << 1);
    iNTT(s, 0), iNTT(b, 0);
    for (int i = 0; i < (l << 1); ++i) {
      s[i] = s[i] * (Mi(2) - b[i] * s[i]);
    }
    iNTT(s, 1);
    s.resize(l);
  }
  return s;
}

求导/积分

\[\frac{\mathrm{d}F(x)}{\mathrm{d}x}=\sum_{i=0}^{n-2}(i+1)[x^{i+1}]F(x)\cdot x^i\\ \int F(x)\mathrm{d}x=C+\sum_{i=1}^{n-1}\frac{[x^{i-1}]F(x)}{i}\cdot x^i \]

Poly Deriv(const Poly &a) {
  int n = a.size();
  Poly s(n);
  for (int i = 0; i < n - 1; ++i) {
    s[i] = a[i + 1] * (i + 1);
  }
  return s;
}
Poly Integ(const Poly &a) {
  int n = a.size();
  Poly s(n);
  for (int i = 1; i < n; ++i) {
    s[i] = a[i - 1] / i;
  }
  return s;
}

牛顿迭代

极其重要,用来解决这种问题:

已知 \(G(x)\)\(G(F(x))=0\),求 \(F(x)\pmod{x^n}\)

假设我们已知 \(F_*(x)\) 使得 \(G(F_*(x))\equiv 0\pmod{x^{n/2}}\),则:

\[F(x)\equiv F_*(x)-\frac{G(F_*(x))}{G'(F_*(x))}\pmod{x^n} \]

计算时的一些注意点:

  • \(F(x)\) 的精度是 \(\pmod{x^n}\) 的。
  • \(F_*(x)\) 的精度只需做到 \(\pmod{x^{n/2}}\)
  • \(G(F_*(x))\) 的最低次项至少是 \(x^{n/2}\),所以 \(G'(F_*(x))\) 的精度只需做到 \(\pmod{x^{n/2}}\)
  • \(G'(F(x))\)\(\dfrac{\mathrm{d}G(F(x))}{\mathrm{d}F(x)}\) 的简写而不是 \(\dfrac{\mathrm{d}G(F(x))}{\mathrm{d}x}\) 的简写。

开根

给定 \(G(F(x))=F(x)^2-A(x)\),求 \(G(F(x))=0\) 的解 \(F(x)\pmod{x^n}\)

发现 \(G'(F(x))=2F(x)\),于是有:

\[F(x)\equiv F_*(x)-\frac{F_*(x)^2-A(x)}{2F_*(x)}\pmod{x^n}\\ F(x)\equiv F_*(x)-(\frac{F_*(x)}{2}-\frac{A(x)}{2F_*(x)})\pmod{x^n}\\ F(x)\equiv \frac{1}{2}(F_*(x)+\frac{A(x)}{F_*(x)})\pmod{x^n}\\ \]

当只有常数项是,开放结果是二次剩余,模板题保证了 \(a_0=1\),故令 \(F_0=1\) 即可。

Poly Sqrt(const Poly &a) {
  int n = a.size();
  Poly s(1, 1);
  Mi iv2 = ~Mi(2);
  for (int l = 2; l <= n; l <<= 1) {
    Poly b = a;
    b.resize(l), b.resize(l << 1);
    Poly _s = ~s;
    _s.resize(l << 1);
    InitBr(l << 1);
    s = (s + b * _s) * iv2;
    s.resize(l);
  }
  return s;
}

带余除法

首先直接求 \(\dfrac{F(x)}{G(x)}\) 的话可以直接对 \(G(x)\) 求逆,然后和 \(F(x)\) 乘起来。

但这样子得到的多项式是 \(n+m\) 次项的,而我们想要的是 \(n-m\) 此项。

考虑原式:

\[F(x)=G(x)\times Q(x)+R(x) \]

我们的 \(\pmod{x^n}\) 可以做到消去 \(x^n\) 及以上的项,但我们此时想要消去的是 \(x^{m}\) 以下的项(即消去 \(R(x)\))。

考虑将多项式反转。我们有 \(F^T(x)=x^nF(x^{-1})\),那么 \(F^T(x)\) 的系数则正好是 \(F(x)\) 的系数反转得到的结果。

那么有:

\[F(x^{-1})=G(x^{-1})\times Q(x^{-1})+R(x^{-1})\\ x^nF(x^{-1})=x^mG(x^{-1})\times x^{n-m}Q(x^{-1})+x^nR(x^{-1})\\ F^T(x)=G^T(x)\times Q^T(x)+x^{n-(m-1)}R^T(x)\\ F^T(x)\equiv G^T(x)\times Q^T(x)\pmod{x^{n-m+1}}\\ Q^T(x)\equiv \frac{F^T(x)}{G^T(x)}\pmod{x^{n-m+1}}\\ \]

此时就能用多项式求逆做了,求出 \(Q^T(x)\) 后再把系数翻回来即可。

我们有:

\[R(x)=F(x)-G(x)\times Q(x) \]

因此 \(R(x)\) 是容易计算的。

posted @ 2023-04-18 18:30  bykem  阅读(101)  评论(0)    收藏  举报