[luoguP4705]玩游戏

最新推荐文章于 2019-02-26 20:44:13 发布

zhouyuheng2003

最新推荐文章于 2019-02-26 20:44:13 发布

阅读量218

点赞数

分类专栏： OI NTT 分治NTT

本文链接：https://blog.csdn.net/zhouyuheng2003/article/details/86695870

版权

OI 同时被 3 个专栏收录

104 篇文章 2 订阅

订阅专栏

NTT

11 篇文章 0 订阅

订阅专栏

分治NTT

3 篇文章 0 订阅

订阅专栏

前言

好的思路题

题目相关

传送门

题意简介

给定长度为 $n$ 的序列 $a$ 和长度为 $m$ 的序列 $b$
对所有 $k\in[1,t]$ ，求在两个序列里分别随机取一个数的和的 $k$ 次方的期望，模 $998244353$

数据范围

$1\le n,m\le10^5,1\le t\le 10^5$

题解

化一下式子，求期望相当于求所有情况的和取平均值
$\begin{aligned} \sum(a+b)^k&=\sum\sum_{i=0}^k\binom{k}{i}a^ib^{k-i}\\ &=\sum_{i=0}^k\binom{k}{i}(\sum a^i)(\sum b^{k-i})\\ &=\sum_{i=0}^k\binom{k}{i}(\sum a^i)(\sum b^{k-i})\\ &=\sum_{i=0}^k\frac{k!}{i!(k-i)!}(\sum a^i)(\sum b^{k-i})\\ &=\sum_{i=0}^kk!\frac{\sum a^i}{i!}\frac{\sum b^{k-i}}{(k-i)!}\\ &=k!\sum_{i=0}^k\frac{\sum a^i}{i!}\frac{\sum b^{k-i}}{(k-i)!}\\ \end{aligned}$
我们如果求出所有 $\sum a^i$ 和 $\sum b^i$
考虑生成函数
$F_i(x)=\sum_{j=0}^{\infty}a_i^jx^j$
根据等比数列求和得到
$F_i(x)=\frac{1}{1-a_ix}$
我们对所有生成函数求和
$\sum_{i=1}^nF_i(x)=\sum_{i=1}^n\frac{1}{1-a_ix}$
然后分治求分子、分母，用多项式求逆将分母乘进分子里
上面的式子算完后别忘了由于是期望所以要除以n*m
算法总复杂度 $\mathcal O(nlog^2n)$

代码

这份代码优化了DFT与IDFT的次数，速度大大提高

#include<cstdio>
#include<cctype>
#include<cstring>
#include<cstdlib>
#include<vector>
namespace fast_IO
{
    const int IN_LEN=10000000,OUT_LEN=10000000;
    char ibuf[IN_LEN],obuf[OUT_LEN],*ih=ibuf+IN_LEN,*oh=obuf,*lastin=ibuf+IN_LEN,*lastout=obuf+OUT_LEN-1;
    inline char getchar_(){return (ih==lastin)&&(lastin=(ih=ibuf)+fread(ibuf,1,IN_LEN,stdin),ih==lastin)?EOF:*ih++;}
    inline void putchar_(const char x){if(oh==lastout)fwrite(obuf,1,oh-obuf,stdout),oh=obuf;*oh++=x;}
    inline void flush(){fwrite(obuf,1,oh-obuf,stdout);}
}
using namespace fast_IO;
#define getchar() getchar_()
#define putchar(x) putchar_((x))
typedef long long ll;
#define rg register
template <typename T> inline T max(const T a,const T b){return a>b?a:b;}
template <typename T> inline T min(const T a,const T b){return a<b?a:b;}
template <typename T> inline T mind(T&a,const T b){a=a<b?a:b;}
template <typename T> inline T maxd(T&a,const T b){a=a>b?a:b;}
template <typename T> inline T abs(const T a){return a>0?a:-a;}
template <typename T> inline void swap(T&a,T&b){T c=a;a=b;b=c;}
template <typename T> inline void swap(T*a,T*b){T c=a;a=b;b=c;}
template <typename T> inline T gcd(const T a,const T b){if(!b)return a;return gcd(b,a%b);}
template <typename T> inline T square(const T x){return x*x;};
template <typename T> inline void read(T&x)
{
    char cu=getchar();x=0;bool fla=0;
    while(!isdigit(cu)){if(cu=='-')fla=1;cu=getchar();}
    while(isdigit(cu))x=x*10+cu-'0',cu=getchar();
    if(fla)x=-x;  
}
template <typename T> void printe(const T x)
{
    if(x>=10)printe(x/10);
    putchar(x%10+'0');
}
template <typename T> inline void print(const T x)
{
    if(x<0)putchar('-'),printe(-x);
    else printe(x);
}
const int maxn=524288,mod=998244353;
inline int Md(const int x){return x>=mod?x-mod:x;}
template<typename T>
inline int pow(int x,T y)
{
    rg int res=1;x%=mod;
    for(;y;y>>=1,x=(ll)x*x%mod)if(y&1)res=(ll)res*x%mod;
    return res;
}
namespace Poly///namespace of Poly
{
int W_[maxn],FW_[maxn],ha[maxn],hb[maxn],Inv[maxn];
inline void init(const int x)
{
    rg int tim=0,lenth=1;
    while(lenth<x)lenth<<=1,tim++;
    for(rg int i=1;i<lenth;i++)
    {
        W_[i]=pow(3,(mod-1)/i/2);
        FW_[i]=pow(W_[i],mod-2);
    }
//	Inv[0]=Inv[1]=1;for(rg int i=2;i<x;i++)Inv[i]=(ll)(mod-mod/i)*Inv[mod%i]%mod;
}
int L;
inline void NTT(int*A,const int fla)//prepare:init L 
{
    for(rg int i=0,j=0;i<L;i++)
    {
        if(i>j)swap(A[i],A[j]);
        for(rg int k=L>>1;(j^=k)<k;k>>=1);
    }
    for(rg int i=1;i<L;i<<=1)
    {
        const int w=fla==-1?FW_[i]:W_[i];
        for(rg int j=0,J=i<<1;j<L;j+=J)
        {
            int K=1;
            for(rg int k=0;k<i;k++,K=(ll)K*w%mod)
            {
                const int x=A[j+k],y=(ll)A[j+k+i]*K%mod;
                A[j+k]=Md(x+y),A[j+k+i]=Md(mod+x-y);
            }
        }
    }
}
inline int Quadratic_residue(const int a)
{
    if(a==0)return 0;
    int b=(rand()<<14^rand())%mod;
    while(pow(b,(mod-1)>>1)!=mod-1)b=(rand()<<14^rand())%mod;
    int s=mod-1,t=0,x,inv=pow(a,mod-2),f=1;
    while(!(s&1))s>>=1,t++,f<<=1;
    t--,x=pow(a,(s+1)>>1),f>>=1;
    while(t)
    {
        f>>=1;
        if(pow((int)((ll)inv*x%mod*x%mod),f)!=1)x=(ll)x*pow(b,s)%mod;
        t--,s<<=1;
    }
    return min(x,mod-x);
}
struct poly
{
    std::vector<int>A;
    poly(){A.resize(0);}
    poly(const int x){A.resize(1),A[0]=x;}
    inline int&operator[](const int x){return A[x];}
    inline int operator[](const int x)const{return A[x];}
    inline void clear(){A.clear();}
    inline unsigned int size()const{return A.size();}
    inline void resize(const unsigned int x){A.resize(x);}
    void RE(const int x)
    {
        A.resize(x);
        for(rg int i=0;i<x;i++)A[i]=0; 
    }
    void readin(const int MAX)
    {
        A.resize(MAX);
        for(rg int i=0;i<MAX;i++)read(A[i]);
    }
    void putout()const
    {
        for(rg unsigned int i=0;i<A.size();i++)print(A[i]),putchar(' ');
    }
    inline poly operator +(const poly b)const
    {
        poly RES;
        RES.resize(max(size(),b.size()));
        for(rg unsigned int i=0;i<RES.size();i++)RES[i]=Md((i<size()?A[i]:0)+(i<b.size()?b[i]:0));
        return RES;
    }
    inline poly operator -(const poly b)const
    {
        poly RES;
        RES.resize(max(size(),b.size()));
        for(rg unsigned int i=0;i<RES.size();i++)RES[i]=Md((i<size()?A[i]:0)+mod-(i<b.size()?b[i]:0));
        return RES;
    }
    inline poly operator *(const int b)const
    {
        poly RES=*this;
        for(rg unsigned int i=0;i<RES.size();i++)RES[i]=(ll)RES[i]*b%mod;
        return RES;
    }
    inline poly operator /(const int b)const
    {
        poly RES=(*this)*pow(b,mod-2);
    	return RES;
    }
    inline poly operator *(const poly b)const
    {
        const int RES=A.size()+b.size()+1;
        L=1;while(L<RES)L<<=1;
        poly c;c.A.resize(RES);
        memset(ha,0,sizeof(int)*L);
        memset(hb,0,sizeof(int)*L);
        for(rg unsigned int i=0;i<A.size();i++)ha[i]=A[i];
        for(rg unsigned int i=0;i<b.A.size();i++)hb[i]=b.A[i];
        NTT(ha,1),NTT(hb,1);
        for(rg int i=0;i<L;i++)ha[i]=(ll)ha[i]*hb[i]%mod;
        NTT(ha,-1);
        const int inv=pow(L,mod-2);
        for(rg int i=0;i<RES;i++)c.A[i]=(ll)ha[i]*inv%mod;
        return c;
    }
    inline poly inv()const
    {
        poly C;
        if(A.size()==1){C=*this;C[0]=pow(C[0],mod-2);return C;}
        C.resize((A.size()+1)>>1);
        for(rg unsigned int i=0;i<C.size();i++)C[i]=A[i];
        C=C.inv();
        L=1;while(L<(int)size()*2-1)L<<=1;
        for(rg unsigned int i=0;i<A.size();i++)ha[i]=A[i];
        for(rg unsigned int i=0;i<C.size();i++)hb[i]=C[i];
        memset(ha+A.size(),0,sizeof(int)*(L-A.size()));
        memset(hb+C.size(),0,sizeof(int)*(L-C.size()));
        NTT(ha,1),NTT(hb,1);
        for(rg int i=0;i<L;i++)ha[i]=(ll)(2-(ll)hb[i]*ha[i]%mod+mod)*hb[i]%mod;
        NTT(ha,-1);
        const int inv=pow(L,mod-2);
        C.resize(size());
        for(rg unsigned int i=0;i<size();i++)C[i]=(ll)ha[i]*inv%mod;
        return C;
    }
/*    inline poly inv()const
    {
        poly C;
        if(A.size()==1){C=*this;C[0]=pow(C[0],mod-2);return C;}
        C.resize((A.size()+1)>>1);
        for(rg unsigned int i=0;i<C.size();i++)C[i]=A[i];
        C=C.inv();
        poly D=(poly)2-*this*C;
        D.resize(size());
        D=D*C;
        D.resize(size());
        return D;
    }*///大常数版本 
    inline void Reverse(const int n)
    {
    	A.resize(n);
    	for(rg int i=0,j=n-1;i<j;i++,j--)swap(A[i],A[j]);
    }
    inline poly operator /(const poly B)const
    {
        poly a=*this,b=B;a.Reverse(size()),b.Reverse(B.size());
        b.resize(size()-B.size()+1);
        b=b.inv();
        b=b*a;
        b.Reverse(size()-B.size()+1);
        return b;
    }
    inline poly operator %(const poly B)const
    {
        poly C=(*this)-(*this)/B*B;C.resize(B.size()-1);
        return C;
    }
    inline poly diff()const
    {
        poly C;C.resize(size()-1);
        for(rg unsigned int i=1;i<size();i++)C[i-1]=(ll)A[i]*i%mod;
        return C;
    }
    inline poly inte()const
    {
        poly C;C.resize(size()+1);
        for(rg unsigned int i=0;i<size();i++)C[i+1]=(ll)A[i]*Inv[i+1]%mod;
        C[0]=0;
        return C;
    }
    inline poly ln()const
    {
        poly C=(diff()*inv()).inte();C.resize(size());
        return C;
    }
    inline poly exp()const
    {
        poly C;
        if(size()==1){C=*this;C[0]=1;return C;}
        C.resize((size()+1)>>1);
        for(rg unsigned int i=0;i<C.size();i++)C[i]=A[i];
        C=C.exp();C.resize(size());
        poly D=(poly)1-C.ln()+*this;
        D=D*C;
        D.resize(size());
        return D;
    }
    inline poly sqrt()const
    {
        poly C;
        if(size()==1)
        {
            C=*this;C[0]=Quadratic_residue(C[0]);
            return C;
        }
        C.resize((size()+1)>>1);
        for(rg unsigned int i=0;i<C.size();i++)C[i]=A[i];
        C=C.sqrt();C.resize(size());
        C=(C+*this*C.inv())*(int)499122177;
        C.resize(size());
        return C;
    }
    inline poly operator >>(const unsigned int x)const
    {
    	poly C;if(size()<x){C.resize(0);return C;}
        C.resize(size()-x);
    	for(rg unsigned int i=0;i<C.size();i++)C[i]=A[i+x];
    	return C;
    }
    inline poly operator <<(const unsigned int x)const
    {
    	poly C;C.RE(size()+x);
    	for(rg unsigned int i=0;i<size();i++)C[i+x]=A[i];
    	return C;
    }
    inline poly Pow(const unsigned int x)const
    {
    	for(rg unsigned int i=0;i<size();i++)
            if(A[i])
            {
                poly C=((((*this/A[i])>>i).ln()*x).exp()*pow(A[i],x))<<(min(i*x,size()));
                C.resize(size());
                return C;
            }
    	return *this;
    }
    inline void cheng(const poly&B)
    {
        for(rg unsigned int i=0;i<size();i++)A[i]=(ll)A[i]*B[i]%mod; 
    }
    inline void jia(const poly&B)
    {
        for(rg unsigned int i=0;i<size();i++)A[i]=Md(A[i]+B[i]); 
    }
    inline void DFT()
    {
        memset(ha,0,sizeof(int)*L);
        for(rg unsigned int i=0;i<A.size();i++)ha[i]=A[i];
        NTT(ha,1);
        resize(L);
        for(rg int i=0;i<L;i++)A[i]=ha[i];
    }
    inline void IDFT()
    {
        memset(ha,0,sizeof(int)*L);
        for(rg unsigned int i=0;i<A.size();i++)ha[i]=A[i];
        NTT(ha,-1);
        const int inv=pow(L,mod-2);
        for(rg int i=0;i<L;i++)A[i]=(ll)ha[i]*inv%mod;
        while(size()&&!A[size()-1])A.pop_back();
    }
};
}///namespace of Poly
Poly::poly a[maxn],za[maxn],b[maxn],zb[maxn],A,B;
int n,m,t;
void fz(Poly::poly*mu,Poly::poly*zi,const int l,const int r)
{
    if(l==r)return;
    const int mid=(l+r)>>1;
    fz(mu,zi,l,mid),fz(mu,zi,mid+1,r);
    const int RES=mu[l].size()<<1|1;
    Poly::L=1;while(Poly::L<RES)Poly::L<<=1;
    const int ll=l,rr=mid+1;
    zi[ll].DFT(),mu[ll].DFT(),zi[rr].DFT(),mu[rr].DFT();
    zi[ll].cheng(mu[rr]);
    zi[rr].cheng(mu[ll]);
    zi[ll].jia(zi[rr]);
    zi[ll].IDFT();
    mu[ll].cheng(mu[rr]);
    mu[ll].IDFT();
}
int fac[maxn+1],inv[maxn+1];
int main()
{
    fac[0]=1;for(rg int i=1;i<=maxn;i++)fac[i]=(ll)fac[i-1]*i%mod;
    inv[maxn]=pow(fac[maxn],mod-2);for(rg int i=maxn;i>=1;i--)inv[i-1]=(ll)inv[i]*i%mod;
    Poly::init(maxn);///namespace of Poly
    read(n),read(m);
    for(rg int i=1;i<=n;i++)
        a[i].resize(2),a[i][0]=1,read(a[i][1]),a[i][1]=Md(mod-a[i][1]),
        za[i].resize(1),za[i][0]=1;
    for(rg int i=1;i<=m;i++)
        b[i].resize(2),b[i][0]=1,read(b[i][1]),b[i][1]=Md(mod-b[i][1]),
        zb[i].resize(1),zb[i][0]=1;
    for(rg int i=1;i<=m;i++);
    read(t);
    fz(a,za,1,n);
    fz(b,zb,1,m);
    a[1].resize(100001);
    b[1].resize(100001);
    A=za[1]*a[1].inv();
    B=zb[1]*b[1].inv();
    A.resize(100001);
    B.resize(100001);
    for(rg int i=1;i<=100000;i++)A[i]=(ll)A[i]*inv[i]%mod,B[i]=(ll)B[i]*inv[i]%mod;
    A=A*B;
    const int INV=pow((ll)n*m%mod,mod-2);
    for(rg int i=1;i<=t;i++)print((ll)A[i]*fac[i]%mod*INV%mod),putchar('\n');
    return flush(),0;
}

题解2

这题有更清真的做法
求所有 $\sum a^i$ 和 $\sum b^i$ 的时候能用ln优化
考虑生成函数
$F_i(x)=\sum_{j=0}^{\infty}a_i^jx^j$
根据等比数列求和得到
$F_i(x)=\frac{1}{1-a_ix}$
我们对所有生成函数求和，设
$f(x)=\sum_{i=1}^nF_i(x)=\sum_{i=1}^n\frac{1}{1-a_ix}$
直接求常数很大
我们考虑求导
众所周知 $(ln(x))'=\frac1x,(ax)'=a$
根据复合函数求导法则
$(ln(1-a_ix))'=\frac{-a_i}{1-a_ix}$
设 $g(x)=\sum_{i=1}^n\frac{-a_i}{1-a_ix}$
我们发现 $f (x) + g (x) x = n$
移项得 $f (x) = n - g (x) x$
那么我们只要求出 $g (x)$ 就能求出 $f (x)$
$\begin{aligned} g(x)&=\sum_{i=1}^n\frac{-a_i}{1-a_ix}\\ &=\sum_{i=1}^n(ln(1-a_ix))'\\ &=(\sum_{i=1}^nln(1-a_ix))'\\ &=(ln(\prod_{i=1}^n(1-a_ix)))'\\ \end{aligned}$
这个算法总复杂度也是 $\mathcal O(nlog^2n)$
不贴代码了，我写的ln常数太大了，过不了
懒得卡常了

总结

一道经典模型好题

zhouyuheng2003

关注

0
点赞
踩
1

收藏

觉得还不错? 一键收藏
2
评论
[luoguP4705]玩游戏

前言好的思路题题目相关传送门题意简介给定长度为nnn的序列aaa和长度为mmm的序列bbb对所有k∈[1,t]k\in[1,t]k∈[1,t]，求在两个序列里分别随机取一个数的和的kkk次方的期望，模998244353998244353998244353数据范围1≤n,m≤105,1≤t≤1051\le n,m\le10^5,1\le t\le 10^51≤n,m≤105,1≤t≤...
复制链接

扫一扫

专栏目录