【清华夏令营模拟2019.5.22】连续段（析合树+多项式牛顿迭代）

Description：

n<=1e5，P是NTT模数

题解：

析合树见WC2019LCA营员交流讲稿

我们考虑把一个序列划分成本源连续段。

怎么划分呢？就是极大划分，假设划分成了x段。

这x段需要满足任取一个区间的段，要么可以全部可以拼起来(1)，要么除了最大的那一个其它都不行(2)。

仔细思考这样就能表示所有的段了。

(1)就是析点，(2)就是合点。

考虑x=2、3时只能是合点

设F表示答案序列的生成函数 所以不难列出这样的方程： $({\sum }_{i=2}^{\infty }{F}^i)\ast 2-F^2-F^3=F-x$ $\frac{F^2}{1-F}-F^2-F^3=F-x$ $F^4+2\ast F^2-(1+x)F+x=0$

那么用牛顿迭代去解方程就行了

Code：

#include<bits/stdc++.h> #define fo(i, x, y) for(int i = x, B = y; i <= B; i ++) #define ff(i, x, y) for(int i = x, B = y; i < B; i ++) #define fd(i, x, y) for(int i = x, B = y; i >= B; i --) #define ul unsigned long long #define ll long long #define pp printf using namespace std; int n, mo, g; ll ksm(ll x, ll y) { ll s = 1; for(; y; y /= 2, x = x * x % mo) if(y & 1) s = s * x % mo; return s; } const int N = 4e5 + 5; int r[N]; void dft(ll *a, int n, int F) { ff(i, 0, n) { r[i] = r[i / 2] / 2 + (i & 1) * (n / 2); if(i < r[i]) swap(a[i], a[r[i]]); } for(int h = 1; h < n; h *= 2) { ll wn = ksm(ksm(g, (mo - 1) / 2 / h), F == 1 ? 1 : mo - 2); for(int j = 0; j < n; j += 2 * h) { ll w = 1, b, *l = a + j, *r = a + j + h; ff(i, 0, h) { b = *r * w, *r = (*l - b) % mo, *l = (*l + b) % mo; w = w * wn % mo, l ++, r ++; } } } if(F == -1) { ll v = ksm(n, mo - 2); ff(i, 0, n) a[i] = (a[i] + mo) * v % mo; } } ll a0[N], a1[N]; typedef vector<ll> V; #define pb push_back #define si size() V operator +(V a, V b) { a.resize(max(a.si, b.si)); ff(i, 0, b.si) a[i] = (a[i] + b[i]) % mo; return a; } V operator -(V a, V b) { a.resize(max(a.si, b.si)); ff(i, 0, b.si) a[i] = (a[i] - b[i] + mo) % mo; return a; } V operator *(V a, ll b) { ff(i, 0, a.si) a[i] = a[i] * b % mo; return a; } V operator *(V a, V b) { int n0 = a.si + b.si - 1, n = 1; while(n < n0) n *= 2; ff(i, 0, n) a0[i] = a1[i] = 0; ff(i, 0, a.si) a0[i] = a[i]; ff(i, 0, b.si) a1[i] = b[i]; dft(a0, n, 1); dft(a1, n, 1); ff(i, 0, n) a0[i] = a0[i] * a1[i] % mo; dft(a0, n, -1); a.resize(n0); ff(i, 0, n0) a[i] = a0[i]; return a; } void dft(V &a, int F) { ff(i, 0, a.si) a0[i] = a[i]; dft(a0, a.si, F); ff(i, 0, a.si) a[i] = a0[i]; } V a, b; V qni(V a) { int n0 = 1; while(n0 < a.si) n0 *= 2; V b; b.resize(1); b[0] = ksm(a[0], mo - 2); for(int n = 2; n <= n0; n *= 2) { V c = a; c.resize(n); c.resize(2 * n); b.resize(2 * n); dft(c, 1); dft(b, 1); ff(i, 0, 2 * n) b[i] = (2 * b[i] - c[i] * b[i] % mo * b[i]) % mo; dft(b, -1); b.resize(n); } b.resize(a.si); return b; } V yy(V a) { fd(i, a.si - 1, 1) a[i] = a[i - 1]; a[0] = 0; return a; } V dd(int n0) { V a; a.resize(1); a[0] = 0; for(int n = 2; n <= n0; n *= 2) { V c; c.resize(2); c[0] = 0; c[1] = 1; V d; d.resize(2); d[0] = -1; d[1] = -1; V b = a; b.resize(n); c = c - b - yy(b); d = d + b * 4; b = b * a; b.resize(n); c = c + b * 2; b = b * a; b.resize(n); d = d + b * 4; b = b * a; b.resize(n); c = c + b; c = c * qni(d); c.resize(n); a.resize(n); a = a - c; } return a; } int main() { freopen("b.in", "r", stdin); freopen("b.out", "w", stdout); scanf("%d %d", &n, &mo); for(g = 2; ; g ++) { if(ksm(g, (mo - 1) / 2) != 1) break; } /* n = 10; a.resize(n + 1); fo(i, 0, n) a[i] = i + 1; b = qni(a); a = a * b; ff(i, 0, n) pp("%lld ", a[i]); pp("\n"); return 0;*/ int n0 = 1; while(n0 <= n) n0 *= 2; a = dd(n0); fo(i, 1, n) pp("%lld\n", a[i]); }