多项式

算法

一些常用函数导数:
\((sin(x))'\ =\ cos(x)\)
\((cos(x))'\ =\ -sin(x)\)
\((a^x)'\ =\ a^x\ ln(a)\)
\(log_a(x)\ =\ \frac{1}{x\ ln(a)}\)
\((uv)'\ =\ u'v\ +\ uv'\)
\((\frac{u}{v})'\ =\ \frac{u'v\ -\ uv'}{v^2}\)
\(f(g(x))'\ =\ f'(g(x))\ g'(x)\)
求逆:
\(A(x)B(x)\ =\ 0\ \bmod\ x^n\)
\(B(x)\ -\ B_0(x)\ =\ 0\ \bmod\ x^{\frac{n}{2}}\)
\(B(x)\ =\ 2B_0(x)\ -\ B_0^2(x)\ A(x)\ \bmod\ x^n\)
除法:
\(A(x)\ =\ B(x)D(x)\ +\ R(x)\)
\(x^{n\ -\ 1}A(x^{-1})\ =\ x^{m\ -\ 1}\ B(x^{-1})x^{n\ -\ m}\ D(x^{-1})\ +\ x^{m\ -\ 2}\ R(x^{-1})\ x^{n\ -\ m\ +\ 1}\)
\(rA(x)\ =\ rB(x)\ =\ rD(x)\ +\ rR(x)\ x^{n\ -\ m\ +\ 1},\ C(x)\ rB(x)\ =\ 1\ \bmod\ x^{n\ -\ m\ +\ 1}\)
\(rD(x)\ =\ C(x)\ rA(x)\ \bmod\ x^{n\ -\ m\ +\ 1}\)
ln:
\(ln(1\ -\ A(x))\ =\ -\sum_{i\ \geq\ 1}\ \frac{A^i(x)}{i}\)
\(B(x)\ =\ ln(A(x))\)
\(B(x)'\ =\ (ln(A(x)))'\ =\ ln'(A(x))\ A(x)\ =\ \frac{A'(x)}{A(x)}\)
牛顿迭代:
\(F(x)\ =\ F_0(x)\ -\ \frac{G(F_0(x))}{G'(F_0(x))}\ \bmod\ x^n\)
exp:
\(B(x)\ =\ e^A(x)\ \bmod\ x^n\)
\(ln(B(x))\ =\ A(x)\ \bmod\ x^n\)
\(B(x)\ =\ B_0(x)(1\ -\ ln(B_0(x))\ +\ A(x))\ \bmod\ x^n\)

应用

例1.BZOJ3456 对于一种划分集合的方案,另 \(G(x)\) 为被划分的集合,\(F(x)\) 为划分的集合,有 \(e^{F(x)}\ =\ G(x)\)
例2.多项式乘方 \(A(x)^N\ =\ e^{N\ ln(A(x))}\),注意要使得 \([x^0]A(x)\ =\ 1\)

代码

ntt版

fft版

ntt版

const int maxn = (1 << 19) | 10; // deal with n <= 1e5
const int mod = 998244353, proot = 3;

inline int pow_mod(int x, int n) { int y = 1; while(n) { if(n & 1) y = (long long) y * x % mod; x = (long long) x * x % mod; n >>= 1; } return y; }

namespace Poly {
	struct poly {
		vector<int> p;

		poly() { p.clear(); }
		poly(const vector<int> &q): p(q) {}
		poly(int n, int *a) { p.resize(n); for (int i = 0; i < n; ++i) p[i] = a[i]; }

		inline int size() const { return p.size(); }
		inline void resize(int n) { p.resize(n); return; }
	};

	void dft(int n, int *a, bool rev) {
		for (int i = 0, j = 0; i < n; ++i) {
			if(i < j) swap(a[i], a[j]);
			for (int k = n >> 1; (j ^= k) < k; k >>= 1);
		}
		for (int hl = 1, l = 2; l <= n; hl = l, l <<= 1) {
			int wn = pow_mod(proot, (mod - 1) / l); if(rev) wn = pow_mod(wn, mod - 2);
			for (int i = 0; i < n; i += l) for (int j = 0, *x = a + i, *y = x + hl, w = 1; j < hl; ++j, ++x, ++y, w = (long long) w * wn % mod) {
				int t = (long long) *y * w % mod; *y = (*x - t) % mod; *x = (*x + t) % mod;
			}
		}
		if(rev) { int inv = pow_mod(n, mod - 2); for (int i = 0; i < n; ++i) a[i] = (long long) a[i] * inv % mod; }
		return;
	}

	poly operator * (const poly &A, const poly &B) {
		static int a[maxn], b[maxn];
		int n = A.size(), m = B.size();
		if(n < 10 || m < 10 || n + m - 1 < 80) {
			int N = n + m - 1; static int c[maxn]; memset(c, 0, sizeof(c[0]) * N);
			for (int i = 0; i < n; ++i) for (int x = A.p[i], j = 0; j < m; ++j) c[i + j] = ((long long) x * B.p[j] + c[i + j]) % mod;
			return poly(N, c);
		}
		int N = 1; while(N < n + m - 1) N <<= 1;
		for (int i = 0; i < N; ++i) a[i] = i < n ? A.p[i] : 0, b[i] = i < m ? B.p[i] : 0;
		dft(N, a, 0); dft(N, b, 0); for (int i = 0; i < N; ++i) a[i] = (long long) a[i] * b[i] % mod; dft(N, a, 1);
		return poly(n + m - 1, a);
	}

	poly operator * (const int &a, const poly &B) {
		static int b[maxn]; int n = B.size();
		for (int i = 0; i < n; ++i) b[i] = (long long) a * B.p[i] % mod;
		return poly(n, b);
	}

	poly operator - (const poly &A, const poly &B) {
		static int a[maxn]; int n = A.size(), m = B.size(), N = max(n, m);
		for (int i = 0; i < N; ++i) a[i] = ((i < n ? A.p[i] : 0) - (i < m ? B.p[i] : 0)) % mod;
		return poly(N, a);
	}

	poly operator + (const poly &A, const poly &B) {
		static int a[maxn]; int n = A.size(), m = B.size(), N = max(n, m);
		for (int i = 0; i < N; ++i) a[i] = ((i < n ? A.p[i] : 0) + (i < m ? B.p[i] : 0)) % mod;
		return poly(N, a);
	}

	poly inv(int n, const poly &A) { // A(x)^-1 mod x^n
		if(n == 1) { return poly(vector<int>{pow_mod(A.p[0], mod - 2)}); }
		static poly B0, B, tA;
		B0 = inv(n + 1 >> 1, A);
		if(A.size() < n) tA = A, tA.resize(n);
		else tA.p.clear(), tA.p.insert(tA.p.end(), A.p.begin(), A.p.begin() + n);
		B = 2 * B0 - B0 * B0 * tA;
		B.resize(n);
		return B;
	}

	poly rev(const poly &A) { static poly B; B = A; reverse(B.p.begin(), B.p.end()); return B; }

	poly operator / (const poly &A, const poly &B) {
		static poly rA, rB, C, D; int n = A.size(), m = B.size();
		rA = rev(A); rB = rev(B); C = inv(n - m + 1, rB); D = C * rA; D.resize(n - m + 1);
		return rev(D);
	}

	poly operator % (const poly &A, const poly &B) { static poly D, ret; D = A / B; ret = A - B * D; ret.resize(B.size() - 1); return ret; }

	poly Rem(const poly &A, const poly &B, const poly &D) { static poly ret; ret = A - B * D; ret.resize(B.size() - 1); return ret; }

	poly dao(const poly &A) {
		static int a[maxn]; int n = A.size() - 1;
		for (int i = 0; i < n; ++i) a[i] = (long long) A.p[i + 1] * (i + 1) % mod;
		return poly(n, a);
	}

	poly ji(const poly &A) {
		static int a[maxn]; int n = A.size();
		for (int i = 1; i <= n; ++i) a[i] = (long long) A.p[i - 1] * pow_mod(i, mod - 2) % mod; a[0] = 0;
		return poly(n + 1, a);
	}

	poly ln(int n, const poly &A) { // ln(A(x)) mod x^n
		static poly B, C;
		C = dao(A); B = inv(n, A);
		B = B * C; B = ji(B); B.resize(n);
		return B;
	}

	poly exp(int n, const poly &A) { // e^A(x) mod x^n
		if(n == 1) { return poly(vector<int>{1}); }
		static poly B0, B, C, tA;
		B0 = exp(n + 1 >> 1, A);
		if(A.size() < n) tA = A, tA.resize(n);
		else tA.p.clear(), tA.p.insert(tA.p.end(), A.p.begin(), A.p.begin() + n);
		C = ln(n, B0);
		B = B0 * (poly(vector<int>{1}) - C + tA);
		B.resize(n);
		return B;
	}

	poly pow_mod(int n, const poly &A, int N) { // A(x)^N mod x^n
		static poly B;
		int na = A.size(), i = 0;
		while(i < na && A.p[i] == 0) ++i;
		if(i == na || (long long) N * i >= n) return poly(vector<int>(n, 0));
		if(i) {
			B.p.clear();
			for (int j = i; j < na; ++j) B.p.push_back(A.p[j]);
			static poly C;
			C = pow_mod(n - N * i, B, N);
			B.p = vector<int>(N * i, 0);
			B.p.insert(B.p.end(), C.p.begin(), C.p.end());
			return B;
		}
		if(A.p[0] != 1) {
			int t = ::pow_mod(A.p[0], mod - 2), s = ::pow_mod(A.p[0], N); B.resize(na);
			for (int i = 0; i < na; ++i) B.p[i] = (long long) A.p[i] * t % mod;
			B = pow_mod(n, B, N);
			for (int i = 0; i < n; ++i) B.p[i] = (long long) s * B.p[i] % mod;
			return B;
		}
		B = ln(n, A);
		for (int i = 0; i < n; ++i) B.p[i] = (long long) B.p[i] * N % mod;
		return exp(n, B);
	}
}

fft版

const int maxn = (1 << 19) | 10, lgn = 20; // deal with n <= 1e5
const int mod = 1e9 + 7;

inline int pow_mod(int x, int n) { int y = 1; while(n) { if(n & 1) y = (long long) y * x % mod; x = (long long) x * x % mod; n >>= 1; } return y; }

namespace Poly { // don't forget to init
	struct comp {
		double r, i;

		comp() {}
		comp(double r, double i): r(r), i(i) {}
	} w[lgn][maxn];

	comp operator * (const comp &a, const comp &b) { return comp(a.r * b.r - a.i * b.i, a.r * b.i + a.i * b.r); }
	comp operator + (const comp &a, const comp &b) { return comp(a.r + b.r, a.i + b.i); }
	comp operator - (const comp &a, const comp &b) { return comp(a.r - b.r, a.i - b.i); }
	comp conj(const comp &a) { return comp(a.r, -a.i); }

	void init() {
		static const double pi = acos(-1);
		for (int i = 0; i < lgn; ++i) {
			int n = 1 << i;
			if(n >= maxn) break;
			for (int j = 0; j < n; ++j) {
				if(i && (~j & 1)) w[i][j] = w[i - 1][j >> 1];
				else {
					double ang = 2. * pi * j / (double) n;
					w[i][j] = comp(cos(ang), sin(ang));
				}
			}
		}
		return;
	}

	inline int toint(double x, double invn) { return llround(x * invn) % mod; }

	void dft(int n, comp *a) {
		for (int i = 0, j = 0; i < n; ++i) {
			if(i < j) swap(a[i], a[j]);
			for (int k = n >> 1; (j ^= k) < k; k >>= 1);
		}
		for (int hl = 1, l = 2, wi = 1; l <= n; hl = l, l <<= 1, ++wi) for (int i = 0; i < n; i += l) {
			comp *x = a + i, *y = x + hl, *z = w[wi];
			for (int j = 0; j < hl; ++j, ++x, ++y, ++z) {
				comp t = *y * *z; *y = *x - t; *x = *x + t;
			}
		}
		return;
	}

	struct poly {
		vector<int> p;

		poly() { p.clear(); }
		poly(const vector<int> &q): p(q) {}
		poly(int n, int *a) { p.resize(n); for (int i = 0; i < n; ++i) p[i] = a[i]; }

		inline int size() const { return p.size(); }
		inline void resize(int n) { p.resize(n); return; }
	};

	poly operator * (const poly &A, const poly &B) {
		int n = A.size(), m = B.size();
		if(n < 10 || m < 10 || n + m - 1 < 80) {
			int N = n + m - 1; static int c[maxn]; memset(c, 0, sizeof(c[0]) * N);
			for (int i = 0; i < n; ++i) for (int x = A.p[i], j = 0; j < m; ++j) c[i + j] = ((long long) x * B.p[j] + c[i + j]) % mod;
			return poly(N, c);
		}
		int N = 1; while(N < n + m - 1) N <<= 1;
		static comp a[maxn], b[maxn], da[maxn], db[maxn], dc[maxn], dd[maxn];
		for (int i = 0; i < n; ++i) { int x = (A.p[i] + mod) % mod; a[i] = comp(x & 32767, x >> 15); } for (int i = n; i < N; ++i) a[i] = comp(0, 0);
		for (int i = 0; i < m; ++i) { int x = (B.p[i] + mod) % mod; b[i] = comp(x & 32767, x >> 15); } for (int i = m; i < N; ++i) b[i] = comp(0, 0);
		dft(N, a); dft(N, b);
		for (int i = 0; i < N; ++i) {
			int j = N - i & N - 1; static comp DA, DB, DC, DD;
			DA = (a[i] + conj(a[j])) * comp(.5, 0); DB = (a[i] - conj(a[j])) * comp(0, -.5);
			DC = (b[i] + conj(b[j])) * comp(.5, 0); DD = (b[i] - conj(b[j])) * comp(0, -.5);
			da[j] = DA * DC; db[j] = DA * DD;
			dc[j] = DB * DC; dd[j] = DB * DD;
		}
		for (int i = 0; i < N; ++i) a[i] = da[i] + db[i] * comp(0, 1), b[i] = dc[i] + dd[i] * comp(0, 1);
		dft(N, a); dft(N, b); double iN = 1.0 / N; N = n + m - 1; static int c[maxn];
		for (int i = 0; i < N; ++i) {
			static int da, db, dc, dd;
			da = toint(a[i].r, iN), db = toint(a[i].i, iN), dc = toint(b[i].r, iN), dd = toint(b[i].i, iN);
			c[i] = (((long long) dd << 30) + ((long long) db + dc << 15) + da) % mod;
		}
		return poly(N, c);
	}

	poly operator * (const int &a, const poly &B) {
		static int b[maxn]; int n = B.size();
		for (int i = 0; i < n; ++i) b[i] = (long long) a * B.p[i] % mod;
		return poly(n, b);
	}

	poly operator - (const poly &A, const poly &B) {
		static int a[maxn]; int n = A.size(), m = B.size(), N = max(n, m);
		for (int i = 0; i < N; ++i) a[i] = ((i < n ? A.p[i] : 0) - (i < m ? B.p[i] : 0)) % mod;
		return poly(N, a);
	}

	poly operator + (const poly &A, const poly &B) {
		static int a[maxn]; int n = A.size(), m = B.size(), N = max(n, m);
		for (int i = 0; i < N; ++i) a[i] = ((i < n ? A.p[i] : 0) + (i < m ? B.p[i] : 0)) % mod;
		return poly(N, a);
	}

	poly inv(int n, const poly &A) { // A(x)^-1 mod x^n
		if(n == 1) { return poly(vector<int>{pow_mod(A.p[0], mod - 2)}); }
		static poly B0, B, tA;
		B0 = inv(n + 1 >> 1, A);
		if(A.size() < n) tA = A, tA.resize(n);
		else tA.p.clear(), tA.p.insert(tA.p.end(), A.p.begin(), A.p.begin() + n);
		B = 2 * B0 - B0 * B0 * tA;
		B.resize(n);
		return B;
	}

	poly rev(const poly &A) { static poly B; B = A; reverse(B.p.begin(), B.p.end()); return B; }

	poly operator / (const poly &A, const poly &B) {
		static poly rA, rB, C, D; int n = A.size(), m = B.size();
		rA = rev(A); rB = rev(B); C = inv(n - m + 1, rB); D = C * rA; D.resize(n - m + 1);
		return rev(D);
	}

	poly operator % (const poly &A, const poly &B) { static poly D, ret; D = A / B; ret = A - B * D; ret.resize(B.size() - 1); return ret; }

	poly Rem(const poly &A, const poly &B, const poly &D) { static poly ret; ret = A - B * D; ret.resize(B.size() - 1); return ret; }

	poly dao(const poly &A) {
		static int a[maxn]; int n = A.size() - 1;
		for (int i = 0; i < n; ++i) a[i] = (long long) A.p[i + 1] * (i + 1) % mod;
		return poly(n, a);
	}

	poly ji(const poly &A) {
		static int a[maxn]; int n = A.size();
		for (int i = 1; i <= n; ++i) a[i] = (long long) A.p[i - 1] * pow_mod(i, mod - 2) % mod; a[0] = 0;
		return poly(n + 1, a);
	}

	poly ln(int n, const poly &A) { // ln(A(x)) mod x^n
		static poly B, C;
		C = dao(A); B = inv(n, A);
		B = B * C; B = ji(B); B.resize(n);
		return B;
	}

	poly exp(int n, const poly &A) { // e^A(x) mod x^n
		if(n == 1) { return poly(vector<int>{1}); }
		static poly B0, B, C, tA;
		B0 = exp(n + 1 >> 1, A);
		if(A.size() < n) tA = A, tA.resize(n);
		else tA.p.clear(), tA.p.insert(tA.p.end(), A.p.begin(), A.p.begin() + n);
		C = ln(n, B0);
		B = B0 * (poly(vector<int>{1}) - C + tA);
		B.resize(n);
		return B;
	}

	poly pow_mod(int n, const poly &A, int N) { // A(x)^N mod x^n
		static poly B;
		int na = A.size(), i = 0;
		while(i < na && A.p[i] == 0) ++i;
		if(i == na || (long long) N * i >= n) return poly(vector<int>(n, 0));
		if(i) {
			B.p.clear();
			for (int j = i; j < na; ++j) B.p.push_back(A.p[j]);
			static poly C;
			C = pow_mod(n - N * i, B, N);
			B.p = vector<int>(N * i, 0);
			B.p.insert(B.p.end(), C.p.begin(), C.p.end());
			return B;
		}
		if(A.p[0] != 1) {
			int t = ::pow_mod(A.p[0], mod - 2), s = ::pow_mod(A.p[0], N); B.resize(na);
			for (int i = 0; i < na; ++i) B.p[i] = (long long) A.p[i] * t % mod;
			B = pow_mod(n, B, N);
			for (int i = 0; i < n; ++i) B.p[i] = (long long) s * B.p[i] % mod;
			return B;
		}
		B = ln(n, A);
		for (int i = 0; i < n; ++i) B.p[i] = (long long) B.p[i] * N % mod;
		return exp(n, B);
	}
}
posted @ 2018-07-09 10:50  King_George  阅读(257)  评论(0编辑  收藏  举报