洛谷 P3726 [AH2017/HNOI2017]抛硬币 exgcd+扩展Lucas定理

最新推荐文章于 2021-10-07 20:30:50 发布

Amber_lylovely

最新推荐文章于 2021-10-07 20:30:50 发布

阅读量656

点赞数

CC 4.0 BY-SA版权

分类专栏：数论

本文链接：https://blog.youkuaiyun.com/liangzihao1/article/details/87867959

数论专栏收录该内容

44 篇文章

订阅专栏

本文探讨了一种基于游戏概率的算法挑战，通过分析硬币抛掷游戏中的胜利条件，介绍了如何计算玩家A在不同抛掷次数下战胜玩家B的概率。算法涉及组合数学、概率论以及大数处理技巧，包括模数运算、中国剩余定理的应用，以及如何高效求解大数阶乘和组合数。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

题目描述
小 A 和小 B 是一对好朋友，他们经常一起愉快的玩耍。最近小 B 沉迷于**师手游，天天刷本，根本无心搞学习。但是已经入坑了几个月，却一次都没有抽到 SSR，让他非常怀疑人生。勤勉的小 A 为了劝说小 B 早日脱坑，认真学习，决定以抛硬币的形式让小 B 明白他是一个彻彻底底的非洲人，从而对这个游戏绝望。两个人同时抛 b 次硬币，如果小 A 的正面朝上的次数大于小 B 正面朝上的次数，则小 A 获胜。

但事实上，小 A 也曾经沉迷过拉拉游戏，而且他一次 UR 也没有抽到过，所以他对于自己的运气也没有太大把握。所以他决定在小 B 没注意的时候作弊，悄悄地多抛几次硬币，当然，为了不让小 B 怀疑，他不会抛太多次。现在小 A 想问你，在多少种可能的情况下，他能够胜过小 B 呢？由于答案可能太大，所以你只需要输出答案在十进制表示下的最后 k 位即可。

输入输出格式
输入格式：
有多组数据，对于每组数据输入三个数a,b,k，分别代表小A抛硬币的次数，小B抛硬币的次数，以及最终答案保留多少位整数。

输出格式：
对于每组数据，输出一个数，表示最终答案的最后 k 位为多少，若不足 k 位以 0 补全。

输入输出样例
输入样例#1：
2 1 9
3 2 1
输出样例#1：
000000004
6

说明
对于第一组数据，当小A抛2次硬币，小B抛1次硬币时，共有4种方案使得小A正面朝上的次数比小B多。
(01,0), (10,0), (11,0), (11,1)
对于第二组数据，当小A抛3次硬币，小B抛2次硬币时，共有16种方案使得小A正面朝上的次数比小B多。
(001,00), (010,00), (100,00), (011,00), (101,00), (110,00), (111,00), (011,01), (101,01), (110,01),(111,01), (011,10), (101,10), (110,10), (111,10), (111,11).

数据范围
10%的数据满足a,b≤20；
30%的数据满足a,b≤100；
70%的数据满足a,b≤100000，其中有20%的数据满足a=b；
100%的数据满足 $1≤a,b≤1015,b≤a≤b+10000,1≤k≤91\le a,b\le 10^{15},b\le a\le b+10000,1\le k\le 9$ ，数据组数小于等于10。

分析：
假设一种方案为 $s$ ，另一种方案 $s^{'}$ 为 $s$ 的每一位取反。
假设 $a = b$ ，那么 $s$ 为胜利状态，那么 $s^{'}$ 必然为失败状态，只要减去平局状态再除以2即可。即 $2a+b−(a+ba)2\frac{2^{a+b}-\binom{a+b}{a}}{2}$ 。
对于 $a > b$ ，假设 $s$ 为平局或失败状态， $s^{'}$ 必然为胜利状态。还有一些胜利状态 $s$ 对应胜利状态 $s^{'}$ 。
因为 $s$ 和 $s^{'}$ 一一对应，显然 $s^{'}$ 所代表的胜利状态恰好为 $2^{a+b-1}$ 种，剩下只要加上 $s$ 中的胜利状态个数。

现在这一类中 $1$ 的个数满足，
$n u m a > n u m b$
$a - n u m a > b - n u m b$
所以
$a - b > n u m a - n u m b$
枚举其中一个的 $1$ 的个数，那么这一个部分的答案为
$ans=∑i=0b∑j=1a−b−1(bi)∗(ai+j)=∑i=0b∑j=1a−b−1(bb−i)∗(ai+j)ans=\sum_{i=0}^{b}\sum_{j=1}^{a-b-1}\binom{b}{i}*\binom{a}{i+j}=\sum_{i=0}^{b}\sum_{j=1}^{a-b-1}\binom{b}{b-i}*\binom{a}{i+j}$
后面相当于直接从 $a + b$ 中选 $j$ 个，枚举 $i$ 只是为了判断一个子集选多少个。
也就是求
$∑i=1a−b−1(a+bb+j)2\frac{\sum_{i=1}^{a-b-1}\binom{a+b}{b+j}}{2}$
总的方案数是 $ans=2a+b+∑i=1a−b−1(a+bb+j)2ans=\frac{2^{a+b}+\sum_{i=1}^{a-b-1}\binom{a+b}{b+j}}{2}$

我们可以直接把模数设为 $10^9$ ，后面只需要保留相应位数即可。
因为模数不是质数，我们可以先把模数因式分解成 $mod=p_1^{a_1}*p_2^{a_2}*...*p_m^{a_m}$ 的形式，对于每一个 $p_i^{a_i}$ 次方分别求解，然后再中国剩余定理合并答案。

考虑怎样求 $pk\binom{n}{m}\ mod\ p^k$ ，其中 $p$ 为质数。
我们把 $n!$ 拆成 $a*p^b$ 的形式，然后再计算。
其中 $b$ 为 $1$ 到 $n$ 所有数包含 $p$ 这个素因子的个数的和，即 $b=np+np2+np3+....b=\frac{n}{p}+\frac{n}{p^2}+\frac{n}{p^ 3}+....$
显然 $n! = 1 * 2 * 3 * . . . . * n$
我们把是 $p$ 的倍数的放在一边，就是
$n!=(∏i=1p−1i)∗(∏i=p+12∗p−1i)∗...∗p∗(2∗p)∗(3∗p)∗....n!=(\prod_{i=1}^{p-1} i)*(\prod_{i=p+1}^{2*p-1} i)*...*p*(2*p)*(3*p)*....$
$=(∏i=1p−1i)∗(∏i=p+12∗p−1i)∗...∗pd∗(n/p)!=(\prod_{i=1}^{p-1} i)*(\prod_{i=p+1}^{2*p-1} i)*...*p^d*(n/p)!$

显然 $pk)\prod_{i=1}^{p-1} i=\prod_{i=p^k+1}^{p^k+p-1} i(mod\ p^k)$
也就是存在循环节，我们需要把一个循环节算 $n/(p^k)$ 次方，再把多出来的乘上就可以了， $p^d$ 已经被统计进 $b$ 中了。
然后阶乘要预处理，这里的阶乘是除掉 $p$ 的倍数的，后面就递归下去即可。
然后再计算组合数时， $p$ 的指数直接减掉，而 $a$ 的部分可以求逆元，显然 $a$ 与 $p^k$ 互质，考虑exgcd。
这样就算出组合数了。
最后要把结果除以 $2$ ，这些组合数可以发现是两两对称的，我们只需要统计一半。如果中间多了一个，这个一定是 $(a+b(a+b)/2)\binom{a+b}{(a+b)/2}$ ，这个数一定是 $2$ 的倍数。

代码：

// luogu-judger-enable-o2
#include <iostream>
#include <cstdio>
#include <cmath>
#define LL long long

const int maxn=2e6+7;
const LL mod[2]={512,1953125};
const LL M=1e9;

using namespace std;

LL n,m,k,ans;
LL jc[maxn][2];

void prework() //预处理阶乘，注意不含p的倍数
{
    jc[0][0]=1;
    for (LL i=1;i<mod[0];i++)
    {
        if (i%2) jc[i][0]=jc[i-1][0]*(LL)i%mod[0];
            else jc[i][0]=jc[i-1][0];
    }
    jc[0][1]=1;
    for (LL i=1;i<mod[1];i++)
    {
        if (i%5) jc[i][1]=jc[i-1][1]*(LL)i%mod[1];
            else jc[i][1]=jc[i-1][1];
    }
}

LL ksm(LL x,LL y,LL k)
{
    if (y==0) return 1;
    LL c=ksm(x,y/2,k);
    c=(c*c)%k;
    if (y&1) c=(c*x)%k;
    return c;
}

LL exgcd(LL a,LL b,LL &x,LL &y)
{
    if (b==0)
    {
        x=1;
        y=0;
        return a;
    }
    LL d=exgcd(b,a%b,x,y);
    LL z=x;
    x=y,y=z-a/b*y;
    return d;
}

LL inv(LL a,LL mod)
{
    LL x,y;
    exgcd(a,mod,x,y);
    return (x%mod+mod)%mod;
}

LL multi(LL n,LL pi,LL pk) //求出不含p的部分的积
{
    if (!n) return 1; //递归结束条件
    LL c=ksm(jc[pk-1][pi!=2],n/pk,pk)*jc[n%pk][pi!=2]%pk; //循环节有n/pk个，每个都是满的，直接乘方后再乘上多出来的n%pk个
    return c*multi(n/pi,pi,pk)%pk; //递归处理
}

LL binom(LL n,LL m,LL pi,LL pk,bool flag)
{
    LL p=0;
    for (LL i=pi;i<=n;i*=pi) p+=n/i;
    for (LL i=pi;i<=m;i*=pi) p-=m/i;
    for (LL i=pi;i<=n-m;i*=pi) p-=(n-m)/i; //算出p的幂次
    if ((flag) && (pi==2)) p--; //要对组合数除以2的情况
    if (p>=9) return 0;
    LL a=multi(n,pi,pk),b=multi(m,pi,pk),c=multi(n-m,pi,pk);
    LL cs=a*inv(b,pk)%pk*inv(c,pk)%pk*ksm(pi,p,pk)%pk;
    if ((flag) && (pi==5)) cs=cs*inv(2,pk)%pk; 
    return cs;
}

LL getans(LL n,LL m,bool flag) //中国剩余定理合并，因为打了exgcd，可以利用上
{
    LL a=binom(n,m,2,mod[0],flag); 
    LL b=binom(n,m,5,mod[1],flag);
    LL x,y;
    exgcd(mod[0],mod[1],x,y);
    x=(x*(b-a)%M*mod[0]%M+a+M)%M;	
    return x;
}

void prin(LL x)
{
    int b[20],num=0;
    for (int i=1;i<=k;i++) b[i]=0;
    while (x)
    {
        b[++num]=x%10;
        x/=10;
    }
    for (int i=k;i>0;i--) printf("%d",b[i]);
    printf("\n");
}

int main()
{
    prework();
    while (scanf("%lld%lld%lld",&n,&m,&k)!=EOF)
    {		
        if (n==m) ans=(ksm(2,n+m-1,M)+M-getans(n+m,n,1))%M;
        else
        {
            ans=ksm(2,n+m-1,M);			
            for (LL i=(n+m)/2+1;i<=n-1;i++) ans=(ans+getans(n+m,i,0))%M;
            if ((n+m)%2==0) ans=(ans+getans(n+m,(n+m)/2,1))%M;
        }
        prin(ans);
    }
}