码迷,mamicode.com
首页 > 其他好文 > 详细

【模板】多项式全家桶_缺斤少两

时间:2020-05-07 01:11:06      阅读:79      评论:0      收藏:0      [点我收藏+]

标签:while   dig   牛顿迭代   emc   fan   memcpy   exp   typename   str   

本文自用,不保证代码可读

功能目录

ntt, intt, 乘法, 积分, 求导(微分), ln, exp, 开方(牛顿迭代)

代码

#include<bits/stdc++.h>
#define LL long long
#define MAXN 100010
#define MOD 998244353
#define YG 3
using namespace std;
template<typename T>void Read(T &cn)
{
	char c;int sig = 1;
	while(!isdigit(c = getchar()))if(c == ‘-‘)sig = -1; cn = c-48;
	while(isdigit(c = getchar()))cn = cn*10+c-48; cn*=sig;
}
template<typename T>void Write(T cn)
{
	if(cn < 0) {putchar(‘-‘); cn = 0-cn; }
	int wei = 0; T cm = 0; int cx = cn%10; cn/=10;
	while(cn)cm = cm*10+cn%10,cn/=10,wei++;
	while(wei--)putchar(cm%10+48),cm/=10;
	putchar(cx+48);
}
const int MAXNTT = MAXN*4+1;
LL jie[MAXN+1], nij[MAXN+1];
LL omg[MAXNTT], inv[MAXNTT], fan[MAXNTT], Mn;
map<int, int> M;
int erwei(int cn) {int guo = 0; while(cn) guo++, cn>>=1; return guo; }
LL ksm(LL cn, LL cm) {LL ans = 1; while(cm) ans = ans*(1+(cn-1)*(cm&1))%MOD, cn = cn*cn%MOD, cm>>=1; return ans; }
int BSGS(int cn)
{
	int Blo = sqrt(MOD)+1;
	int xian = 1, yuan = 3; M.clear();
	for(int i = 1;i<=Blo;i++) M[(int)(1ll*cn*xian%MOD)] = i-1, xian = 1ll*xian*yuan%MOD;
	yuan = xian; xian = 1;
	for(int i = 1;i<=Blo;i++)
	{
		xian = 1ll*xian * yuan%MOD;
		if(M.find(xian) != M.end()) return Blo * i - M[xian];
	}
}
void yuchu_omg(int cn)
{
	Mn = cn = 1<<erwei(cn); omg[0] = inv[0] = 1; omg[1] = ksm(YG, MOD/cn); inv[1] = ksm(omg[1], MOD-2);
	for(int i = 2;i<cn;i++) omg[i] = omg[i-1]*omg[1]%MOD, inv[i] = inv[i-1]*inv[1]%MOD;
}
void yuchu_jie(int cn)
{
	jie[0] = 1;
	for(int i = 1;i<=cn;i++) jie[i] = jie[i-1]*i%MOD;
	nij[cn] = ksm(jie[cn], MOD-2);
	for(int i = cn-1;i>=0;i--) nij[i] = nij[i+1]*(i+1)%MOD;
}
struct Poly{
	int n;
	LL a[MAXNTT];
	void yuchu_fan(int cn) {int lin = erwei(cn)-2; fan[0] = 0; for(int i = 1;i<cn;i++) fan[i] = (fan[i>>1]>>1) | ((i&1)<<lin); }
	void qing(int cn) {memset(a+n, 0, sizeof(a[0])*(cn-n)); }
	void kong() {memset(a,0,sizeof(a)); }//for safety
	void do_ntt(LL omg[], int cn)
	{
		yuchu_fan(cn);
		for(int i = 0;i<cn;i++) if(fan[i] > i) swap(a[i], a[fan[i]]);
		for(int i = 2, m = 1;i<=cn;i = (m = i)<<1)
		for(int j = 0;j<cn;j+=i)
		for(int k = 0;k<m;k++)
		{
			LL lin1 = a[j+k], lin2 = 1ll*a[j+k+m] * omg[Mn/i*k]%MOD;
			a[j+k] = lin1 + lin2 < MOD ? lin1+lin2 : lin1+lin2-MOD;
			a[j+k+m] = lin1-lin2 < 0 ? lin1-lin2+MOD : lin1-lin2;
		}
	}
	void enlarge(int cn) {qing(cn); n = cn; }
	void set_len(int cn) {n = cn; }
	void copy(Poly &A, int cn) {int lin = min(cn, A.n); memcpy(a, A.a, lin*sizeof(a[0])); if(lin < cn-1) memset(a+lin, 0, (cn-lin)*sizeof(a[0])); n = cn; }
	void ntt(int cn) {qing(cn); do_ntt(omg, cn); }
	void intt(int cn, int cm) {do_ntt(inv, cn); int lin = ksm(cn, MOD-2); for(int i = 0;i<cm;i++) a[i] = 1ll*a[i]*lin%MOD; set_len(cm); qing(cn); }
	void wf() {for(int i = 0;i<n-1;i++) a[i] = 1ll*a[i+1]*(i+1)%MOD; a[n-1] = 0; n--;  }
	void jf() {for(int i = n;i>=1;i--) a[i] = 1ll*a[i-1]*jie[i-1]%MOD*nij[i]%MOD; a[0] = 0; n++; }
	void getit(int cn) {n = cn; for(int i = 0;i<n;i++) Read(a[i]); }
	void outit() {for(int i = 0;i<n;i++) Write(a[i]), putchar(‘ ‘); puts(""); }
	void outit(int cn) {for(int i = 0;i<cn;i++) Write(a[i]), putchar(‘ ‘); puts(""); }
	inline int r(int cn) {return cn < n ? a[cn] : 0; }//get an element safely
};
Poly A, B, G, nG, nB, nA, cA, cB;
void cheng(Poly &A, Poly &B, Poly &C)
{
	int lin = 1<<erwei(A.n + B.n-1);
	cA.copy(A, A.n); cB.copy(B, B.n);
	cA.ntt(lin); cB.ntt(lin); 
	for(int i = 0;i<lin;i++) C.a[i] = 1ll*cA.a[i]*cB.a[i]%MOD;
	C.intt(lin, A.n+B.n-1);
}
void qiu_ni(Poly &A, Poly &B, Poly &G, int n)
{
	if(n == 1) {B.a[0] = ksm(A.a[0], MOD-2); B.set_len(1); return; }
	int m = (n+1)>>1, lin = 1<<erwei(m*2+n); qiu_ni(A, G, B, m); G.enlarge(n);
	nG.copy(G, m); B.copy(A, n); 
	B.ntt(lin); nG.ntt(lin);
	for(int i = 0;i<lin;i++) B.a[i] = 1ll*B.a[i]*nG.a[i]%MOD*nG.a[i]%MOD;
	B.intt(lin, n); 
	for(int i = 0;i<n;i++) B.a[i] = (2*G.a[i] - B.a[i]+MOD)%MOD;
}
void qiu_ln(Poly &A, Poly &B, Poly &G, int n)
{
	nA.copy(A, n);
	qiu_ni(nA, G, B, n); nA.wf(); 
	cheng(nA, G, B); B.set_len(n-1); B.jf(); 
}
void qiu_exp(Poly &A, Poly &B, Poly &G, int n)
{
	if(n == 1) {B.a[0] = 1; B.set_len(1); return; }
	int m = (n+1)>>1; qiu_exp(A, G, B, m); G.enlarge(n);
	qiu_ln(G, nB, B, n);
	for(int i = 0;i<n;i++) nB.a[i] = (A.a[i] - nB.a[i] + MOD)%MOD; nB.a[0]++;
	cheng(nB, G, B); B.set_len(n);
}
void qiu_sqrt(Poly &A, Poly &B, Poly &G, int n)
{
	if(n == 1) {int lin = BSGS(A.a[0]); B.a[0] = min(ksm(3, lin/2), ksm(3, (lin+MOD-1)/2)); B.set_len(1); return; }
	int m = (n+1)>>1; qiu_sqrt(A, G, B, m); G.enlarge(n);
	qiu_ni(G, nB, B, n); nA.copy(A, n); cheng(nB, nA, B); 
	for(int i = 0;i<n;i++) B.a[i] = 1ll*nij[2]*(G.a[i] + B.a[i])%MOD;
	B.set_len(n);
}
int n;
int main()
{
	Read(n); A.getit(n); yuchu_omg(n*2+1); yuchu_jie(n+1);
	qiu_sqrt(A,B,G,n);
	B.outit();
	return 0;
}

【模板】多项式全家桶_缺斤少两

标签:while   dig   牛顿迭代   emc   fan   memcpy   exp   typename   str   

原文地址:https://www.cnblogs.com/czyarl/p/12840002.html

(0)
(0)
   
举报
评论 一句话评论(0
登录后才能评论!
© 2014 mamicode.com 版权所有  联系我们:gaon5@hotmail.com
迷上了代码!