普通多項式筆記
\(\textrm{Newton's Method}\) （牛頓迭代）
應用於解決已知 \(g(x)\) 的情況下，求出 \(g(f(x))\equiv 0\mod x^n\)。
首先透過列出方程顯然，\(f(x) \mod x^n\) 在此時是唯一的。
那麼我們假設已知 \(g(f_0(x))\equiv 0\mod x^{n/2}\)，顯然此時 \(f_0(x)\mod x^{n/2}\) 也是唯一的，且 \(\deg f_0(x)=\frac{n}{2},f(x)\equiv f_0(x) \mod x^{n/2}\)。
這時我們把 \(g(f(x))\) 在 \(f_0(x)\) 處展開，可以得到：
\[g(f(x))= \sum_{k=0}^{+\infty}\frac{g^{(k)}(f_0(x))}{k!}(f(x)-f_0(x))^k
\]
因此有：
\[\sum_{k=0}^{+\infty}\frac{g^{(k)}(f_0(x))}{k!}(f(x)-f_0(x))^k\equiv 0 \mod x^n
\]
因為 \(f(x)-f_0(x)\equiv 0 \mod x^{n/2}\)，所以有： \((f(x)-f_0(x))^2 \equiv 0 \mod x^n\)。
當然 \(n=1\) 時應該要特殊地求出 \(f(x)\)。
所以：
\[\sum_{k=0}^{+\infty}\frac{g^{(k)}(f_0(x))}{k!}(f(x)-f_0(x))^k\equiv g(f_0(x))+g'(f_0(x))(f(x)-f_0(x))\equiv 0 \mod x^n
\\
f(x)\equiv f_0(x)-\frac{g(f_0(x))}{g'(f_0(x))} \mod x^n
\]
最終的結果看起來沒什麼令人驚奇的地方，但是實際上應用卻是非常多的。
多項式 \(\textrm{inverse}\)，逆
已知 \(h(x)\)，求 \(f(x)=h^{-1}(x)\)。
考慮使用 \(\textrm{Newton's Method}\) 解決這個問題。
構造 \(g(f(x))=f^{-1}(x)-h(x)\)，此時我們把 \(h(x)\) 當作一個常數。
因此有：
\[\begin{aligned}
f(x)&\equiv f_0(x)-\frac{f^{-1}_0(x)-h(x)}{-f_0^{-2}(x)} &\mod x^n \\
&\equiv f_0(x)_+f_0^2(x)\left(f_0^{-1}(x)-h(x)\right) &\mod x^n
\end{aligned}
\]
時間複雜度 \(O(n\log_2 n)\)。
多項式 \(\textrm{division}\)，帶餘除法
已知 \(A(x),B(x)\)，求 \(Q(x),R(x)\) 使得 \(A(x)=B(x)\cdot Q(x)+R(x)\)。
記 \(\deg A(x)=n,\deg B(x)=m,n\ge m\)，那麼欽定 \(\deg Q(x)=n-m,\deg R(x)=m-1\)。
記對多項式 \(f(x)\) 的 \(\textrm{Reverse}\) 操作為 \(f^R(x)=x^{\deg f(x)}f\left(\frac1 x\right)\)。
因此有：
\[A(\frac{1}{x})=B(\frac{1}{x})\cdot Q(\frac{1}{x}) +R(\frac{1}{x})
\\
x^nA(\frac{1}{x})=x^{m}B(\frac{1}{x})\cdot x^{n-m}Q(\frac{1}{x}) +x^{n-m+1}\cdot x^{m-1}R(\frac{1}{x})
\\
A^R(x)=B^R(x)\cdot Q^R(x)+x^{n-m+1} R^R(x)
\\
A^R(x)\equiv B^R(x)\cdot Q^R(x) \mod x^{n-m+1}
\]
可以求出 \(Q^R(x)\) 進而求出 \(Q(x),R(x)\)。
多項式 \(\textrm{sqrt}\)，開方
已知 \(h(x)\)，求 \(f(x)=\sqrt {h(x)}\)。
考慮使用 \(\textrm{Newton's Method}\) 解決這個問題。
構造 \(g(f(x))=f^2(x)-h(x)\)，此時我們把 \(h(x)\) 當作一個常數。
因此有：
\[\begin{aligned}
f(x)\equiv f_0(x)-\frac{f_0^2(x)-h(x)}{2f_0(x)} \mod x^n
\end{aligned}
\]
時間複雜度 \(O(n\log_2n)\)。
多項式 \(\ln\)，自然對數
已知 \(f(x)\)，求 \(\ln f(x)\)。
考慮求導：
\[(\ln f(x))'=\frac{f'(x)}{f(x)}
\]
積分還原後：
\[\ln f(x)=\int \frac{f'(x)}{f(x)}
\]
有一個問題是最終答案的常數項。
但是我們可以證明若 \(\ln f(x)\) 存在，其充要條件是 \(f(x)\) 常數項為 \(1\)。
令 \(g(x) = 1-f(x)\)
考慮：
\[\ln f(x)=\ln (1-g(x))=-\sum_{k=0}^{+\infty} \frac{g^k(x)}{k}
\]
若 \([x^0]g(x)\ne 0\)，則 \(\ln(1-g(x))=-\sum_{k=0}^{+\infty}\frac{g^k(x)}{k}\) 的常數項可能不收斂，如果僅考慮整數的話則一定不收斂。可以求出當且僅當 \([x^0]g(x)\in (0,2)\) 時才收斂。
由此我們可以推得 \([x^0]\ln f(x)=0\)。
多項式 \(\exp\)，自然常數的冪
已知 \(h(x)\)，求 \(f(x)=\exp h(x)\)。
考慮使用 \(\textrm{Newton's Method}\) 解決這個問題。
構造 \(g(f(x))=\ln f(x)-h(x)\)，此時我們把 \(h(x)\) 當作一個常數。
因此有：
\[\begin{aligned}
f(x)&\equiv f_0(x)-\frac{\ln f_0(x)-h(x)}{f_0^{-1}(x)} &\mod x^n \\
&\equiv f_0(x)(1-\ln f_0(x)+h(x)) &\mod x^n
\end{aligned}
\]
時間複雜度 \(O(n\log_2 n)\)。
需要注意的是 \(h(x)\) 的常數項一定為 \(0\)，否則如果使用快速數論變換時會出現問題，因為我們無法得到 \(\exp [x^0]h(x)\) 在取模意義下的取值。
多項式 \(\textrm{power}\)，冪
已知 \(h(x)\)，求 \(h^m(x)\)。
顯然有：
\[h^m(x)=\exp (m \ln h(x))
\]
當然，如果出現 \([x^0]h(x)=0\) 或者 \([x^0]h(x)\ne 1,0\) 的情況要先對 \(h(x)\) 做處理，最後再乘回來。
時間複雜度 \(O(n\log_2 n)\)。
全家桶程式碼
// Not afraid to dark.
#include <bits/stdc++.h>
using namespace std;
clock_t start_time, end_time;
#define GET_START start_time = clock ();
#define GET_END end_time = clock (); fprintf (stderr, "TIME COSSEMED : %0.3lf\n", 1.0 * (end_time - start_time) / CLOCKS_PER_SEC);
#define int long long
namespace io {
int read_pos, read_dt; char read_char;
__inline__ __attribute__ ((always_inline)) int read (int &p = read_pos){
p = 0, read_dt = 1; read_char = getchar ();
while (! isdigit (read_char)){
if (read_char == '-')
read_dt = - 1;
read_char = getchar ();
}
while (isdigit (read_char)){
p = (p << 1) + (p << 3) + read_char - 48;
read_char = getchar ();
}
return p = p * read_dt;
}
int write_sta[65], write_top;
__inline__ __attribute__ ((always_inline)) void write (int x){
if (x < 0)
putchar ('-'), x = - x;
write_top = 0;
do
write_sta[write_top ++] = x % 10, x /= 10;
while (x);
while (write_top)
putchar (write_sta[-- write_top] + 48);
}
int llen;
__inline__ __attribute__ ((always_inline)) int get_string (char c[], int &len = llen){
len = 0;
read_char = getchar ();
while (read_char == ' ' || read_char == '\n' || read_char == '\r')
read_char = getchar ();
while (read_char != ' ' && read_char != '\n' && read_char != '\r'){
c[++ len] = read_char;
read_char = getchar ();
}
return len;
}
}
#define log2_(x) ((x) ? (31 ^ __builtin_clz (x)) : 0)
namespace Region {
const int mod = 998244353, G = 3, Gi = 332748118;
__inline__ __attribute__ ((always_inline)) int power (int a, int p){
int res = 1;
while (p > 0){
if (p & 1)
res = res * a % mod;
a = a * a % mod;
p >>= 1;
}
return res;
}
__inline__ __attribute__ ((always_inline)) int Inverse (int x){
return power (x, mod - 2);
}
}
namespace Polynomial {
// ordinary Number-Theoretic Tranform
using namespace Region;
const int N = 5e5;
__inline__ __attribute__ ((always_inline)) void Reverse (vector < int > &a, int len){
static int rev[N + 5];
rev[0] = 0;
for (int i = 1;i < len;++ i){
rev[i] = rev[i >> 1] >> 1;
if (i & 1)
rev[i] |= len >> 1;
if (i < rev[i])
swap (a[i], a[rev[i]]);
}
}
struct pnm {
int len;
vector < int > a;
pnm (int _len_ = 0, int _c_ = 0){
len = _len_;
a.resize (len, 0);
if (len)
a[0] = _c_;
}
__inline__ __attribute__ ((always_inline)) int & operator [] (const int x) {
return a[x];
}
__inline__ __attribute__ ((always_inline)) int operator [] (const int x) const {
return a[x];
}
__inline__ __attribute__ ((always_inline)) int deg (){
for (int i = len - 1;i > 0;-- i)
if (a[i] != 0)
return i;
return 0;
}
__inline__ __attribute__ ((always_inline)) void flush (int x = 0){
if (x < len)
x = len;
for (int i = 2 << log2_ (x - 1);i > len;-- i)
a.push_back (0);
len = (int) a.size ();
}
__inline__ __attribute__ ((always_inline)) void release (){
len = 0, a.clear ();
}
__inline__ __attribute__ ((always_inline)) pnm Cut (int l) const {
pnm res (l);
for (int i = min (l, len) - 1;i >= 0;-- i)
res[i] = a[i];
return res;
}
__inline__ __attribute__ ((always_inline)) void Delete (int l){
while (len > l){
-- len;
a.pop_back ();
}
}
__inline__ __attribute__ ((always_inline)) void NTT (int swi){
Reverse (a, len);
int om, mid, x, y;
for (int l = 2;l <= len;l <<= 1){
om = power ((swi == - 1) ? Gi : G, (mod - 1) / l), mid = l >> 1;
for (int i = 0;i < len;i += l)
for (int j = i, pw = 1;j < i + mid;++ j, pw = pw * om % mod){
x = a[j], y = a[j + mid] * pw % mod;
a[j] = (x + y) % mod;
a[j + mid] = (x - y + mod) % mod;
}
}
if (swi == - 1){
int iv = Inverse (len);
for (int i = 0;i < len;++ i)
a[i] = a[i] * iv % mod;
}
}
};
__inline__ __attribute__ ((always_inline)) pnm operator + (const pnm a, const pnm b){
pnm res (max (a.len, b.len));
for (int i = 0;i < a.len;++ i)
res[i] = a[i];
for (int i = 0;i < b.len;++ i)
res[i] = (res[i] + b[i]) % mod;
return res;
}
__inline__ __attribute__ ((always_inline)) pnm operator - (const pnm a, const pnm b){
pnm res (max (a.len, b.len));
for (int i = 0;i < a.len;++ i)
res[i] = a[i];
for (int i = 0;i < b.len;++ i)
res[i] = (res[i] - b[i] + mod) % mod;
return res;
}
__inline__ __attribute__ ((always_inline)) pnm operator ^ (const pnm a, const pnm b){
// assert (a.len == b.len);
pnm res = a;
for (int i = 0;i < b.len;++ i)
res[i] = res[i] * b[i];
return res;
}
__inline__ __attribute__ ((always_inline)) pnm operator * (pnm a, pnm b){
static int l = 0;
l = a.len + b.len - 1;
a.flush (l), b.flush (l);
l = a.len;
a.NTT (1), b.NTT (1);
a = a ^ b;
a.NTT (- 1);
return a;
}
__inline__ __attribute__ ((always_inline)) pnm inv (pnm f, int le = 0){
// assert (f[0] != 0);
f.flush ();
if (! le)
le = f.len;
pnm iv (1, Inverse (f[0]));
for (int l = 2;l <= le;l <<= 1)
iv = iv * (pnm (1, 2) - (f.Cut (l) * iv).Cut (l)), iv.Delete (l);
return iv;
}
__inline__ __attribute__ ((always_inline)) pnm operator / (const pnm a, pnm b){
// assert (b[0] != 0);
return a * inv (b);
}
__inline__ __attribute__ ((always_inline)) void Rev (pnm &f, int l = 0){
if (! l)
l = f.len;
f.flush (l);
for (int i = 0;i < l - i - 1;++ i)
swap (f[i], f[l - i - 1]);
}
__inline__ __attribute__ ((always_inline)) pnm operator | (pnm a, pnm b){
int nn = a.deg (), mm = b.deg ();
Rev (a, nn + 1), Rev (b, mm + 1);
pnm c = a.Cut (nn - mm + 1) * inv (b.Cut (nn - mm + 1));
c.Delete (nn - mm + 1);
Rev (c, nn - mm + 1);
return c;
}
__inline__ __attribute__ ((always_inline)) pnm operator % (const pnm a, const pnm b){
return a - ((a | b) * b);
}
__inline__ __attribute__ ((always_inline)) pnm inte (pnm f){
if (f.len == 0)
return f;
for (int i = 1;i < f.len;++ i)
f[i - 1] = f[i] * i % mod;
-- f.len, f.a.pop_back ();
return f;
}
__inline__ __attribute__ ((always_inline)) pnm deri (pnm f, int _c = 0){
static int jc[N + 5], inj;
jc[0] = 1;
for (int i = 1;i <= f.len;++ i)
jc[i] = jc[i - 1] * i % mod;
inj = Inverse (jc[f.len]);
++ f.len;
f.a.push_back (0);
for (int i = f.len - 1;i > 0;inj = inj * (i --) % mod)
f[i] = f[i - 1] * jc[i - 1] % mod * inj % mod;
f[0] = _c;
return f;
}
__inline__ __attribute__ ((always_inline)) pnm LN (pnm f, int m = 0){
// assert (f[0] == 1);
f.flush ();
if (! m)
m = f.len;
return deri (inte (f) * inv (f, m)).Cut (m);
}
__inline__ __attribute__ ((always_inline)) pnm EXP (pnm f, int m = 0){
// assert (f[0] == 0);
f.flush ();
static int p;
pnm g (1, 1);
p = 1;
if (! m)
m = f.len;
while (p < m){
p <<= 1;
g = g * (pnm (1, 1) - LN (g, p) + f.Cut (p));
g.Delete (p);
}
g.Delete (m);
return g;
}
__inline__ __attribute__ ((always_inline)) pnm SQRT (pnm f, int m = 0){
// assert (f[0] == 1)
f.flush ();
static int p;
pnm g (1, 1);
p = 1;
if (! m)
m = f.len;
while (p < m){
p <<= 1;
g = g - ((g * g).Cut (p) - f.Cut (p)) * pnm (1, Inverse (2)) * inv (g, p);
g.Delete (p);
}
g.Delete (m);
return g;
}
__inline__ __attribute__ ((always_inline)) pnm POWER (pnm f, int m, int n = 0){
f.flush ();
if (! n)
n = f.len;
int lowx = 0;
while (lowx < f.len && ! f[lowx])
++ lowx;
if (lowx == f.len)
return pnm (n);
if (lowx && (m >= mod || m * lowx > n))
return pnm (n);
int v = f[lowx];
int iv = Inverse (v);
for (int i = 0;i < f.len;++ i){
if (i + lowx < f.len)
f[i] = f[i + lowx] * iv % mod;
else
f[i] = 0;
}
f = EXP (pnm (1, m % mod) * LN (f), n).Cut (n);
if (lowx){
pnm g (m % mod * lowx + 1);
g[m % mod * lowx] = power (v, m % (mod - 1));
return f * g;
}
return f * pnm (1, power (v, m % (mod - 1)));
}
}
using Polynomial::EXP;
using Polynomial::LN;
using Polynomial::POWER;
using Polynomial::SQRT;
signed main (){
// GET_START
// GET_END
return 0;
}