BZOJ 4403:浅谈Lucas定理应用及组合数建模

来源:互联网 发布:数据挖掘的技术基础 编辑:程序博客网 时间:2024/06/05 00:26

这里写图片描述
世界真的很大
在做数论题时,组合数是一个不得不讨论的问题,因为这可能是数论基础中的基础,组合数的计算公式不必多说,但是很多题看似和组合数有很大关系但是硬生生没法建立组合数的模型,套用公式。很多题看起来很难但如果对于组合数建立模型成功的话,几乎马上就会变成直接计算组合数的水题。所以以这道题为例,浅谈一下数论题中组合数的建模思路。
先来描述一下题意:

给定三个正整数N、L和R,统计长度在1N之间,元素大小都在L到R之间的单调不降序列的数量。输出答案对10^6+3取模的结果。

input:

输入第一行包含一个整数T,表示数据组数。第2到第T+1行每行包含三个整数N、L和R,N、L和R的意义如题所述。1N,L,R≤10^91≤T≤100,输入数据保证L≤R。

output:

输出包含T行,每行有一个数字,表示你所求出的答案对10^6+3取模的结果。

Sample input:

21 4 52 4 5

Sample output:

25//【样例说明】满足条件的2个序列为[4]和[5]。

乍一看的确是组合数,从区间里任意选出数字,排个序,就能构成单调不降序列,无关乎选择顺序,就是求l~r区间这么多个数里面,选出i(1<=i<=n)个数的组合数,如果是这样的建模就没有必要浅谈了,,
并没有这样简单,选数其实是可以重复的,题目要求的是单调不降序列,如:4,4,4,4,也是单调不降的,这里就有问题了,组合数没法处理重复的情况,组合数算出来的数量里不包含4,4之类的组合。
我们尝试模板化建模过程:
我们想要得到的,如样例,4-5区间里的单调不降序列的个数,即选i长度的序列的个数,即从l-r选i个数的包含重复的组合数,比如4-5中选长度为2的,有【4,4】,【4,5】,【5,5】三种组合数C(2,2)算出的仅有【4,5】,无法算出重复的情况。所以我们考虑将4,4分开,把4,5分成4,4,5,5,计算C(4,2),但是又有新的问题,【4,5】区间被计算了两次,这是由于组合数计算时把两个4当成了独立的个体,并不知道这两个4其实是相同的。考虑投影元素(如4)失败了,于是乎想到,因为最后我们需要的是【4,4】,【4,5】,【5,5】三种组合,想办法将组合分开,给予每个数新的“映射”,将组合投影,不为是一种选择。
考虑【4,4】,保持序列第一个不变,将第二个数+1,得到(4,5),将【4,5】同样处理得到(4,6),处理【5,5】得到(5,6),原序列的各种合法组合通过映射的方式得到了新的组合,而由新的组合,反映出了新的区间(4,6),其中的每一种新的组合对应了原来的一种合法组合,且避免了【4,4】一类组合数无法处理的情况。在枚举长度为i时,因为依次+,最大得数r就加了i-1,那区间长度就变成了i+m-1,(m=r-l+1),组合数就变成了C(i,m+i-1),等价于C(m-1,m+i-1),接下来就枚举i从1到n求和,就是答案。
设M=R−L+1
长度为i,元素大小在1…M之间的单调不降序列的数量有CM−1i+M−1个
故答案为
∑Ni=1CM−1i+M−1
=(∑Ni=1CM−1i+M−1)+CMM−1
=(∑Ni=2CM−1i+M−1)+CMM+1−1
=(∑Ni=3CM−1i+M−1)+CMM+2−1

=CMN+M−1
这里用到了组合数的一个公式:

C(M-1N-1)+C(M-1N)=C(M,N

相当于我们直接计算C(m+n,n)-1就行了,即我们成功将其模板化了(注意理解模板化的思路步骤)
但是数据范围是1e9,无法直接计算那么大的组合数,所以再次用数论知识来优化它,即Lucas定理:

C(a,b)=C(a/mod,b/mod)*C(a%mod,b%mod) (mod mod)

mod是一个质数这道题里就是1e6+3.
首先预处理出阶乘saber数组:

void init(int n){    register int i;     saber[0]=1;    for(i=1;i<=n;i++)        saber[i]=(dnt) i*saber[i-1]%mod;}

然后使用递归的方式调用Lucas定理,因为阶乘一旦超过mod后,就全是0了,因为mod mod嘛,所以用Lucas定理将范围缩小到mod以内。代码:

dnt Misaka(int a,int b){    if(a<b) return 0;    if(a<mod) return saber[a]*Railgun(saber[b]*saber[a-b]%mod)%mod;    return (dnt) Misaka(a/mod,b/mod)*Misaka(a%mod,b%mod)%mod;}

这里的Railgun函数是求逆元的过程,因为saber数组都是对mod取过mod的,用除法处理会爆精度,所以把除法改成乘以他的逆元,即在mod mod意义下saber[b]*saber[a-b]的倒数。
逆元我们在线求,用扩展欧几里得的方法,代码:

int exgcd( int a, int b, int &x, int &y ) {    if(b==0)     {        x=1;        y=0;        return a;    }       int x0,y0;    int cd=exgcd(b,a%b,x0,y0);    x=y0;    y=x0-(a/b)*y0;    return cd;}dnt Railgun( int a ) {    int x, y;    exgcd(a,mod,x,y);    return (dnt) (x%mod+mod)%mod;}

嗯,就是这样,完整代码:

#include<stdio.h>#include<algorithm>using namespace std;typedef long long dnt;int n,m,mod=1e6+3,tt,l,r;dnt saber[1000010];void init(int n){    register int i;     saber[0]=1;    for(i=1;i<=n;i++)        saber[i]=(dnt) i*saber[i-1]%mod;}int exgcd( int a, int b, int &x, int &y ) {    if(b==0)     {        x=1;        y=0;        return a;    }       int x0,y0;    int cd=exgcd(b,a%b,x0,y0);    x=y0;    y=x0-(a/b)*y0;    return cd;}dnt Railgun( int a ) {    int x, y;    exgcd(a,mod,x,y);    return (dnt) (x%mod+mod)%mod;}dnt Misaka(int a,int b){    if(a<b) return 0;    if(a<mod) return saber[a]*Railgun(saber[b]*saber[a-b]%mod)%mod;    return (dnt) Misaka(a/mod,b/mod)*Misaka(a%mod,b%mod)%mod;}int main(){    init(1e6+10);    scanf("%d",&tt);    while(tt--)    {        scanf("%d%d%d",&n,&l,&r);        m=r-l+1;        printf("%lld\n",(Misaka(m+n,n)%mod-1+mod)%mod);    }    return 0;}