首页 > 其他分享 >我的多项式全家桶

我的多项式全家桶

时间:2024-01-22 09:02:09浏览次数:29  
标签:return lc int 多项式 全家 mid poly rc

第一个自己实现的多项式全家桶。打比赛时终于有板子了!

代码是从之前学转置原理的那篇博客里升级来的。

但是功能很残?效率很逊?接口很怪?评价是能用就行。

目前封装了:乘、逆、除(取模)、开方(无二次剩余,不会 qwq)、对数、指数、多点求值、快速插值、常系数齐次线性递推、Berlekamp-Massey 算法。

为什么突然开始多项式?原因是在尝试改鸡的时候发现了 sjy 的究极神秘题解,感觉 BM 算法竟然有意想不到的用处,于是学了 BM 算法,然后递归到了常系数其次线性递推,然后递归到了除法,然后升级了多项式板子。

不过话说 BM 算法其实与多项式没有很大关系不是吗?

#include <cstdio>
#include <vector>
#include <cassert>
#include <algorithm>
#define lc (p<<1)
#define rc (p<<1|1)
#define ALL(p) p.begin(),p.end()
#define IL inline
// #define getchar() (p1==p2&&(p2=(p1=buf)+fread(buf,1,1<<22,stdin)),p1==p2?EOF:*p1++)
using namespace std;
char buf[1<<22],*p1=buf,*p2=buf;
int read(){
	char c=getchar();int x=0;bool f=0;
	while(c<48||c>57) f|=(c=='-'),c=getchar();
	do x=(x<<1)+(x<<3)+(c^48),c=getchar();
	while(c>=48&&c<=57);
	if(f) return -x;
	return x;
}
typedef vector<int> vi;
const int P=998244353;
typedef long long ll;
IL int qp(int a,int b=P-2){
	int res=1;
	while(b){
		if(b&1) res=(ll)res*a%P;
		a=(ll)a*a%P;b>>=1;
	}
	return res;
}
int len,ilen,bt;
int rev[1<<20],cw[1<<20|1];
IL void init(int _len){ // mod x^len
	len=1;bt=-1;
	while(len<_len) len<<=1,++bt;
	int w=qp(3,(P-1)>>(bt+1));
	cw[0]=cw[len]=1;
	for(int i=1;i<len;++i){
		cw[i]=(ll)cw[i-1]*w%P;
		rev[i]=(rev[i>>1]>>1)|((i&1)<<bt);
	}
	ilen=qp(len);
}
struct poly{
	vi f;
	poly():f(){}
	poly(int Len):f(Len){}
	poly(initializer_list<int> Init):f(Init){}
	IL int& operator[](const int &x){return f[x];}
	IL const int& operator[](const int &x)const{return f[x];}
	IL void NTT(){
		f.resize(len,0);
		for(int i=1;i<len;++i) if(rev[i]<i) swap(f[rev[i]],f[i]);
		for(int i=1,tt=len>>1;i<len;i<<=1,tt>>=1)
			for(int j=0;j<len;j+=(i<<1))
				for(int k=j,t=0;k<(j|i);++k,t+=tt){
					int x=f[k],y=(ll)f[k|i]*cw[t]%P;
					if((f[k]=x+y)>=P) f[k]-=P;
					if((f[k|i]=x-y)<0) f[k|i]+=P;
				}
	}
	IL void INTT(){
		for(int i=1;i<len;++i) if(rev[i]<i) swap(f[rev[i]],f[i]);
		for(int i=1,tt=len>>1;i<len;i<<=1,tt>>=1)
			for(int j=0;j<len;j+=(i<<1))
				for(int k=j,t=len;k<(j|i);++k,t-=tt){
					int x=f[k],y=(ll)f[k|i]*cw[t]%P;
					if((f[k]=x+y)>=P) f[k]-=P;
					if((f[k|i]=x-y)<0) f[k|i]+=P;
				}
		for(int i=0;i<len;++i) f[i]=(ll)f[i]*ilen%P;
	}
	IL int size(){return f.size();}
	IL void reduce(){while(!f.empty()&&!f.back()) f.pop_back();}
	IL void rever(){reverse(ALL(f));}
	IL void show(){
		for(int x:f) printf("%d ",x);
		putchar('\n');
	}
	IL void trunc(int lim){
		if(lim<int(f.size())) f.erase(f.begin()+lim,f.end());
	}
	IL poly inv(int lim){ // mod x^lim
		assert(f[0]);
		poly cur({qp(f[0])});
		for(int t=1;t<lim;t<<=1){
			poly ff(t<<2);
			copy(f.begin(),f.begin()+min(t<<1,size()),ff.f.begin());
			init(t<<2);ff.NTT();cur.NTT();
			poly tmp(len);
			for(int i=0;i<len;++i){
				tmp[i]=(2ll-(ll)cur[i]*ff[i])%P*cur[i]%P;
				if(tmp[i]<0) tmp[i]+=P;
			}
			tmp.INTT();
			cur.f.swap(tmp.f);
			cur.trunc(t<<1);
		}
		cur.trunc(lim);
		return cur;
	}
	friend poly operator+(poly A,poly B){
		int mx=max(A.size(),B.size());
		A.f.resize(mx,0);B.f.resize(mx,0);
		poly C(mx);
		for(int i=0;i<mx;++i){
			C[i]=A[i]+B[i];
			if(C[i]>=P) C[i]-=P;
		}
		return C;
	}
	friend poly operator-(poly A,poly B){
		int mx=max(A.size(),B.size());
		A.f.resize(mx,0);B.f.resize(mx,0);
		poly C(mx);
		for(int i=0;i<mx;++i){
			C[i]=A[i]-B[i];
			if(C[i]<0) C[i]+=P;
		}
		C.reduce();
		return C;
	}
	friend poly operator*(poly A,poly B){
		init(A.size()+B.size()-1);A.NTT();B.NTT();
		poly C(len);
		for(int i=0;i<len;++i) C[i]=(ll)A[i]*B[i]%P;
		C.INTT();C.reduce();
		return C;
	}
	poly operator*(int X){
		int n=size();
		poly F(n);
		for(int i=0;i<n;++i) F[i]=(ll)f[i]*X%P;
		return F;
	}
	friend poly operator%(poly F,poly G){
		int n=F.size()-1,m=G.size()-1;
		if(n<m) return F;
		F.rever();G.rever();
		poly Q=G.inv(n-m+1);
		Q.prod(Q,F);Q.f.resize(n-m+1,0);
		F.rever();G.rever();Q.rever();
		poly R=F-Q*G;
		R.reduce();
		return R;
	}
	IL void prod(poly A,poly B){
		init(A.size()+B.size()-1);A.NTT();B.NTT();
		f.resize(len);
		for(int i=0;i<len;++i) f[i]=(ll)A[i]*B[i]%P;
		INTT();reduce();
	}
	IL void prodT(poly A,poly B){
		int an=A.size()-1,bn=B.size()-1;
		reverse(ALL(B.f));prod(A,B);
		for(int i=0;i<=an;++i) f[i]=f[i+bn];
		trunc(an+1);
	}
	IL int calc(int t){
		int pw=1,res=0;
		for(int x:f){
			res=(res+(ll)pw*x)%P;
			pw=(ll)pw*t%P;
		}
		return res;
	}
	IL poly deriv(){
		int n=f.size();
		poly D(n-1);
		for(int i=1;i<n;++i) D[i-1]=(ll)f[i]*i%P;
		return D;
	}
	IL poly integ(){
		int n=f.size();
		poly D(n+1),inv(n+1);
		for(int i=1;i<=n;++i){
			if(i>1) inv[i]=(ll)inv[P%i]*(P-P/i)%P;
			else inv[1]=1;
			D[i]=(ll)f[i-1]*inv[i]%P;
		}
		return D;
	}
	IL poly getln(int lim){ // mod x^lim
		assert(f[0]==1);
		poly F=deriv()*inv(lim);
		F.f.resize(lim-1,0);
		return F.integ();
	}
	IL poly getexp(int lim){ // mod x^lim
		assert(f[0]==0);
		poly F({1});
		for(int t=1;t<lim;t<<=1){
			poly ff(t<<1);
			copy(f.begin(),f.begin()+min(t<<1,size()),ff.f.begin());
			F=F*(poly({1})-F.getln(t<<1)+ff);
			F.f.resize(t<<1,0);
		}
		F.trunc(lim);
		return F;
	}
	IL poly sqrt(int lim){ // without quadratic residue , mod x^lim
		assert(f[0]==1);
		poly F({1});
		for(int t=1;t<lim;t<<=1){
			poly ff(t<<1);
			copy(f.begin(),f.begin()+min(t<<1,size()),ff.f.begin());
			F=(F+F.inv(t<<1)*ff)*((P+1)>>1);
			F.f.resize(t<<1,0);
		}
		F.trunc(lim);
		return F;
	}
}F;
namespace extend{
	const int N=200003;
	int n;
	int px[N],py[N],res[N];
	namespace multieval{
		poly G[N<<2],H[N<<2];
		void calc(int p,int l,int r){
			if(l==r){G[p]=poly({1,(P-px[r])%P});return;}
			int mid=(l+r)>>1;
			calc(lc,l,mid);
			calc(rc,mid+1,r);
			G[p].prod(G[lc],G[rc]);
		}
		void solve(int p,int l,int r){
			H[p].trunc(r-l+1);
			if(l==r){res[r]=H[p][0];return;}
			int mid=(l+r)>>1;
			G[p].f.clear();G[p].f.shrink_to_fit();
			H[lc].prodT(H[p],G[rc]);
			H[rc].prodT(H[p],G[lc]);
			H[p].f.clear();H[p].f.shrink_to_fit();
			solve(lc,l,mid);
			solve(rc,mid+1,r);
		}
		void sol(){
			// init n, px, py first
			calc(1,1,n);
			H[1].prodT(F,G[1].inv(n+1));
			solve(1,1,n);
		}
	}
	namespace interplot{
		poly T[N<<2];
		void eval(int p,int l,int r){
			if(l==r){T[p]=poly({(P-px[r])%P,1});return;}
			int mid=(l+r)>>1;
			eval(lc,l,mid);
			eval(rc,mid+1,r);
			T[p].prod(T[lc],T[rc]);
		}
		poly proc(int p,int l,int r){
			if(l==r) return poly({py[r]});
			int mid=(l+r)>>1;
			poly L=proc(lc,l,mid);
			poly R=proc(rc,mid+1,r);
			L.prod(L,T[rc]);
			R.prod(R,T[lc]);
			return L+R;
		}
		poly solve(){
			// init n, px, py first
			eval(1,1,n);
			F=T[1].deriv();
			multieval::sol();
			for(int i=1;i<=n;++i) py[i]=(ll)py[i]*qp(res[i])%P;
			return proc(1,1,n);
		}
	}
}
namespace berlekamp_massey{
	const int N=200003;
	int n,k,kk,las,pos,ex;
	int s[N],ss[N],a[N],tmp[N];
	int solve(){
		n=read();ex=read();
		for(int i=0;i<n;++i) a[i]=read();
		bool fl=1;
		for(int i=0;i<n;++i){
			int cur=a[i];
			for(int j=1;j<=k;++j) cur=(cur-(ll)a[i-j]*s[j])%P;
			if(cur<0) cur+=P;
			if(!cur) continue;
			if(fl){
				k=i+1;kk=0;pos=i;las=cur;fl=0;
				continue;
			}
			int nk=max(k,i-pos+kk);
			int val=(ll)cur*qp(las)%P;
			for(int j=1;j<i-pos;++j) tmp[j]=0;
			tmp[i-pos]=val;
			for(int j=1;j<=kk;++j) tmp[i-pos+j]=(ll)(P-val)*ss[j]%P;
			if(k-i<kk-pos){
				kk=k;pos=i;las=cur;
				for(int j=1;j<=k;++j) ss[j]=s[j];
			}
			k=nk;
			for(int j=1;j<=k;++j){
				s[j]+=tmp[j];
				if(s[j]>=P) s[j]-=P;
			}
		}
		poly M(k+1),F(2),G({1});
		M[k]=1;
		for(int i=1;i<=k;++i)
			if(s[i]) M[k-i]=P-s[i];
		F[1]=1;
		while(ex){
			if(ex&1) G=G*F%M;
			F=F*F%M;ex>>=1;
		}
		G.f.resize(k,-1);
		int ans=0;
		for(int i=0;i<k;++i) ans=(ans+(ll)G[i]*a[i])%P;
		if(ans<0) ans+=P;
		return ans;
	}
}
int main(){
	int n=read();
	poly F(n);
	for(int i=0;i<n;++i) F[i]=read();
	F=F.sqrt(n);
	F.show();
	return 0;
}

标签:return,lc,int,多项式,全家,mid,poly,rc
From: https://www.cnblogs.com/yyyyxh/p/17979236/poly

相关文章

  • 无涯教程-MATLAB - 多项式(Polynomials)
    MATLAB将多项式表示为行向量,其中包含按降序排序的系数。例如,方程P(x)=x4+7x3-5x+9可以表示为-p=[170-59];判断多项式polyval函数用于以指定值判断多项式。例如,要判断我们先前的多项式p,在x=4处,键入-p=[170-59];polyval(p,4)MATLAB执行上述语句并返......
  • 洛谷题单指南-模拟和高精度-P1067 [NOIP2009 普及组] 多项式输出
    原题链接:https://www.luogu.com.cn/problem/P1067题意解读:模拟法依次输出多项式内容即可,但是需要考虑的周全,主要有以下关键点:1、系数为0时不输出多项式2、第一个符号,只有负号才输出3、次数不为0时,不输出为1的系数;次数为0时,输出所有系数4、次数为1时,不输出次数;次数为0时不输......
  • 一些多项式常用操作的公式推导
    怕我以后忘记了看不懂自己的板子(((牛顿迭代用于求解函数零点的近似值。设函数\(f(x)\)的零点近似值为\(x_0\),过点\((x_0,f(x_0))\)作\(f(x)\)的切线,切线与\(x\)轴交点的横坐标即为新的近似值。切线解析式为\(y=f'(x_0)(x-x_0)+f(x_0)\),当\(y=0\)时\(x=x_0-\dfrac{f......
  • 多项式相关
    板子:1constintN=4e5+5;23structNTT{4constintmod=998244353,G=3,invG=332748118;5intrev[N],res[N],Inv[N],Res[N],eres[N],Eres[N];67voidpreinv(intn){8Inv[0]=Inv[1]=1;9......
  • 多项式计数杂谈
    多项式计数杂谈-command_block的博客-洛谷博客command_block的博客导航切换首页文章command_block的博客多项式计数杂谈postedon2020-02-1118:16:01|under算法|Ox00.前面的话&前置芝士本文Markd......
  • 多项式求值软件下载Polynomial evaluation software mus 2025 download
    本软件是Windows下64位软件。本软件能计算如a0+a1*x+a2*x^2+......+an*x^n的式子的对b1的求值结果。具体的方法就是在多多项式系数区输入a0到an的值,然后点击计算多项式的结果即可在结果栏算出结果。最大项数为1000项。多项式系数输入时1项1行,从上到下是a0到an,中间不能空行。T......
  • 切比雪夫多项式
    切比雪夫多项式通常我们使用切比雪夫多项式时都在范围[-1,1]之间。定义切比雪夫多项式在[-1,1]上的定义是:\(T_n(x)=cos(narccos(x)),-1\leqx\leq1\),其中,T_n(x)是阶数为n的切比雪夫多项式。性质\(T_n(x)\)是n阶多项式。\(T_n(x)\)的奇偶性和n的奇偶性一致。\(T_n(x)\)在区......
  • 手机全家桶买哪个品牌好?
    手机全家桶买哪个品牌好?选择手机全家桶品牌时,通常要考虑多个因素,包括品牌的声誉、产品质量、性能、价格、生态系统和客户支持等。以下是几个备受认可的手机全家桶品牌:1.Apple(苹果):苹果提供一系列无缝集成的产品,包括iPhone手机、iPad平板电脑、Mac电脑和AppleWatch智能手表。它......
  • 多项式基础
    FFT/NTT板子,就不说了。分治FFT给定序列\(g_{1\dotsn-1}\),求序列\(f_{0\dotsn-1}\)。其中\(f_i=\sum_{j=1}^if_{i-j}g_j\),边界为\(f_0=1\)。Solution1朴素地做是\(\mathcal{O}(n^2)\)的。观察到后面的项依赖前面的项。于是我们可以考虑用cdq分治实现。具体......
  • 多项式定积分计算软件2025 64位WIN版下载Polynomial definite integral calculation s
    本软件功能强大,价格实惠,欢迎试用本软件的WIN64位版本。本软件能计算如a0+a1*x+a2*x^2+......+an*a^n的式子的对b1和b2的积分的结果。具体的方法就是在多多项式系数区输入a0到an的值,然后点击计算积分结果即可在结果栏算出结果。最大项数为1000项。多项式系数输入时1项1行,从上......