文档章节

多项式各种操作

o
 osc_fmg49rzg
发布于 2019/03/20 12:06
字数 2513
阅读 7
收藏 0

精选30+云产品,助力企业轻松上云!>>>

多项式各种操作

这里的操作还不是特别完备,可能以后会更新??gugugu

多项式求逆

给定一个多项式$A(x)$,求多项式$F(x)$满足: $$ F(x)A(x)\equiv 1\pmod{x^n} $$ 当然如果是良心出题人会让你模$998244353$之类的模数。


首先要把$n$补成$2$的次幂,好理解也好写一些。

假设我们现在已经求出了$F_0(x)$,满足: $$ F_0(x)A(x)\equiv 1\pmod{ x^{n/2}} $$ 现在要求$F(x)$,满足: $$ F(x)A(x)\equiv 1\pmod{x^n} $$ 显然也满足: $$ F(x)A(x)\equiv 1\pmod{x^{n/2}} $$ 那么两式相减可得: $$ F(x)-F_0(x)\equiv0\pmod{x^{n/2}} $$ 两边平方: $$ F^2(x)-2F(x)F_0(x)+F_0^2(x)\equiv0\pmod{x^{n}} $$ 注意这里模数也平方了,因为左边$<n/2$的每一项都是$0$,平方之后,对于第$i$项,$i\in [n/2,n)$,可得:$a_i=\sum_{j=0}^i a_ja_{i-j}$,注意这里$j$和$i-j$必然有一个$<n/2$,也就是说这一项也为$0$,所以是正确的。

那么两边都乘上一个$A(x)$: $$ F(x)\equiv2F_0(x)-A(x)F_0^2(x)\pmod{x^n} $$ 那么我们只要递归做就好了。

复杂度$T(n)=T(n/2)+O(n\log n)$,所以$T(n)=O(n\log n)$。常数可想而知

代码大概长这个样子:

void get_inv(int *r,int *t,int m) {
	if(m==1) return t[0]=qpow(r[0],mod-2),void();
	get_inv(r,t,(m+1)>>1);
	N=1,bit=0;while(N<(m<<1)) N<<=1,bit++;
	for(int i=0;i<N;i++) pos[i]=pos[i>>1]>>1|((i&1)<<(bit-1));
	for(int i=0;i<m;i++) tmp[i]=r[i];
	for(int i=m;i<N;i++) tmp[i]=0;
	ntt(tmp,1),ntt(t,1);
	for(int i=0;i<N;i++) t[i]=1ll*t[i]*(2%mod-1ll*tmp[i]*t[i]%mod+mod)%mod;
	ntt(t,-1);for(int i=m;i<N;i++) t[i]=0;
}

完整代码会贴在下面。

多项式求ln

给定一个多项式$A(x)$,求$\ln A(x)$,保证$[0]A(x)=1$。


多项式的对数...大概可以如下定义: $$ \ln(1-F(x))=\sum_{n=1}^\infty -\frac{F^n(x)}{n} $$ 其实就是把这玩意泰勒展开了一下。


那么计算其实很简单了,令: $$ F(x)=\ln A(x) $$ 两边对$x$求导: $$ F'(x)=\frac{A'(x)}{A(x)} $$ 然后注意到后面是可以算出来的,直接多项式求逆就好了。

然后两边积分: $$ F(x)=\int \frac{A'(x)}{A(x)} {\rm{d}}x $$ 注意到求导和积分都是可以$O(n)$算的,所以总复杂度$O(n\log n)$。

牛顿迭代

给定一个函数$G(x)$,求多项式$F(x)$,使得$G(F(x))=0$。


假定我们现在已经求出了$F_0(x)$,满足:

$$ G(F_0(x))\equiv0\pmod{x^{n/2}} $$

考虑对$G(F(x))$在$F_0(x)$的泰勒展开: $$ G(F(x))\equiv\sum_{n=0}^\infty G^{(n)}(F_0(x))(F(x)-F_0(x))^n\pmod{x^n} $$ 这里$G^{(n)}$表示$n$阶导,注意这里是对$F(x)$求导,也就是: $$ \frac{{\rm{d}}^nG(F(x))}{{\rm{d}}F(x)^n} $$ 注意到$F(x)$和$F_0(x)$的前$n/2$项相同,那么可以发现: $$ (F(x)-F_0(x))^2\equiv0\pmod{x^n} $$ 所以后面的项都消掉了,剩下前两项: $$ G(F(x))=G(F_0(x))+G'(F_0(x))(F(x)-F_0(x))=0 $$ 然后可以把$F(x)$解出来: $$ F(x)=F_0(x)-\frac{G(F_0(x))}{G'(F_0(x))} $$ 然后就做完了。

下面举两个栗子看看这玩意有啥用。

多项式求exp

给定一个多项式$A(x)$,求$\exp(A(x))$,保证$[0]A(x)=0$。


同样的,用麦克劳林级数定义如下: $$ e^{A(x)}=\sum_{n=0}^\infty \frac{F^n(x)}{n!} $$


现在是要求$F(x)$满足: $$ e^{A(x)}=F(x) $$ 两边取$\ln$: $$ A(x)=\ln F(x)\ \ln F(x)-A(x)=0 $$ 构造函数$G(x)$如下: $$ G(F(x))=\ln F(x)-A(x) $$ 那么: $$ G'(F(x))=\frac{1}{F(x)} $$ 带入先前的牛顿迭代式子: $$ F(x)=F_0(x)-(\ln F_0(x)-A(x))F_0(x)\ =F_0(x)(1-\ln F_0(x)+A(x)) $$ 带进去算就好了,复杂度$O(n\log n)$,这里面还套了个多项式求逆,常数可想而知

多项式开根

给定一个多项式$A(x)$,求$F(x)$,满足$F^2(x)\equiv A(x)\pmod{x^n}$。


同样构造函数$G(x)$: $$ G(F(x))=F^2(x)-A(x)=0 $$ 那么带进去: $$ \begin{align} F(x)&=F_0(x)-\frac{F_0^2(x)-A(x)}{2F_0(x)}=\frac{F_0^2(x)+A(x)}{2F_0(x)}\ &=\frac{1}{2}(F_0(x)+\frac{A(x)}{F_0(x)}) \end{align} $$ 复杂度同样$O(n\log n)$。

多项式快速幂

给定一个多项式$A(x)$,求$F(x)$,满足$F(x)\equiv A^k(x)\pmod{x^n}$。


两边取$\ln$: $$ \ln F(x)=k\ln A(x) $$ 然后$\exp$回来: $$ F(x)=\exp(k\ln A(x)) $$ 复杂度$O(n\log n)$。

这个也可以实现多项式开根,只是跑的慢一些而已。

模板

多项式$\exp$:题目来自洛谷模板

#include<bits/stdc++.h>
using namespace std;
 
void read(int &x) {
    x=0;int f=1;char ch=getchar();
    for(;!isdigit(ch);ch=getchar()) if(ch=='-') f=-f;
    for(;isdigit(ch);ch=getchar()) x=x*10+ch-'0';x*=f;
}
 
void print(int x) {
    if(x<0) putchar('-'),x=-x;
    if(!x) return ;print(x/10),putchar(x%10+48);
}
void write(int x) {if(!x) putchar('0');else print(x);putchar('\n');}

const int maxn = 1e6+10;
const int mod = 998244353;

int s[maxn],pos[maxn],N,bit,n,s2[maxn];

int qpow(int a,int x) {
	int res=1;
	for(;x;x>>=1,a=1ll*a*a%mod) if(x&1) res=1ll*res*a%mod;
	return res;
}

void ntt(int *r,int op) {
	for(int i=0;i<N;i++) if(pos[i]>i) swap(r[i],r[pos[i]]);
	for(int i=1;i<N;i<<=1) {
		int wn=qpow(op==1?3:qpow(3,mod-2),(mod-1)/(i<<1));
		for(int j=0;j<N;j+=(i<<1)) 
			for(int k=0,w=1;k<i;k++,w=1ll*w*wn%mod) {
				int x=r[j+k],y=1ll*w*r[i+j+k]%mod;
				r[j+k]=(x+y)%mod,r[i+j+k]=(x-y+mod)%mod;
			}
	}
	if(op==-1) {
		int inv=qpow(N,mod-2);
		for(int i=0;i<N;i++) r[i]=1ll*r[i]*inv%mod;
	}
}

int tmp[maxn],tmp2[maxn],tmp3[maxn];

void get_inv(int *r,int *t,int m) {
	if(m==1) return t[0]=qpow(r[0],mod-2),void();
	get_inv(r,t,(m+1)>>1);
	N=1,bit=0;while(N<(m<<1)) N<<=1,bit++;
	for(int i=0;i<N;i++) pos[i]=pos[i>>1]>>1|((i&1)<<(bit-1));
	for(int i=0;i<m;i++) tmp[i]=r[i];
	for(int i=m;i<N;i++) tmp[i]=0;
	ntt(tmp,1),ntt(t,1);
	for(int i=0;i<N;i++) t[i]=1ll*t[i]*(2%mod-1ll*tmp[i]*t[i]%mod+mod)%mod;
	ntt(t,-1);for(int i=m;i<N;i++) t[i]=0;
}

void integral(int *r,int *t,int m) {
	N=1;while(N<(m<<1)) N<<=1;t[0]=0;
	for(int i=0;i<m;i++) t[i+1]=1ll*r[i]*qpow(i+1,mod-2)%mod;
	for(int i=m+1;i<N;i++) t[i]=0;
}

void derivative(int *r,int *t,int m) {
	N=1;while(N<(m<<1)) N<<=1;
	for(int i=1;i<m;i++) t[i-1]=1ll*r[i]*i%mod;
	for(int i=m-1;i<N;i++) t[i]=0;
}

void get_ln(int *r,int *t,int m) {
	derivative(r,tmp3,m);
	get_inv(r,tmp2,m);
	ntt(tmp2,1),ntt(tmp3,1);
	for(int i=0;i<N;i++) tmp2[i]=1ll*tmp2[i]*tmp3[i]%mod;
	ntt(tmp2,-1);integral(tmp2,t,m);
	for(int i=0;i<=m<<1;i++) tmp2[i]=0;
}

int t2[maxn];

void get_exp(int *r,int *t,int m) {
	if(m==1) return t[0]=1,void();
	get_exp(r,t,(m+1)>>1);	
	N=1,bit=0;while(N<(m<<1)) N<<=1,bit++;
	for(int i=m;i<N;i++) t[i]=0;
	get_ln(t,t2,m);
	for(int i=m;i<N;i++) t2[i]=0;
	for(int i=0;i<m;i++) t2[i]=(-t2[i]+r[i]+mod)%mod;
	t2[0]=(t2[0]+1)%mod;
	for(int i=0;i<N;i++) pos[i]=pos[i>>1]>>1|((i&1)<<(bit-1));
	ntt(t2,1),ntt(t,1);
	for(int i=0;i<N;i++) t[i]=1ll*t[i]*t2[i]%mod;
	ntt(t,-1);for(int i=m;i<N;i++) t[i]=0;
}

int main() {
	read(n);
	for(int i=0;i<n;i++) read(s[i]);
	get_exp(s,s2,n);
	for(int i=0;i<n;i++) printf("%d ",(s2[i]+mod)%mod);puts("");
	return 0;
}

多项式开根:题目来自洛谷模板

#include<bits/stdc++.h>
using namespace std;
 
void read(int &x) {
    x=0;int f=1;char ch=getchar();
    for(;!isdigit(ch);ch=getchar()) if(ch=='-') f=-f;
    for(;isdigit(ch);ch=getchar()) x=x*10+ch-'0';x*=f;
}
 
void print(int x) {
    if(x<0) putchar('-'),x=-x;
    if(!x) return ;print(x/10),putchar(x%10+48);
}
void write(int x) {if(!x) putchar('0');else print(x);putchar('\n');}

#define lf double
#define ll long long 

const int maxn = 6e5+10;
const int inf = 1e9;
const lf eps = 1e-8;
const int mod = 998244353;
const int inv2 = 499122177;

int f[maxn],g[maxn],n,m,mxn,bit,N,w[maxn],rw[maxn],s[maxn],t[maxn],pos[maxn];

int qpow(int a,int x) {
	int res=1;
	for(;x;x>>=1,a=1ll*a*a%mod) if(x&1) res=1ll*res*a%mod;
	return res;
}

void prepare() {
	w[0]=1;w[1]=qpow(3,(mod-1)/mxn);
	for(int i=2;i<=mxn;i++) w[i]=1ll*w[i-1]*w[1]%mod;
	rw[0]=1,rw[1]=qpow(qpow(3,mod-2),(mod-1)/mxn);
	for(int i=2;i<=mxn;i++) rw[i]=1ll*rw[i-1]*rw[1]%mod;
}

void ntt(int *r,int op) {
	for(int i=1;i<N;i++) if(pos[i]>i) swap(r[i],r[pos[i]]);
	for(int i=1,d=mxn>>1;i<N;i<<=1,d>>=1) 
		for(int j=0;j<N;j+=i<<1)
			for(int k=0;k<i;k++) {
				int x=r[j+k],y=1ll*r[i+j+k]*(op==1?w:rw)[k*d]%mod;
				r[j+k]=(x+y)%mod,r[i+j+k]=(x-y+mod)%mod;
				}
	if(op==-1) {
		int inv=qpow(N,mod-2);
		for(int i=0;i<N;i++) r[i]=1ll*r[i]*inv%mod;
	}
}

int tmp1[maxn],tmp2[maxn],tmp3[maxn];

void get_pos(int len) {
	for(bit=0,N=1;N<len;N<<=1,bit++);
	for(int i=1;i<N;i++) pos[i]=pos[i>>1]>>1|((i&1)<<(bit-1));
}

void poly_inv(int *r,int *b,int len) {
	if(len==1) return b[0]=qpow(r[0],mod-2),void();
	poly_inv(r,b,len>>1);
	for(int i=0;i<len;i++) tmp1[i]=b[i],tmp2[i]=r[i];
	get_pos(len<<1);
	ntt(tmp1,1),ntt(tmp2,1);
	for(int i=0;i<N;i++) b[i]=((2ll*tmp1[i]%mod-1ll*tmp2[i]*tmp1[i]%mod*tmp1[i]%mod)%mod+mod)%mod;
	ntt(b,-1);
	for(int i=len;i<N;i++) b[i]=0;
	for(int i=0;i<len<<1;i++) tmp1[i]=tmp2[i]=0;
}

void poly_sqrt(int *r,int *b,int len) {
	if(len==1) return b[0]=r[0],void();
	poly_sqrt(r,b,len>>1);
	poly_inv(b,tmp3,len);
	get_pos(len<<1);
	for(int i=0;i<len;i++) tmp2[i]=r[i];
	ntt(tmp2,1),ntt(tmp3,1);
	for(int i=0;i<N;i++) tmp3[i]=1ll*tmp3[i]*tmp2[i]%mod;
	ntt(tmp3,-1);
	for(int i=0;i<len;i++) b[i]=1ll*inv2*(b[i]+tmp3[i])%mod;
	for(int i=0;i<len<<1;i++) tmp3[i]=tmp2[i]=0;
}

int main() {
	read(n);
	for(int i=0;i<n;i++) read(g[i]);
	for(mxn=1;mxn<=n<<1;mxn<<=1);
	prepare();
	poly_sqrt(g,t,mxn>>1);
	for(int i=0;i<n;i++) printf("%d ",t[i]);puts("");
	return 0;
}

代码是不同时期写的,可能相同的函数有点出入。

o
粉丝 0
博文 500
码字总数 0
作品 0
私信 提问
加载中
请先登录后再评论。
数据结构上机实验(2)

1、顺序表的各种基本运算操作 2、单链表的各种基本运算算法 3、双链表的各种基本运算算法 4、循环单链表的各种基本运算操作 5、循环双链表的各种基本运算操作 6、将单链表按基准划分 7、将两...

osc_f5ujy662
05/24
10
0
数据结构上机实验(2)

1、顺序表的各种基本运算操作 2、单链表的各种基本运算算法 3、双链表的各种基本运算算法 4、循环单链表的各种基本运算操作 5、循环双链表的各种基本运算操作 6、将单链表按基准划分 7、将两...

我在吃大西瓜呢
05/23
0
0
多项式的各种操作

多项式求逆 给定多项式$A(x)$,求一个多项式$B(x)$满足$A(x)B(x)=1left(text{mod} x^nright)$。 假设已求出多项式$C(x)$满足$A(x)C(x)=1left(text{mod} x^{lceilfrac{n}{2}rceil}right)$。 ...

osc_y0xqgfqd
2019/03/14
2
0
第二周(多变量线性回归 +Matlab使用)-【机器学习-Coursera Machine Learning-吴恩达】

目录: 多变量线性回归(模型、梯度下降技巧) 特征选择和多项式回归 正规方程 Matlab学习 1 多变量线性回归1)模型 - 假设函数: - 参数:全部的 theta - 代价函数:和单变量回归一样 - 梯度...

kevinbetterq
2018/03/05
0
0
STARKs, Part II: Thank Goodness It's FRI-day

在本系列的上一篇文章中,我们谈到了,如何能够做出一些非常有意思且简洁的计算证明,比如通过利用多项式复合和除法技术,证明你算出了第一百万个斐波那契数。但是,它依托于一个非常重要的元...

liuchengxu
2017/12/23
0
0

没有更多内容

加载失败,请刷新页面

加载更多

Asp.net core之NLog

NuGet添加 NLog.Web.AspNetCore。 <PackageReference Include="Microsoft.AspNetCore.App" /> 添加配置文件 新建一个文件nlog.config(建议全部小写,linux系统中要注意), 并右键点击其属性......

一介草民Coder
13分钟前
7
0
.NET中的struct和class有什么区别? - What's the difference between struct and class in .NET?

问题: .NET中的struct和class有什么区别? 解决方案: 参考一: https://stackoom.com/question/3OT/NET中的struct和class有什么区别 参考二: https://oldbug.net/q/3OT/What-s-the-differ...

富含淀粉
55分钟前
23
0
android:layout_weight是什么意思? - What does android:layout_weight mean?

问题: I don't understand how to use this attribute. 我不明白如何使用这个属性。 Can anyone tell me more about it? 谁能告诉我更多关于它的事情? 解决方案: 参考一: https://stacko...

javail
今天
17
0
CSS背景不透明度[重复] - CSS Background Opacity [duplicate]

问题: This question already has an answer here: 这个问题已经在这里有了答案: How do I give text or an image a transparent background using CSS? 如何使用CSS为文本或图像提供透明背...

fyin1314
今天
31
0
node http 获取gb2312网页如何转为utf8

最初,我想当然认为是下述做法,但被证明是错误的 const http = require('http'), iconv = require('iconv-lite');const url = 'http://xxx';http.get(url, function(res) { var bo......

高延
今天
24
0

没有更多内容

加载失败,请刷新页面

加载更多

返回顶部
顶部