zoukankan      html  css  js  c++  java
  • 多项式ln、牛顿迭代学习笔记

    这篇博客写得很好,本文大部分参考此博客。

    一些证明的细节在本文中不会提及。

    前置知识:多项式求逆,多项式求导、积分 , 泰勒展开。

    多项式ln

    给定(G(x)),求

    [F(x)=ln G(x) ]

    两边求导,得

    [F'(x)=frac{G'(x)}{G(x)} ]

    [F(x)=int frac{G'(x)}{G(x)} ]

    求逆,计算即可。

    牛顿迭代

    以下我们把一个多项式当做一个数来理解。

    (delta(F(x)))是一个关于多项式(F(x))的一个函数,现在我们要求一个(F(x)),满足

    [delta(F(x))=0 ]

    假设我们已经求出了(F(x))的前(frac{n}{2})项,记为(F_0(x)),现在我们要求它的前(n)

    (delta(F(x)))进行泰勒展开,有

    [delta(F(x))=delta(F_0(x))+frac{delta '(F_0(x))}{1!}(F(x)-F_0(x))+frac{delta ''(F_0(x))}{2!}(F(x)-F_0(x))^2+... ]

    其中(delta '(F(x)))(delta (F(x)))关于(F(x))的导数

    因为(F_0(x))(F(x))(frac{n}{2})项相同,所以((F(x)-F_0(x))^2equiv 0 (mod x^n)),所以上式从第三项开始的值全为0。

    也就是

    [delta(F(x))=delta(F_0(x))+delta '(F_0(x))(F(x)-F_0(x))=0 ]

    整理得

    [F(x)=F_0(x)-frac{delta(F_0(x))}{delta '(F_0(x))} ]

    如果没搞懂的话可以参照下面的例子来理解。

    牛顿迭代——多项式求逆

    给定(G(x)),求(F(x))满足

    [F(x)G(x)equiv1 (mod x^n) ]

    就是要使得(F(x)G(x)-1=0)

    所以设 (delta(F(x))=F(x)G(x)- 1)

    就有

    [F(x)=F_0(x)-frac{delta(F_0(x))}{delta '(F_0(x))} ]

    [=F_0(x)-frac{F_0(x)G(x)-1}{G(x)} ]

    [=F_0(x)-F_0(F_0(x)G(x)-1) ]

    [=2F_0(x)-F_0^2(x)G(x) (mod x^n) ]

    牛顿迭代——多项式开根

    即求(F^2(x)-G(x)=0)

    [F(x)=F_0(x)-frac{F_0^2(x)-G(x)}{2F_0(x)} ]

    [=frac{F_0^2(x)+G(x)}{2F_0(x)} ]

    多项式Exp

    [F(x)=e^ {G(x)} ]

    [ln F(x)=G(x) ]

    [ln F(x)-G(x)=0 ]

    [F(x)=F_0(x)-frac{ln F_0(x)-G(x)}{frac{1}{F_0(x)}} ]

    [=F_0(x)-F_0(x)(ln F_0(x)-G(x)) ]

    [=F_0(x)(1-ln F_0(x)+G(x)) ]

    最后奉上我丑陋的代码模板:

    #include<bits/stdc++.h>
    using namespace std;
    #define N 600007
    #define M 2100007
    #define ll long long
    const ll mod=998244353;
    const int lim=4e5;
    ll inv[N];
    int rev[M],len;
    ll p2(ll x){return x*x%mod;}
    ll pw(ll x,ll p)
    {
    	return p?p2(pw(x,p/2))*(p&1?x:1)%mod:1;
    }
    void getrev()
    {
    	for(int i=0;i<len;i++)
    		rev[i]=rev[i>>1]>>1|(i&1?len>>1:0);
    }
    void getlen(int n)
    {
    	for(len=1;len<=n;len<<=1);
    	getrev();
    }
    void NTT(ll *a,int op)
    {
    	for(int i=0;i<len;i++)if(i<rev[i])swap(a[i],a[rev[i]]);
    	for(int i=1;i<len;i<<=1)
    	{
    		ll nw=pw(3,(mod-1)/(i<<1));
    		for(int j=0;j<len;j+=i<<1)
    		{
    			ll w=1;
    			for(int k=j;k<j+i;k++)
    			{
    				ll x=a[k],y=a[k+i]*w%mod;
    				a[k]=(x+y)%mod,a[k+i]=(x-y+mod)%mod;
    				w=w*nw%mod;
    			}
    		}
    	}
    	if(op<0)
    	{
    		reverse(a+1,a+len);
    		ll Inv=pw(len,mod-2);
    		for(int i=0;i<len;i++)a[i]=a[i]*Inv%mod;
    	}
    }
    ll inv_c[N];
    void getinv(int p,ll *a,ll *b)
    {
    	if(p==1)return a[0]=pw(b[0],mod-2),(void)1;
    	getinv((p+1)/2,a,b);
    	getlen(p<<1);
    	ll *c=inv_c;
    	for(int i=0;i<p;i++)c[i]=b[i];
    	for(int i=p;i<len;i++)c[i]=0;
    	NTT(a,1),NTT(c,1);
    	for(int i=0;i<len;i++)a[i]=a[i]*(2-a[i]*c[i]%mod+mod)%mod;
    	NTT(a,-1);
    	for(int i=p;i<len;i++)a[i]=0;
    }
    ll mul_c[M],mul_d[M];
    void mul(ll *t,ll *a,ll *b)
    {
    	ll *c=mul_c,*d=mul_d;
    	for(int i=0;i<len;i++)c[i]=a[i],d[i]=b[i];
    	NTT(c,1),NTT(d,1);
    	for(int i=0;i<len;i++)c[i]=c[i]*d[i]%mod;
    	NTT(c,-1);
    	for(int i=0;i<len;i++)t[i]=c[i];
    }
    void devir(ll *a,int n)
    {
    	for(int i=1;i<=n;i++)a[i-1]=a[i]*i%mod;
    	a[n]=0;
    }
    void inter(ll *a,int n)
    {
    	for(int i=n;i>=0;i--)a[i+1]=a[i]*inv[i+1]%mod;
    	a[0]=0;
    }
    ll ln_c[N],ln_d[N];
    void getln(int n,ll *a,ll *b)
    {
    	ll *c=ln_c,*d=ln_d;
    	getlen(2*n);
    	for(int i=0;i<len;i++)c[i]=b[i],d[i]=0;
    	getinv(n,d,c),devir(c,n);
    	getlen(2*n);
    	mul(a,c,d);
    	inter(a,n);
    	for(int i=n;i<len;i++)a[i]=0;
    }
    ll exp_c[N];
    void getexp(int p,ll *a,ll *b)
    {
    	if(p==1)return a[0]=1,(void)1;
    	getexp((p+1)/2,a,b);
    	ll *c=exp_c;
    	getlen(2*p);
    	for(int i=0;i<len;i++)c[i]=0;
    	getln(p,c,a);
    	for(int i=0;i<p;i++)c[i]=(b[i]-c[i]+mod)%mod;
    	c[0]=(c[0]+1)%mod;
    	getlen(2*p);
    	mul(a,a,c);
    	for(int i=p;i<len;i++)a[i]=0;
    }
    ll sqrt_c[N],sqrt_d[N];
    void getsqrt(int p,ll *a,ll *b)
    {
    	if(p==1)return a[0]=1,(void)1;
    	getsqrt((p+1)/2,a,b);
    	ll *c=sqrt_c,*d=sqrt_d;
    	getlen(2*p);
    	for(int i=0;i<len;i++)c[i]=d[i]=0;
    	for(int i=0;i<p;i++)d[i]=b[i];
    	getinv(p,c,a);
    	getlen(2*p),mul(c,c,d);
    	for(int i=0;i<p;i++)a[i]=(a[i]+c[i])*inv[2]%mod;
    }
    ll a[M],b[M],c[M];
    int n;
    void work_mul()
    {
    	int n,m;
    	scanf("%d%d",&n,&m);
    	for(int i=0;i<=n;i++)scanf("%lld",&a[i]);
    	for(int i=0;i<=m;i++)scanf("%lld",&b[i]);
    	getlen(n+m);
    	mul(c,a,b);
    	for(int i=0;i<=n+m;i++)printf("%lld ",c[i]);
    }
    void read()
    {
    	scanf("%d",&n);
    	for(int i=0;i<n;i++)scanf("%lld",&a[i]);
    }
    void write(ll *a)
    {
    	for(int i=0;i<n;i++)printf("%lld ",a[i]);
    }
    void work_inv()
    {
    	read();
    	getinv(n,b,a);
    	write(b);
    }
    void work_sqrt()
    {
    	read();
    	getsqrt(n,b,a);
    	write(b);
    }
    void work_ln()
    {
    	read();
    	getln(n,b,a);
    	write(b);
    }
    void work_exp()
    {
    	read();
    	getexp(n,b,a);
    	write(b);
    }
    void Init()
    {
    	inv[1]=1;
    	for(int i=2;i<=lim;i++)
    		inv[i]=mod-mod/i*inv[mod%i]%mod;
    }
    int main()
    {
    	//freopen("data.in","r",stdin);
    	Init();
    	work_exp();
    	return 0;
    }
    
    
    
  • 相关阅读:
    LNMP安装后MYSQL数据库无法远程访问解决
    Ubuntu忘记root密码怎么办?
    composer安装出现proc_open没有开启问题的解决方案
    LNMP搭建环境遇到的N多坑
    lnmp HTTP ERROR 500
    LNMP集成运行(开发)环境的部署
    最新javamail 使用方案,可以异步发送邮件
    vi常用快捷键
    Dom4j解析XML文件
    Multiple markers at this line @Override的解决方法
  • 原文地址:https://www.cnblogs.com/lishuyu2003/p/12123705.html
Copyright © 2011-2022 走看看