多项式全家桶
多项式全家桶(更新至开根)
- 多项式常用复杂度分析 \(T(n) = T(n/2) + O(nlogn) = O(nlogn)\)
FFT快速傅里叶变换
流程: 将多项式\(\Theta (nlog_n)\)转成点值表示形式 进行卷积, 再\(\Theta (nlog_n)\) 转回来
离散傅里叶变换:
朴素转为点值, 需要将一个一个x带入, 而这里傅里叶搞到了几个可以优化的复数根
利用复数, 在复平面上画出一个单位圆, 将单位圆n等分, 每一个等分点为$\omega_n^i $,
接下来将\(\omega^1_n \omega^2_n\cdots\omega^n_n\) n个根带入, 其中\(\omega_n^i = cos(\frac{k}{n}2\pi) + i * sin(\frac{k}{n}2\pi)\)
稍等, 还要几个小性质才可以:
性质一: \(\omega^{2k}_{2n} = \omega_n^k\)
性质二: \(\omega_n^{k+\frac{n}{2}} = - \omega_n^k\)
接下来开始快速傅里叶变换:
设A(x) = \(a_0 + a_1x + a_2x^2 + \cdots + a_{n-1}x^{n-1}\)
利用分治, 将A按x的指数分为奇偶两部分
\(A(x)=(a_0+a_2x^2+\cdots + a_{n-2}x^{n-2})+(a_1x+a_3x^3+\cdots\) $ a_{n-1}x^{n-1})$
设偶部分 \(A_1(x) = a_0 + a_2x + a_4x^2 + \cdots + a_{n-2}x^{\frac{n}{2}-1}\)
奇部分\(A_2(x) = a_1 + a_3x + a_5x^2 + \cdots + a_{n-1}x^{\frac{n}{2}-1}\)
则有: \(A(x) = A_1(x^2) + x A_2(x^2)\)
设k < n/2 带入\(x = \omega^k_n\)
\(A(\omega^k_n) = A_1(\omega^{2k}_n) + \omega_n^kA_2(\omega^{2k}_n)\)
再带入 \(\omega^{k+\frac{n}{2}}_n\)
得\(A(\omega^{k+\frac{n}{2}}_n) = A_1(\omega^{2k}_n) -\omega_n^kA_2(\omega^{2k}_n)\)
你fa♂现了吗, 分治时, 两个求值一次解决
IDFT: 离散傅里叶逆变换
设\((y_0,y_1,\cdots,y_{n-1})\) 为多项式A(x) = \(a_0 + a_1x + a_2x^2 + \cdots + a_{n-1}x^{n-1}\)的离散傅里叶变换
设多项式 B(x) = A(x) = \(y_0 + y_1x + y_2x^2 + \cdots + y_{n-1}x^{n-1}\)
将n个单位根的倒数带入得新的离散傅里叶变换
此处略去推导, 得\(a_i = \frac{z_i}{n}\)
非递归版fft
每个数递归到底层时的二进制表示恰是原来的表示翻转得到的
如 6(011) 到 3(110)
再加上蝴蝶变换就可以啦
#include<iostream>
#include<cstdio>
#include<cstring>
#include<cmath>
using namespace std;
const int N = 3000600;
template <typename T>
void read(T &x) {
x = 0; bool f = 0;
char c = getchar();
for (;!isdigit(c);c=getchar()) if (c=='-') f=1;
for (;isdigit(c);c=getchar()) x=(x<<3)+(x<<1)+(c^48);
if (f) x = -x;
}
struct Complex {
double x, y;
Complex (double xx = 0,double yy = 0) {x=xx,y=yy;}
Complex operator + (Complex &i) const {
return Complex(x + i.x, y + i.y);
}
Complex operator - (Complex &i) const {
return Complex(x - i.x, y - i.y);
}
Complex operator * (Complex &i) const {
return Complex(x * i.x - y * i.y, x * i.y + y * i.x);
}
}A[N], B[N];
const double Pi = acos(-1.0);
int n, m;
int lim = 1, L;
int r[N];
void FFT(Complex *a,int type) {
for (int i = 0;i < lim; i++)
if (r[i] > i) swap(a[r[i]], a[i]);
for (int mid = 1;mid < lim; mid <<= 1) {
Complex T(cos(Pi/mid), type * sin(Pi/mid));
for (int j = 0;j < lim; j += (mid << 1)) {
Complex t(1, 0);
for (int k = 0;k < mid; k++, t = t * T) {
Complex x = a[j + k], y = t * a[mid + j + k];
a[k + j] = x + y;
a[mid + j + k] = x - y;
}
}
}
}
int main() {
read(n), read(m);
for (int i = 0;i <= n; i++) {
int x; read(x);
A[i].x = x;
}
for (int i = 0;i <= m; i++) {
int x; read(x);
B[i].x = x;
}
while (lim <= (n + m)) lim <<= 1, L++;
for (int i = 0;i < lim; i++)
r[i] = (r[i>>1]>>1) | ((i&1) << (L-1));
FFT(A, 1); FFT(B, 1);
for (int i = 0;i <= lim; i++) A[i] = A[i] * B[i];
FFT(A, -1);
for (int i = 0;i <= n + m; i++)
printf ("%d ", (int)(A[i].x / lim + 0.5));
return 0;
}
NTT快速数论变换
直接上代码
#include<iostream>
#include<cstdio>
#include<cstring>
#define ll long long
using namespace std;
const int N = 3e6+6;
const int P = 998244353;
const int G = 3;
const int Gi = (P+1)/G;
inline int read(void) {
int x = 0, f = 1;
char c = getchar();
for (;!isdigit(c);c=getchar()) if (c=='-') f=-1;
for (;isdigit(c);c=getchar()) x=(x<<3)+(x<<1)+c-'0';
return x * f;
}
int n, m, r[N], L;
ll A[N], B[N];
int lim = 1;
ll pow(ll di,ll mi) {
ll ans = 1, a = di;
while (mi) {
if (mi & 1) ans = ans * a % P;
a = a * a % P;
mi >>= 1;
}
return ans;
}
void NTT(ll *a, int tag) {
for (int i = 0;i < lim; i++)
if (i < r[i]) swap(a[i], a[r[i]]);
for (int mid = 1;mid < lim ;mid <<= 1) {
ll wn = pow(tag == 1 ? G : Gi, (P-1) / (mid << 1));
for (int j = 0;j < lim; j += (mid << 1)) {
ll w = 1;
for (int k = 0;k < mid ; k++, w = (w * wn) % P) {
ll x = a[j+k], y = w * a[j+k+mid] % P;
a[j+k] = (x + y) % P;
a[j+k+mid] = (x - y + P) % P;
}
}
}
}
int main() {
n = read(), m = read();
for (int i = 0;i <= n; i++) A[i] = read();
for (int j = 0;j <= m; j++) B[j] = read();
while (lim <= n + m) L++, lim <<= 1;
for (int i = 0;i < lim; i++)
r[i] = (r[i>>1]>>1) | ((i&1) << (L-1));
NTT(A, 1), NTT(B, 1);
for (int i = 0;i <= lim; i++) A[i] = (A[i] * B[i]) % P;
NTT(A, -1);
ll inv = pow(lim, P-2);
for (int i = 0;i <= n + m; i++)
printf ("%d ", A[i] * inv % P);
return 0;
}
多项式求逆
\[\text{求 }g \equiv f^{-1} \mod x^n \\
\text{设 }h*f \equiv 1 \mod x^{\lceil \frac n2 \rceil}\\
(g - h)*f \equiv 0 \mod x^{\lceil \frac n2 \rceil}\\
(g-h)^2 \equiv 0 \mod x^n\\
f*(g-h)^2 \equiv 0 \mod x^n\\
g-2*h+h^2*f \equiv 0 \mod x^n\\
g\equiv h*(2-h*f) \mod x^n\\
\]
#include<iostream>
#include<cstdio>
#include<cstring>
#define ll long long
using namespace std;
const int P = 998244353;
const int G = 3;
const int Gi = (P+1)/G;
const int N = 3e5+5;
template <typename T>
void read(T &x) {
x = 0; int f = 1;
char c = getchar();
for (;!isdigit(c);c=getchar()) if (c=='-') f=-1;
for (;isdigit(c);c=getchar()) x=(x<<3)+(x<<1) + c-'0';
x *= f;
}
ll A[N], n, lim = 1, L;
int r[N];
ll pow(ll di,ll mi) {
ll ans = 1;
while (mi) {
if (mi & 1) ans = (ans * di) % P;
di = (di * di) % P;
mi >>= 1;
}
return ans;
}
void NTT(ll *a, int type) {
for (int i = 0;i < lim; i++)
if (i < r[i]) swap(a[i], a[r[i]]);
for (int mid = 1;mid < lim; mid <<= 1) {
ll wn = pow(type == 1 ? G : Gi, (P-1)/(mid<<1));
for (int j = 0;j < lim;j += (mid << 1)) {
ll w = 1;
for (int k = 0;k < mid; k++, w = (w * wn) % P) {
ll x = a[j+k], y = w * a[j+k+mid] % P;
a[j+k] = (x + y) % P, a[j+k+mid] = (x - y + P) % P;
}
}
}
if (type == 1) return;
ll inv = pow(lim, P-2);
for (int i = 0;i < lim; i++)
a[i] = a[i] * inv % P;
}
ll B[N], C[N], D[N];
void Inv(int deg, ll *a, ll *b) {
b[0] = pow(a[0], P - 2);
int len = 1;
for ( ; len <= (n << 1); len <<= 1) {
lim = len << 1;
for (int i = 1;i < lim; i++)
r[i] = (r[i>>1]>>1) | ((i & 1) ? len : 0);
for (int i = 0;i < len; i++) C[i] = a[i], D[i] = b[i];
NTT(C, 1), NTT(D, 1);
for (int i = 0;i < lim; i++) b[i] = (2 - C[i] * D[i] % P + P) * D[i] % P;
NTT(b, -1);
for (int i = len;i < lim; i++) b[i] = 0;
}
}
int main() {
read(n);
for (int i = 0;i < n; i++) read(A[i]);
Inv(n, A, B);
for (int i = 0;i < n; i++) printf ("%lld ", B[i]);
return 0;
}
分治FFT
应用CDQ分治思想即可
#include<iostream>
#include<cstdio>
#include<cstring>
#define ll long long
using namespace std;
const int P = 998244353;
const int G = 3, Gi = (P+1) / G;
const int N = 405000;
ll fpw(ll di, ll mi) {
ll res = 1;
while (mi) {
if (mi & 1) res = res * di % P;
di = di * di % P;
mi >>= 1;
}
return res;
}
int lim, L;
ll r[N];
void NTT(ll *a,int type) {
for (int i = 0;i < lim; i++)
if (i < r[i]) swap(a[i], a[r[i]]);
for (int mid = 1;mid < lim; mid <<= 1) {
ll wn = fpw(type == 1 ? G : Gi, (P-1) / (mid<<1));
for (int i = 0;i < lim; i += mid << 1) {
ll w = 1;
for (int j = 0;j < mid; j++, w = w * wn % P) {
ll x = a[j + i], y = a[mid + i + j] * w % P;
a[i + j] = (x + y) % P;
a[i + j + mid] = (x - y + P) % P;
}
}
}
if (type == 1) return;
ll inv = fpw(lim, P-2);
for (int i = 0;i < lim; i++) a[i] = a[i] * inv % P;
}
ll a[N], f[N], b[N], g[N];
void solve(int l,int R) {
if (l == R) return;
int mid = (l + R) >> 1;
solve(l, mid); lim = 1, L = 0;
int len = R - l + 1;
for (int i = l;i <= mid; i++) a[i-l] = f[i];
for (int i = 0;i < len; i++) b[i] = g[i];
len += mid - l + 1;
while (lim <= len) lim <<= 1, L++;
for (int i = 0;i < lim; i++)
r[i] = (r[i>>1]>>1) | ((i&1)<<(L-1));
NTT(a, 1), NTT(b, 1);
for (int i = 0;i < lim; i++) a[i] = a[i] * b[i] % P;
NTT(a, -1);
for (int i = mid + 1;i <= R; i++) f[i] = (f[i] + a[i-l]) % P;
for (int i = 0;i < lim; i++) a[i] = b[i] = 0;
solve(mid + 1, R);
}
int read(void) {
int x = 0; char c = getchar();
while (!isdigit(c)) c=getchar();
for (;isdigit(c);c=getchar()) x=(x<<3)+(x<<1)+(c^48);
return x;
}
int n;
int main() {
n = read();
for (int i = 1;i < n; i++) g[i] = read();
f[0] = 1;
solve(0, n - 1);
for (int i = 0;i < n; i++) printf ("%lld ", f[i]);
return 0;
}
多项式ln：
\[G \equiv ln F\\
G' \equiv \frac {F'}{F}
\]
求逆再积分即可
#include<iostream>
#include<cstring>
#include<cstdio>
#define ll long long
using namespace std;
const int N = 505000;
const int P = 998244353;
const int G = 3;
const int Gi = (P+1) / G;
int read(void) {
int x = 0; bool f = 0;
char c = getchar();
for (;!isdigit(c);c=getchar()) if (c=='-') f=1;
for (;isdigit(c);c=getchar()) x=(x<<3)+(x<<1)+(c^48);
return f ? -x : x;
}
ll fpw(ll di,ll mi) {
ll ans = 1;
while (mi) {
if (mi & 1) ans = ans * di % P;
di = di * di % P;
mi >>= 1;
}
return ans;
}
ll A[N], B[N], C[N];
int r[N], n, lim, L;
void NTT(ll *a, int type) {
for (int i = 0;i < lim; i++)
if (i < r[i]) swap(a[i], a[r[i]]);
for (int j = 1;j < lim; j <<= 1) {
ll wn = fpw(type == 1 ? G : Gi, (P-1) / (j<<1));
for (int k = 0;k < lim; k += j << 1) {
ll w = 1;
for (int p = 0;p < j; p++, w = w * wn % P) {
ll x = a[p+k], y = a[p+k+j] * w % P;
a[p+k] = (x + y) % P;
a[p+k+j] = (x - y + P) % P;
}
}
}
if (~type) return;
ll inv = fpw(lim, P-2);
for (int i = 0;i < lim; i++)
a[i] = a[i] * inv % P;
}
void work(int deg, ll *a, ll *b) {
if (deg == 1) {
b[0] = fpw(a[0], P-2);
return;
}
work((deg+1)>>1, a, b);
lim = 1, L = 0;
while (lim < (deg<<1)) lim <<= 1, L++;
for (int i = 0;i < lim; i++)
r[i] = (r[i>>1]>>1) | ((i&1)<<(L-1));
for (int i = 0;i < deg; i++) C[i] = a[i];
for (int i = deg;i < lim; i++) C[i] = 0;
NTT(C, 1); NTT(b, 1);
for (int i = 0;i < lim; i++)
b[i] = ((ll)2 - C[i] * b[i] % P + P) % P * b[i] % P;
NTT(b, -1);
for (int i = deg;i < lim; i++) b[i] = 0;
}
int m;
void qiudoor(void) {
for (int i = 1;i < n; i++) A[i-1] = A[i] * i % P;
A[n-1] = 0;
}
void jinitaimei(void) {
for (int i = n-1;i >= 1; i--) A[i] = A[i-1] * fpw(i, P-2) % P;
A[0] = 0;
}
int main() {
n = read();
for (int i = 0;i < n; i++) A[i] = read();
work(n, A, B); qiudoor();
lim = 1, L = 0; m = n;
while (lim <= n + m) lim <<= 1, L++;
NTT(A, 1); NTT(B, 1);
for (int i = 0;i < lim; i++) A[i] = A[i] * B[i] % P;
NTT(A, -1);
jinitaimei();
for (int i = 0;i < n; i++) printf ("%lld ", A[i] % P);
return 0;
}
多项式exp
\[牛顿迭代\\
F(G_0) = 0 \mod {x^{\lceil \frac n2 \rceil}}\\
F(G)=F(G_0) + \frac{F'(G_0)}{1!}(G - G_0) + \frac{F''(G_0)}{2!}(G-G_0)^2+\cdots \mod {x^n}\\
易知G - G_0 = 0\mod {x^{\lceil \frac n2 \rceil}}\\
因为后面二次三次方的系数已经大于x^n, 所以可以忽略\\
又 \because F(G) = 0 \mod x^n\\
\therefore G = G_0 - \frac{F(G_0)}{F'(G_0)}
\]
下面回到本题
\[e^F = H\\
F = lnH\\
ln H-F = 0\\
设G = ln H-F, 即求G的零点\\
以H为变量,F为常量求导\\
G'=\frac {1}H\\
H = H_0-\frac{lnH_0-F}{\frac {1}H_0}=H_0(1-lnH_0+F)
\]
#include <iostream>
#include <cstdio>
#include <cstring>
#define ll long long
using namespace std;
const int P = 998244353;
const int G = 3, Gi = (P+1) / G;
const int N = 400005;
template <typename T>
void read(T &x) {
x = 0; int f = 0; char c = getchar();
for (;!isdigit(c);c=getchar()) if (c=='-') f=1;
for (;isdigit(c);c=getchar()) x=(x<<3)+(x<<1)+(c^48);
if (f) x=-x;
}
ll fpw(ll di, ll mi) {
ll res = 1;
while (mi) {
if (mi & 1) res = res * di % P;
di = di * di % P;
mi >>= 1;
}
return res;
}
int lim, r[N], n, L;
ll A[N], B[N], C[N], D[N], E[N];
void NTT(ll *a, int type = 1) {
for (int i = 0;i < lim; i++)
if (r[i] > i) swap(a[i], a[r[i]]);
for (int mid = 1;mid < lim; mid <<= 1) {
ll p = mid << 1, wn = fpw(type == 1 ? G : Gi, (P-1) / (mid << 1));
for (int j = 0;j < lim; j += p) {
ll w = 1;
for (int k = 0;k < mid; k++, w = (w * wn) % P) {
ll x = a[j+k], y = a[j+k+mid] * w % P;
a[j+k] = (x + y) % P, a[j+k+mid] = (x - y + P) % P;
}
}
}
if (type == 1) return;
ll inv = fpw(lim, P - 2);
for (int i = 0;i < lim; i++) a[i] = a[i] * inv % P;
}
void Inv(int deg, ll *a, ll *b) {
b[0] = fpw(a[0], P-2); int len;
for (len = 1;len < (deg << 1); len <<= 1) {
lim = len << 1;
for (int i = 0;i < len; i++) A[i] = a[i], B[i] = b[i];
for (int i = 0;i < lim; i++)
r[i] = (r[i>>1]>>1) | ((i&1)?len:0);
NTT(A, 1); NTT(B, 1);
for (int i = 0;i < lim; i++)
b[i] = (2 - A[i] * B[i] % P + P) * B[i] % P;
NTT(b, -1);
for (int i = len;i < lim; i++) b[i] = 0;
}
for (int i = 0;i < len; i++) A[i] = B[i] = 0;
}
ll inv[N];
void init(void) {
inv[1] = 1;
for (int i = 2;i <= 200080; i++)
inv[i] = P - P / i * inv[P%i] % P;
}
void jinitaimei(ll *a, int deg) {
for (int i = deg - 1;i >= 1; i--)
a[i] = a[i-1] * inv[i] % P;
a[0] = 0;
}
void qiudoor(ll *a, int deg) {
for (int i = 1;i < deg; i++) a[i-1] = a[i] * i % P;
a[deg-1] = 0;
}
void get_ln(ll *a, ll deg) {
Inv(deg, a, C); qiudoor(a, deg);
NTT(a, 1); NTT(C, 1);
for (int i = 0;i < lim; i++) a[i] = a[i] * C[i] % P;
NTT(a, -1); jinitaimei(a, deg);
}
void get_exp(int deg, ll *a, ll *b) {
b[0] = 1; int len;
for (len = 2;len < (deg << 1); len <<= 1) {
lim = len << 1;
for (int i = 0;i < len; i++) E[i] = b[i];
get_ln(E, len);
E[0] = (a[0] + 1 - E[0] + P) % P;
for (int i = 1;i < len; i++)
E[i] = (a[i] - E[i] + P) % P;
NTT(E, 1); NTT(b, 1);
for (int i = 0;i < lim; i++) b[i] = b[i] * E[i] % P;
NTT(b, -1);
for (int i = len;i < lim; i++) b[i] = 0;
}
}
ll f[N], g[N];
int main() {
read(n); init();
for (int i = 0;i < n; i++) read(f[i]);
get_exp(n, f, g);
for (int i = 0;i < n; i++) printf ("%lld ", g[i]);
return 0;
}
多项式开根
方法一(大常数) :
\[G^2 = F\\
2*lnG=lnF\\
G = e^{\frac {lnF}2}
\]
方法二(牛顿迭代):
\[G^2-F=0\\
G=G_0-\frac{(G_0^2-F)}{2G_0}=\frac{G_0^2+F}{2G_0}
\]
方法三(与牛顿迭代形式相同):
\[G_0^2=F \mod {x^{\lceil \frac n2 \rceil}}\\
G^2=F \mod x^n\\
G^2-G_0^2=0 \mod {x^{\lceil \frac n2 \rceil}}\\
\because G + G_0!= 0\\
\therefore G - G_0= 0\\
(G-G_0)^2=0 \mod x^n\\
G^2-2*G*G_0+G_0^2=0\\
F-2*G*G_0+G_0^2=0\\
G=\frac{F+G_0^2}{2*G_0}
\]
#include<iostream>
#include<cstdio>
#include<cstring>
#define ll long long
using namespace std;
const int P = 998244353;
const int G = 3, Gi = (P+1) / G;
const int INV2 = (P+1) / 2;
const int N = 300500;
int read(void) {
int x = 0, f = 0;
char c = getchar();
for (;!isdigit(c);c=getchar()) if (c=='-') f=1;
for (;isdigit(c);c=getchar()) x=(x<<3)+(x<<1)+(c^48);
return f ? -x : x;
}
ll fpw(ll di, ll mi) {
ll res = 1;
while (mi) {
if (mi & 1) res = res * di % P;
di = di * di % P;
mi >>= 1;
}
return res;
}
int lim, r[N], n, L;
ll f[N], g[N], A[N], C[N], D[N];
void NTT(ll *a, int type) {
for (int i = 0;i < lim; i++)
if (i < r[i]) swap(a[i], a[r[i]]);
for (int mid = 1;mid < lim; mid <<= 1) {
ll wn = fpw(type == 1 ? G : Gi, (P-1)/(mid<<1));
for (int j = 0;j < lim;j += (mid << 1)) {
ll w = 1;
for (int k = 0;k < mid; k++, w = (w * wn) % P) {
ll x = a[j+k], y = w * a[j+k+mid] % P;
a[j+k] = (x + y) % P, a[j+k+mid] = (x - y + P) % P;
}
}
}
if (type == 1) return;
ll inv = fpw(lim, P-2);
for (int i = 0;i < lim; i++)
a[i] = a[i] * inv % P;
}
ll B[N];
void Inv(int deg, ll *a, ll *b) {
b[0] = fpw(a[0], P-2);
int len;
for (len = 1;len < (deg<<1);len <<= 1) {
lim = len << 1;
for (int i = 0;i < len; i++) A[i] = a[i], B[i] = b[i];
for (int i = 0;i < lim; i++)
r[i] = (r[i>>1]>>1) | ((i&1)?len:0);
NTT(A, 1), NTT(B, 1);
for (int i = 0;i < lim; i++)
b[i] = (2 - A[i] * B[i] % P + P) * B[i] % P;
NTT(b, -1);
for (int i = len;i < lim; i++) b[i] = 0;
}
for (int i = 0;i < len; i++) A[i] = B[i] = 0;
for (int i = deg;i < len; i++) b[i] = 0;
}
void Sqrt(int deg, ll *a, ll *b) {
b[0] = 1;
ll len;
for (len = 1;len < (deg << 1);len <<= 1) {
lim = len << 1;
for (int i = 0;i < len; i++) C[i] = a[i];
Inv(len, b, D);
for (int i = 0;i < lim; i++)
r[i] = (r[i>>1]>>1) | ((i&1) ? len : 0);
NTT(C, 1), NTT(D, 1);
for (int i = 0;i < lim; i++) C[i] = C[i] * D[i] % P;
NTT(C, -1);
for (int i = 0;i < len; i++) b[i] = (b[i] + C[i]) * INV2 % P;
for (int i = len;i < lim; i++) b[i] = 0;
}
}
int main() {
n = read();
for (int i = 0;i < n; i++) f[i] = read();
Sqrt(n, f, g);
for (int i = 0;i < n; i++) printf ("%lld ", g[i]);
cout << endl;
return 0;
}
多项式除法
\[F=Q*G+R (已知G, F)\\
reverseA(x) = A(\frac 1x)*x^n = rA\\
F(\frac 1x)*x^n=Q(\frac 1x)*G(\frac 1x)*x^n+R(\frac 1x)*x^n\\
rF=rQ*rG+rR*x^{n-m+1} \mod x^{n-m+1}\\
rQ=rF*rG^{-1}\\
R = F-Q*G
\]
注意因为mod的是\(x^{n-m+1}\)大于Q的次数, 所以Q唯一
#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#define ll long long
using namespace std;
const int P = 998244353;
const int G = 3;
const int Gi = (P+1)/G;
const int N = 3e5+5;
template <typename T>
void read(T &x) {
x = 0; int f = 1;
char c = getchar();
for (;!isdigit(c);c=getchar()) if (c=='-') f=-1;
for (;isdigit(c);c=getchar()) x=(x<<3)+(x<<1) + c-'0';
x *= f;
}
ll A[N], n, lim = 1, L;
int r[N];
ll pow(ll di,ll mi) {
ll ans = 1;
while (mi) {
if (mi & 1) ans = (ans * di) % P;
di = (di * di) % P;
mi >>= 1;
}
return ans;
}
void NTT(ll *a, int type) {
for (int i = 0;i < lim; i++)
if (i < r[i]) swap(a[i], a[r[i]]);
for (int mid = 1;mid < lim; mid <<= 1) {
ll wn = pow(type == 1 ? G : Gi, (P-1)/(mid<<1));
for (int j = 0;j < lim;j += (mid << 1)) {
ll w = 1;
for (int k = 0;k < mid; k++, w = (w * wn) % P) {
ll x = a[j+k], y = w * a[j+k+mid] % P;
a[j+k] = (x + y) % P, a[j+k+mid] = (x - y + P) % P;
}
}
}
if (type == 1) return;
ll inv = pow(lim, P-2);
for (int i = 0;i < lim; i++)
a[i] = a[i] * inv % P;
}
ll B[N], C[N];
void work(int deg, ll *a, ll *b) {
if (deg == 1) {
b[0] = pow(a[0], P-2);
return;
}
work((deg+1)>>1, a, b);
lim = 1, L = 0;
while (lim < (deg<<1)) lim <<= 1, L++;
for (int i = 0;i < lim; i++)
r[i] = (r[i>>1]>>1) | ((i&1) << (L-1));
for (int i = 0;i < deg; i++) C[i] = a[i];
for (int i = deg;i < lim; i++) C[i] = 0;
NTT(C, 1), NTT(b, 1);
for (int i = 0;i < lim; i++)
b[i] = ((ll)2 - C[i] * b[i] % P + P) % P * b[i] % P;
NTT(b, -1);
for (int i = deg;i < lim; i++) b[i] = 0;
}
ll f[N], g[N], q[N], ff[N], gg[N];
int m;
int main() {
freopen("hs.in","r",stdin);
freopen("hs.out","w",stdout);
read(n), read(m);
for (int i = n;i >= 0; i--) read(f[i]), ff[n-i] = f[i];
for (int j = m;j >= 0; j--) read(g[j]), gg[m-j] = g[j];
for (int i = n - m + 1;i <= n; i++) f[i] = g[i] = 0;
work(n-m+1, g, B);
NTT(B, 1); NTT(f, 1);
for (int i = 0;i < lim; i++) q[i] = f[i] * B[i] % P;
NTT(q, -1);
reverse(q, q + n - m + 1);
for (int i = 0;i <= n - m; i++) printf ("%lld ", q[i]);
lim = 1, L = 0;
while (lim <= n) lim <<= 1, L++;
for (int i = n - m + 1;i < lim; i++) q[i] = 0;
for (int i = 0;i < lim; i++) r[i] = (r[i>>1]>>1) | ((i&1)<<(L-1));
NTT(gg, 1); NTT(q, 1);
for (int i = 0;i < lim; i++) gg[i] = gg[i] * q[i] % P;
NTT(gg, -1);
putchar('\n');
for (int i = 0;i < m; i++) printf ("%lld ", (ff[i] - gg[i] + P) % P);
return 0;
}
拉格朗日插值
ll fpw(ll x,ll mi) {
ll res = 1;
while (mi) {
if (mi & 1) res = res * x % P;
x = x * x % P;
mi >>= 1;
}
return res;
}
int x[N], y[N], n, k;
int main() {
read(n), read(k);
for (int i = 1;i <= n; i++) read(x[i]), read(y[i]);
ll res = 0;
for (int i = 1;i <= n; i++) {
ll ans1 = 1, ans2 = 1;
for (int j = 1;j <= n; j++) {
if (j == i) continue;
ans1 = ans1 * (k - x[j]) % P;
ans2 = ans2 * (x[i] - x[j]) % P;
}
res += y[i] * ans1 % P * fpw(ans2, P-2) % P;
}
res %= P;
if (res < 0) res += P;
cout << res << endl;
return 0;
}
多项式快速幂
\[F^k=G\\
klnF=lnG\\
e^{klnF}=G
\]
#include<iostream>
#include<cstdio>
#include<cstring>
#define ll long long
using namespace std;
const int P = 998244353;
const int G = 3, Gi = (P+1) / G;
const int N = 1005000;
int read(void) {
int x = 0; char c = getchar();
while (!isdigit(c)) c=getchar();
for (;isdigit(c);c=getchar()) x=(x<<3)+(x<<1)+(c^48);
return x;
}
ll fpw(ll di, ll mi) {
ll res = 1;
while (mi) {
if (mi & 1) res = res * di % P;
di = di * di % P;
mi >>= 1;
}
return res;
}
int lim, L;
int r[N];
void NTT(ll *a,int type) {
for (int i = 0;i < lim; i++)
if (r[i] > i) swap(a[r[i]], a[i]);
for (int mid = 1;mid < lim; mid <<= 1) {
ll wn = fpw(type == 1 ? G : Gi, (P-1) / (mid << 1));
for (int i = 0;i < lim; i += mid << 1) {
ll w = 1;
for (int j = 0;j < mid; j++, w = w * wn % P) {
int x = a[i + j], y = a[i + j + mid] * w % P;
a[i + j] = (x + y) % P;
a[i + j + mid] = (x - y + P) % P;
}
}
}
if (type == 1) return;
ll inv = fpw(lim, P-2);
for (int i = 0;i < lim; i++) a[i] = a[i] * inv % P;
}
ll a[N], b[N], c[N];
void work(int deg, ll *a, ll *b) {
if (deg == 1) {
b[0] = fpw(a[0], P-2);
return;
}
work((deg+1) >> 1, a, b);
lim = 1, L = 0;
while (lim < (deg << 1)) lim <<= 1, L++;
for (int i = 0;i < lim ;i++)
r[i] = (r[i>>1]>>1) | ((i&1) << (L-1));
for (int i = 0;i < deg; i++) c[i] = a[i];
for (int i = deg; i < lim; i++) c[i] = 0;
NTT(b, 1); NTT(c, 1);
for (int i = 0;i < lim; i++)
b[i] = (2 - b[i] * c[i] % P + P) % P * b[i] % P;
NTT(b, -1);
for (int i = deg;i < lim; i++) b[i] = 0;
}
int n;
void qiudoor(ll *a, int deg) {
for (int i = 1;i < n; i++) a[i-1] = a[i] * i % P;
a[n-1] = 0;
}
ll inv[N];
void init(void) {
inv[1] = 1;
for (int i = 2;i <= n * 2; i++)
inv[i] = (P - P / i * inv[P%i] % P) % P;
}
void jinitaimei(ll *a, int deg) {
for (int i = n-1;i >= 1; i--) a[i] = a[i-1] * inv[i] % P;
a[0] = 0;
}
ll A[N], B[N];
void get_ln(ll *a, ll n) {
for (int i = 0;i < n*2; i++) B[i] = 0;
// put(A, n), put(B, n*2);
work(n, a, B);
qiudoor(a, n);
NTT(a, 1); NTT(B, 1);
for (int i = 0;i < lim; i++) a[i] = a[i] * B[i] % P;
NTT(a, -1); jinitaimei(a, n);
}
void solve(int deg, ll *a, ll *b) {
if (deg == 1) return (void) (b[0] = 1);
solve(deg>>1, a, b);
for (int i = 0;i < deg; i++) A[i] = b[i];
get_ln(A, deg);
A[0] = (a[0] + 1 - A[0] + P) % P;
for (int i = 1;i < deg; i++)
A[i] = (a[i] - A[i] + P) % P;
NTT(A, 1); NTT(b, 1);
for (int i = 0;i < lim; i++)
b[i] = b[i] * A[i] % P;
NTT(b, -1);
for (int i = deg;i < lim; i++) b[i] = 0;
}
ll f[N], g[N];
ll get_k(void) {
ll x = 0; char c = getchar();
while (!isdigit(c)) c=getchar();
for (;isdigit(c);c=getchar()) {
x=(x<<3)+(x<<1)+(c^48);
if (x >= P) x %= P;
}
return x;
}
ll k;
int main() {
n = read(); init(); k = get_k();
for (int i = 0;i < n; i++) f[i] = read();
get_ln(f, n);
for (int i = 0;i < n; i++) f[i] = f[i] * k % P;
lim = 1;
while (lim <= n) lim <<= 1;
solve(lim, f, g);
for (int i = 0;i < n; i++)
printf ("%d ", g[i]);
return 0;
}
多项式三角函数: 直接欧拉公式即可
\[e^{i\theta}=sin\theta+cosx\theta
\]
#include <iostream>
#include <cstring>
#include <cstdio>
#define ll long long
using namespace std;
const int P = 998244353;
const int I = 86583718;
const int G = 3, Gi = (P + 1) / G;
const int N = 300500;
template <typename T>
void read(T &x) {
x = 0; int f = 0; char c = getchar();
for (;!isdigit(c);c=getchar()) if (c=='-') f=1;
for (;isdigit(c);c=getchar()) x=(x<<3)+(x<<1)+(c^48);
if (f) x=-x;
}
ll fpw(ll di, ll mi) {
ll res = 1;
while (mi) {
if (mi & 1) res = res * di % P;
di = di * di % P;
mi >>= 1;
}
return res;
}
ll lim, r[N], n, L;
ll A[N], B[N], C[N], D[N], E[N];
void NTT(ll *a, int type = 1) {
for (int i = 0;i < lim; i++)
if (r[i] > i) swap(a[i], a[r[i]]);
for (int mid = 1;mid < lim; mid <<= 1) {
ll p = mid << 1, wn = fpw(type == 1 ? G : Gi, (P-1) / (mid << 1));
for (int j = 0;j < lim; j += p) {
ll w = 1;
for (int k = 0;k < mid; k++, w = (w * wn) % P) {
ll x = a[j+k], y = a[j+k+mid] * w % P;
a[j+k] = (x + y) % P, a[j+k+mid] = (x - y + P) % P;
}
}
}
if (type == 1) return;
ll inv = fpw(lim, P - 2);
for (int i = 0;i < lim; i++) a[i] = a[i] * inv % P;
}
void Inv(int deg, ll *a, ll *b) {
b[0] = fpw(a[0], P-2); int len;
for (len = 1;len < (deg << 1); len <<= 1) {
lim = len << 1;
for (int i = 0;i < len; i++) A[i] = a[i], B[i] = b[i];
for (int i = 0;i < lim; i++)
r[i] = (r[i>>1]>>1) | ((i&1)?len:0);
NTT(A, 1); NTT(B, 1);
for (int i = 0;i < lim; i++)
b[i] = (2 - A[i] * B[i] % P + P) * B[i] % P;
NTT(b, -1);
for (int i = len;i < lim; i++) b[i] = 0;
}
for (int i = 0;i < len; i++) A[i] = B[i] = 0;
}
ll inv[N];
void init(void) {
inv[1] = 1;
for (int i = 2;i <= 100500; i++)
inv[i] = (P - P / i) * inv[P % i] % P;
}
void jinitaimei(ll *a, int deg) {
for (int i = deg - 1;i >= 1; i--)
a[i] = a[i-1] * inv[i] % P;
a[0] = 0;
}
void qiudoor(ll *a, int deg) {
for (int i = 1;i < deg; i++) a[i-1] = a[i] * i % P;
a[deg-1] = 0;
}
void get_ln(ll *a, ll deg) {
Inv(deg, a, C); qiudoor(a, deg);
NTT(a, 1); NTT(C, 1);
for (int i = 0;i < lim; i++) a[i] = a[i] * C[i] % P;
NTT(a, -1); jinitaimei(a, deg);
}
void get_exp(int deg, ll *a, ll *b) {
b[0] = 1; int len;
for (len = 2;len < (deg << 1); len <<= 1) {
lim = len << 1;
for (int i = 0;i < len; i++) E[i] = b[i];
get_ln(E, len);
E[0] = (a[0] + 1 - E[0] + P) % P;
for (int i = 1;i < len; i++)
E[i] = (a[i] - E[i] + P) % P;
NTT(E, 1); NTT(b, 1);
for (int i = 0;i < lim; i++) b[i] = b[i] * E[i] % P;
NTT(b, -1);
for (int i = len;i < lim; i++) b[i] = 0;
}
}
ll f[N], t1[N], t2[N], type;
int main() {
read(n), read(type); init();
for (int i = 0;i < n; i++) read(f[i]), f[i] = f[i] * I % P;
get_exp(n, f, t1); Inv(n, t1, t2);
ll INV = fpw(2 * I % P, P - 2);
for (int i = 0;i < n; i++)
if (type) printf ("%lld ", (t1[i] + t2[i]) * inv[2] % P);
else printf ("%lld ", (t1[i] - t2[i] + P) * INV % P);
return 0;
}
