[NOI2019]斗主地

斗主地

题解

我们可以先来考虑位置 i i i的牌一轮中被洗到位置 j j j的概率。
如果 i ⩽ A i\leqslant A iA,贡献显然是 ( n − j ) ! A ! ( n − A ) ! n ! ( A − i ) ! ( n − A − j + i ) ! ( j − 1 i − 1 ) = ( j − 1 i − 1 ) ( n − j A − i ) ( n A ) \frac{(n-j)!A!(n-A)!}{n!(A-i)!(n-A-j+i)!}\binom{j-1}{i-1}=\frac{\binom{j-1}{i-1}\binom{n-j}{A-i}}{\binom{n}{A}} n!(Ai)!(nAj+i)!(nj)!A!(nA)!(i1j1)=(An)(i1j1)(Ainj),我们可以把每个 i i i j j j的贡献的算出来,这相当于一个矩阵。
初始值相当于是一个向量,它乘上这 m m m个矩阵就相当于进行了 m m m次洗牌。
而所有 A A A都相等的部分相当于是直接矩阵快速幂,这样我们就能拿到 40 p t s 40pts 40pts了。

但让我们理性观察一下后面的点,它们的 n n n都达到了 1 0 7 10^7 107的级别,显然不可能用我们上面矩阵的方法做了,理性考虑一下还有什么可能做的方法。
我们观察一下初始的这些 f ( i ) = i f(i)=i f(i)=i f ( i ) = i 2 f(i)=i^2 f(i)=i2,它们好像都是多项式的形式,理性猜测一下, 我们变化后每个位置的 f ′ ( i ) f'(i) f(i)依旧是一个多项式。
然后你随便打个插值什么的,可以发现, t y p = 1 typ=1 typ=1就是一个一次多项式, t y p = 2 typ=2 typ=2是一个二次多项式。
所以我们可以暴力手算出下一轮三个位置上的值,然后把下一轮的多项式给插出来。
好的,做法这里就讲完了,下面我们来讲讲它为什么是一个多项式。

考虑归纳证明,显然,最开始它肯定是多项式,考虑多项式变化一轮后的样子。
首先,多项式肯定是可以通过斯特林反演变化成组合数的,我们假定我们上一层的多项式为 f ( x ) = ∑ k a k ( x − 1 k ) f(x)=\sum_{k}a_k\binom{x-1}{k} f(x)=kak(kx1)。(这里化成 ( x − 1 k ) \binom{x-1}{k} (kx1)而不是 ( x k ) \binom{x}{k} (kx)是为了后面好做点)
考虑下一层位置 i i i上的值怎么计算,利用我们上面推出来的贡献,这里照样是只算 ⩽ A \leqslant A A部分的贡献值,因为大于 A A A的与小于 A A A的部分形式是基本一样的,如果小的部分是多项式,大的部分也会是多项式。
f i ′ = ∑ k a k ( n A ) ∑ j = 1 A ( j − 1 k ) ( i − 1 j − 1 ) ( n − i A − j ) = ∑ k a k ( n A ) ∑ j = 1 A ( i − 1 k ) ( i − k − 1 j − k − 1 ) ( n − i A − j ) = ∑ k a k ( i − 1 k ) ( n A ) ∑ j = 1 A − k − 1 ( i − k − 1 j ) ( n − i A − j − k − 1 ) = ∑ a k ( n − k − 1 A − k − 1 ) ( n A ) ( i − 1 k ) f'_i=\sum_{k}\frac{a_k}{\binom{n}{A}}\sum_{j=1}^A\binom{j-1}{k}\binom{i-1}{j-1}\binom{n-i}{A-j}\\ =\sum_k\frac{a_k}{\binom{n}{A}}\sum_{j=1}^A\binom{i-1}{k}\binom{i-k-1}{j-k-1}\binom{n-i}{A-j}\\ =\sum_{k}\frac{a_k\binom{i-1}{k}}{\binom{n}{A}}\sum_{j=1}^{A-k-1}\binom{i-k-1}{j}\binom{n-i}{A-j-k-1}\\ =\sum\frac{a_k\binom{n-k-1}{A-k-1}}{\binom{n}{A}}\binom{i-1}{k} fi=k(An)akj=1A(kj1)(j1i1)(Ajni)=k(An)akj=1A(ki1)(jk1ik1)(Ajni)=k(An)ak(ki1)j=1Ak1(jik1)(Ajk1ni)=(An)ak(Ak1nk1)(ki1)这不又是一个多项式吗?
根据这个式子我们可以发现,上一个多项式里每一项系数乘上一个只与 n n n A A A相关的常数就会变成下一个多项式。当然,我们这里的多项式是组合数形式下的
组合数当然可以重新拆成次数项的多项式,所以我们也就能还原出来当前的多项式。并且,多项式的次数也没有增加。
我们就证明了我们每一次洗牌操作后每个位置的期望权值仍然是一个多项式。
那么我们之后可以做的方法就多了,可以像上面说的那样插值,也可以直接维护组合数形式下的多项式,不过这样要预处理出来 n n n以下的阶乘。

时间复杂度 O ( m log ⁡ P   or   m + n ) O\left(m\log P\, \text{or} \,m+n\right) O(mlogPorm+n)
当然,由于多项式只有常数项,当然不会写进时间复杂度里,不过两种实现常数肯定会有一定区别。
另外,我觉得正常人类应该不大可能去猜这东西是一个多项式

源码

#include<bits/stdc++.h>
using namespace std;
typedef long long LL;
typedef unsigned long long uLL;
typedef pair<LL,int> pii;
#define MAXN 500005
#define pb push_back
#define mkpr make_pair
#define fir first
#define sec second
#define lowbit(x) (x&-x)
const int mo=998244353;
const int inv2=5e8+4;
const int jzm=2333;
const int zero=15;
const LL INF=0x3f3f3f3f3f3f3f3f;
const double Pi=acos(-1.0);
const double eps=1e-9;
const int lim=1000000;
const int orG=3,ivG=332748118;
template<typename _T>
_T Fabs(_T x){return x<0?-x:x;}
template<typename _T>
void read(_T &x){
    _T f=1;x=0;char s=getchar();
    while(s<'0'||s>'9'){if(s=='-')f=-1;s=getchar();}
    while('0'<=s&&s<='9'){x=(x<<3)+(x<<1)+(s^48);s=getchar();}
    x*=f;
}
int add(int x,int y,int p){return x+y<p?x+y:x+y-p;}
void Add(int &x,int y,int p){x=add(x,y,p);}
int qkpow(int a,int s,int p){int t=1;while(s){if(s&1)t=1ll*a*t%p;a=1ll*a*a%p;s>>=1;}return t;}
int N,m,typ,Q;
struct poly{
    int a[3],n;poly(){}
    void clear(){a[0]=a[1]=a[2]=0;}
    poly operator + (const poly &rhs){
        poly res;res.clear();res.n=max(n,rhs.n);
        for(int i=0;i<=n;i++)res.a[i]=add(a[i],rhs.a[i],mo);
        return res;
    }
    poly operator * (const int &rhs){
        poly res;res.clear();res.n=n;
        for(int i=0;i<=n;i++)res.a[i]=1ll*rhs*a[i]%mo;
        return res;
    }
    poly operator ^ (const int &rhs){
        poly res;res.clear();res.n=n+1;
        for(int i=0;i<=n;i++)res.a[i]=mo-1ll*rhs*a[i]%mo;
        for(int i=0;i<=n;i++)Add(res.a[i+1],a[i],mo);
        return res;
    }
    int ask(const int &x){
        int res=0,now=1;if(x<1||x>N)return 0;
        for(int i=0;i<=n;i++,now=1ll*now*x%mo)
            Add(res,1ll*now*a[i]%mo,mo);
        return res;
    }
}P,B;
int main(){
    read(N);read(m);read(typ);P.clear();
    int invn=qkpow(N,mo-2,mo),inv1=qkpow(N-1,mo-2,mo);
    if(typ==1)P.n=P.a[1]=1;else P.n=2,P.a[2]=1;
    for(int i=1;i<=m;i++){
        int A;read(A);if(A==N)continue;
        if(typ==1){
            int t1=1ll*(1ll*A*P.ask(1)+1ll*(N-A)*P.ask(A+1))%mo*invn%mo;
            int t2=1ll*(1ll*A*P.ask(A)+1ll*(N-A)*P.ask(N))%mo*invn%mo;
            P.clear();P.n=0;P.a[0]=t1;
            B.clear();B.n=1;B.a[1]=1;B.a[0]=mo-1;
            int tmp=1ll*add(t2,mo-P.ask(N),mo)*qkpow(B.ask(N),mo-2,mo)%mo;P=B*tmp+P;
        }
        else{
            int t1=1ll*(1ll*A*P.ask(1)+1ll*(N-A)*P.ask(A+1))%mo*invn%mo;
            int t2=1ll*(1ll*A*P.ask(A)+1ll*(N-A)*P.ask(N))%mo*invn%mo;
            int t3=(1ll*(1ll*(A-1)*P.ask(2)+1ll*(N-A)*P.ask(A+1))%mo*invn%mo*inv1%mo*A+
                1ll*(1ll*A*P.ask(1)+1ll*(N-A-1)*P.ask(A+2))%mo*invn%mo*inv1%mo*(N-A))%mo;
            P.clear();P.n=0;P.a[0]=t3;
            B.clear();B.n=1;B.a[1]=1;B.a[0]=mo-2;
            int tmp1=1ll*add(t1,mo-P.ask(1),mo)*qkpow(B.ask(1),mo-2,mo)%mo;P=B*tmp1+P;B=B^1;
            int tmp2=1ll*add(t2,mo-P.ask(N),mo)*qkpow(B.ask(N),mo-2,mo)%mo;P=B*tmp2+P;
        }
    }
    read(Q);
    for(int i=1,x;i<=Q;i++)
        read(x),printf("%d\n",P.ask(x));
    return 0;
}

谢谢

  • 2
    点赞
  • 1
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值