多项式开方

2021-01-27 15:35:23 阅读：155 来源： 互联网

\[\begin{align} A(x) &\equiv B'(x)^2 &\mod x^{n/2}\\ A(x) &\equiv B(x)^2 &\mod x^{n}\\ A(x) &\equiv B(x)^2 &\mod x^{n/2}\\ B(x)^2 - B'(x)^2 &\equiv 0 &\mod x^{n/2}\\ B^4(x)+B'(x)^4-2B(x)^2B'(x)^2 &\equiv 0 &\mod x^n\\ B^4(x)+B'(x)^4 + 2B(x)^2B'(x)^2 &\equiv 4B(x)^2B'(x)^2 &\mod x^n\\ B^2(x)+B'(x)^2 &\equiv 2B(x)B'(x) &\mod x^n\\ A(x)+B'(x)^2 &\equiv 2B(x)B'(x) &\mod x^n\\ B(x) &\equiv \dfrac{A(x)+B'(x)^2}{2B'(x)} &\mod x^n \end{align} \]

跟求逆用的那种倍增方法极度相似。

就是要用 O(n log² n) 的时间来求了。真的要学牛顿迭代了吗（

但是这个好像可以在倍增的过程中维护 B 的逆，有点可以，可以单 log，以后再补。

#include<bits/stdc++.h>
typedef long long LL;
using namespace std;

const int N = 4e5 + 233, mo = 998244353;
LL ksm(LL a, LL b) {
	LL res = 1ll;
	for(; b; b >>= 1, a = a * a % mo)
		if(b & 1) res = res * a % mo;
	return res;
}
const LL g = 3, ig = ksm(g, mo - 2), inv_2 = ksm(2, mo - 2);

int rv[N];
void NTT(LL *a, int n, int type) {
	for(int i = 0; i < n; ++i) if(i < rv[i]) swap(a[i], a[rv[i]]);
	for(int m = 2; m <= n; m = m << 1) {
		LL w = ksm(type == 1 ? g : ig, (mo - 1) / m);
		for(int i = 0; i < n; i += m) {
			LL tmp = 1ll;
			for(int j = 0; j < (m >> 1); ++j) {
				LL p = a[i + j], q = tmp * a[i + j + (m >> 1)] % mo;
				a[i + j] = (p + q) % mo, a[i + j + (m >> 1)] = (p - q + mo) % mo;
				tmp = tmp * w % mo;
			}
		}
	}
	if(type == -1) {
		LL Inv = ksm(n, mo - 2);
		for(int i = 0; i < n; ++i) a[i] = a[i] * Inv % mo;
	}
}

LL t[N];
void poly_inv(int deg, LL *a, LL *b) {
	if(deg == 1) { b[0] = ksm(a[0], mo - 2); return; }
	poly_inv((deg + 1) >> 1, a, b);
	int len = 1; while(len < (deg << 1)) len = len << 1;
	for(int i = 0; i < deg; ++i) t[i] = a[i];
	for(int i = deg; i < len; ++i) b[i] = t[i] = 0ll;
	for(int i = 1; i < len; ++i) rv[i] = (rv[i>>1]>>1) | (i&1?len>>1:0);
	NTT(b, len, 1), NTT(t, len, 1);
	for(int i = 0; i < len; ++i) b[i] = b[i] * (2ll - t[i] * b[i] % mo) % mo;
	NTT(b, len, -1);
	for(int i = deg; i < len; ++i) b[i] = 0ll;
}

LL binv[N], c[N];
void poly_sqrt(int deg, LL *a, LL *b) {
	if(deg == 1) { b[0] = 1ll; return; }
	poly_sqrt((deg + 1) >> 1, a, b);
	for(int i = 0; i < deg; ++i) c[i] = 2ll * b[i] % mo;
	poly_inv(deg, c, binv);
	int len = 1; while(len < (deg << 1)) len = len << 1;
	for(int i = 1; i < len; ++i) rv[i] = (rv[i>>1]>>1) | (i&1?len>>1:0);
	for(int i = 0; i < deg; ++i) c[i] = a[i];
	for(int i = deg; i < len; ++i) b[i] = binv[i] = c[i] = 0ll;
	NTT(b, len, 1), NTT(binv, len, 1), NTT(c, len, 1);
	for(int i = 0; i < len; ++i) b[i] = (c[i] + b[i] * b[i] % mo ) % mo * binv[i] % mo;
	NTT(b, len, -1);
	for(int i = deg; i < len; ++i) b[i] = 0ll;
}

int n;
LL a[N], b[N];
int main() {
	scanf("%d", &n);
	for(int i = 0; i < n; ++i) scanf("%lld", &a[i]);
	poly_sqrt(n, a, b);
	for(int i = 0; i < n; ++i) cout << (b[i]%mo + mo) % mo << ' ';
	return 0;
}

标签：int,多项式,mo,len,++,开方,LL,deg
来源： https://www.cnblogs.com/tztqwq/p/14335243.html

本站声明： 1. iCode9 技术分享网（下文简称本站）提供的所有内容，仅供技术学习、探讨和分享；
2. 关于本站的所有留言、评论、转载及引用，纯属内容发起人的个人观点，与本站观点和立场无关；
3. 关于本站的所有言论和文字，纯属内容发起人的个人观点，与本站观点和立场无关；
4. 本站文章均是网友提供，不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属；如您发现该文章侵犯了您的权益，可联系我们第一时间进行删除；
5. 本站为非盈利性的个人网站，所有内容不会用来进行牟利，也不会利用任何形式的广告来间接获益，纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。

ICode9

多项式开方