求组合数 II($a<=10^5,b<=10^5$，预处理的方法)

• 时间复杂度为$O(nlogn)$
• 肯定不能用组合数I的递推公式了,不然时间复杂度为$O(10^{10})$,肯定TLE
• 组合数II的转换公式为：$C_{a}^{b}=\frac{a!}{(a-b)!*b!}=fact[a]\times infact[a-b]\times infact[b]$
• 上述公式中($fact[a]$表示$a!$,$infact[b]$表示$b!$的逆元)
• 为什么可以这么转换，去搜下取模的运算法则就ok
• 这题为什么可以用费马小定理：因为mod=1e9+7是质数，只要$i!$不等于mod，两者肯定互质

我从这道题学到了什么

• 性质1：如果x是a的逆元,y是b的逆元,即$ax\equiv1,by\equiv1$，所以$abxy\equiv1$,即xy是ab的逆元，说明要求ab的逆元,直接将a的逆元和b的逆元相乘即可
• infact[i]=(LL)infact[i-1]*qmi(i,MOD-2,MOD)%MOD;这行代码就是根据性质1来的
• 证明如下：
• $(i!)^{-1}\equiv ((i-1)!)^{-1}\times (i)^{-1}\equiv ((i-1)!)^{-1}\times (i)^{p-2} $就是上面那行代码
• 其实这行代码infact[i]=(LL)infact[i-1]*qmi(i,MOD-2,MOD)%MOD;
• 我觉得这样写更好理解infact[i]=(LL)qmi(fact[i],MOD-2,MOD)%MOD;，因为$(i!)^{-1}\equiv (i!)^{p-2}$

#include <iostream>
#include <cstring>
#include <algorithm>
using namespace std;

typedef long long LL;

const int N = 100010,MOD=1e9+7;
LL fact[N],infact[N];//fact[i]代表i的阶乘,ifact[i]代表i的阶乘的逆元
int n;

int qmi(int a,int b,int p)
{
    LL res=1;
    while(b)
    {
        if(b&1) res=(LL)res*a%p;
        a=(LL)a*a%p;
        b>>=1;
    }
    return res;
}

int main()
{
    cin>>n;
    fact[0]=infact[0]=1;
    for(int i=1;i<N;i++)
    {
        fact[i]=(LL)fact[i-1]*i%MOD;//乘的时候记得都要转成LL,不然会爆
        infact[i]=(LL)infact[i-1]*qmi(i,MOD-2,MOD)%MOD;//打卡会详细说明
    }
    while (n -- )
    {
        int a,b;
        cin>>a>>b;
        cout<<(LL)fact[a]*infact[a-b]%MOD*infact[b]%MOD<<endl;//中间取一次模是为了防止爆longlong,三个数直接相乘再取模可能会爆
    }
}

求组合数 II($a<=10^5,b<=10^5$，预处理的方法)

• 时间复杂度为$O(nlogn)$
• 肯定不能用组合数I的递推公式了,不然时间复杂度为$O(10^{10})$,肯定TLE
• 组合数II的转换公式为：$C_{a}^{b}=\frac{a!}{(a-b)!*b!}=fact[a]\times infact[a-b]\times infact[b]$
• 上述公式中($fact[a]$表示$a!$,$infact[b]$表示$b!$的逆元)
• 为什么可以这么转换，去搜下取模的运算法则就ok
• 这题为什么可以用费马小定理：因为mod=1e9+7是质数，只要$i!$不等于mod，两者肯定互质

我从这道题学到了什么

• 性质1：如果x是a的逆元,y是b的逆元,即$ax\equiv1,by\equiv1$，所以$abxy\equiv1$,即xy是ab的逆元，说明要求ab的逆元,直接将a的逆元和b的逆元相乘即可
• infact[i]=(LL)infact[i-1]*qmi(i,MOD-2,MOD)%MOD;这行代码就是根据性质1来的
• 证明如下：
• $(i!)^{-1}\equiv ((i-1)!)^{-1}\times (i)^{-1}\equiv ((i-1)!)^{-1}\times (i)^{p-2} $就是上面那行代码
• 其实这行代码infact[i]=(LL)infact[i-1]*qmi(i,MOD-2,MOD)%MOD;
• 我觉得这样写更好理解infact[i]=(LL)qmi(fact[i],MOD-2,MOD)%MOD;，因为$(i!)^{-1}\equiv (i!)^{p-2}$

#include <iostream>
#include <cstring>
#include <algorithm>
using namespace std;

typedef long long LL;

const int N = 100010,MOD=1e9+7;
LL fact[N],infact[N];//fact[i]代表i的阶乘,ifact[i]代表i的阶乘的逆元
int n;

int qmi(int a,int b,int p)
{
    LL res=1;
    while(b)
    {
        if(b&1) res=(LL)res*a%p;
        a=(LL)a*a%p;
        b>>=1;
    }
    return res;
}

int main()
{
    cin>>n;
    fact[0]=infact[0]=1;
    for(int i=1;i<N;i++)
    {
        fact[i]=(LL)fact[i-1]*i%MOD;//乘的时候记得都要转成LL,不然会爆
        infact[i]=(LL)infact[i-1]*qmi(i,MOD-2,MOD)%MOD;//打卡会详细说明
    }
    while (n -- )
    {
        int a,b;
        cin>>a>>b;
        cout<<(LL)fact[a]*infact[a-b]%MOD*infact[b]%MOD<<endl;//中间取一次模是为了防止爆longlong,三个数直接相乘再取模可能会爆
    }
}

算法1：对1-n中的每一个数进行质因数分解,统计每个质因数之和,时间复杂度为$O(n\sqrt n)$,肯定会TLE

算法2:

• 先用欧拉筛,筛出1000000内的所有质数
• 对每个质数p,怎么统计1-n中p的总共次数是多少呢？
• $p$的总共次数=$p$的倍数个数+$p^2$的倍数个数+…+$p^k$的倍数个数，其中$p^k$是小于$n$中最大的$p$次幂

证明如下：

• 举个例子来说：比如n=8,p=2
• 其中2的总共次数为8次,因为1-8中有2、4、6、8里含有2,即1+2+2+3=8
• 套入公式中：2的总共次数=2的倍个数数+$2^2$的倍数个数+$2^3$的倍数个数=$n/2+n/2^2+n/2^3$=8

#include <iostream>
#include <cstring>
using namespace std;

typedef long long LL;
const int N = 1000010;
int primes[N],cnt;
bool st[N];
int n;

void get_primes(int n)
{
    for(int i=2;i<=n;i++)
    {
        if(!st[i]) primes[cnt++]=i;
        for(int j=0;primes[j]<=n/i;j++)
        {
            st[i*primes[j]]=true;
            if(i%primes[j]==0) break; 
        }
    }
}

int main () 
{
    cin>>n;
    get_primes(N-1);//初始化筛出所有质数
    for(int i=0;i<cnt;i++)
    {
        int p=primes[i];
        int s=0;
        for(int j=n;j;j/=p) s+=j/p;
        if(s) cout<<p<<" "<<s<<endl;
    }
}

算法1：对1-n中的每一个数进行质因数分解,统计每个质因数之和,时间复杂度为$O(n\sqrt n)$,肯定会TLE

算法2:

• 先用欧拉筛,筛出1000000内的所有质数
• 对每个质数p,怎么统计1-n中p的总共次数是多少呢？
• $p$的总共次数=$p$的倍数个数+$p^2$的倍数个数+…+$p^k$的倍数个数，其中$p^k$是小于$n$中最大的$p$次幂

证明如下：

• 举个例子来说：比如n=8,p=2
• 其中2的总共次数为8次,因为1-8中有2、4、6、8里含有2,即1+2+2+3=8
• 套入公式中：2的总共次数=2的倍个数数+$2^2$的倍数个数+$2^3$的倍数个数=$n/2+n/2^2+n/2^3$=8

#include <iostream>
#include <cstring>
using namespace std;

typedef long long LL;
const int N = 1000010;
int primes[N],cnt;
bool st[N];
int n;

void get_primes(int n)
{
    for(int i=2;i<=n;i++)
    {
        if(!st[i]) primes[cnt++]=i;
        for(int j=0;primes[j]<=n/i;j++)
        {
            st[i*primes[j]]=true;
            if(i%primes[j]==0) break; 
        }
    }
}

int main () 
{
    cin>>n;
    get_primes(N-1);//初始化筛出所有质数
    for(int i=0;i<cnt;i++)
    {
        int p=primes[i];
        int s=0;
        for(int j=n;j;j/=p) s+=j/p;
        if(s) cout<<p<<" "<<s<<endl;
    }
}

算法过程

• 先用欧拉筛预处理出2-1000000中的所有质数
• 然后分析一下数字性质,由于n是偶数且要求组成n的两个数为奇素数
• 那么根据偶数=奇数+奇数的性质,只要有一个是奇数,另外一个必是奇数
• 并且在所有的素数集里除了2是偶数外,其他所有质数都为奇数,所以我们不需要去判断i素数是否为奇数
• 我们只需要去判断n-primes[i]是否是素数即可，利用st[n-primes[i]]判断是否是素数

#include <iostream>
#include <cstring>
#include <algorithm>
using namespace std;

const int N = 1000010;
int primes[N];
bool st[N];
int cnt;
int n;

void get_primes(int n)
{
    for(int i=2;i<=n;i++)
    {
        if(!st[i]) primes[cnt++]=i;
        for(int j=0;primes[j]<=n/i;j++)
        {
            st[i*primes[j]]=true;
            if(i%primes[j]==0) break;
        }
    }
}

int main()
{
    ios::sync_with_stdio(0),cin.tie(0),cout.tie(0);
    get_primes(1000000);
    int n;
    while(cin>>n,n)
    {
        int flag=0;
        for(int i=0;i<cnt;i++)
        {
            if(!st[n-primes[i]])
            {
                printf("%d = %d + %d\n",n,primes[i],n-primes[i]);
                flag=1;
                break;
            }

        }
        if(!flag) cout<<"Goldbach's conjectrue is wrong."<<endl;
    }
}

算法过程

• 先用欧拉筛预处理出2-1000000中的所有质数
• 然后分析一下数字性质,由于n是偶数且要求组成n的两个数为奇素数
• 那么根据偶数=奇数+奇数的性质,只要有一个是奇数,另外一个必是奇数
• 并且在所有的素数集里除了2是偶数外,其他所有质数都为奇数,所以我们不需要去判断i素数是否为奇数
• 我们只需要去判断n-primes[i]是否是素数即可，利用st[n-primes[i]]判断是否是素数

#include <iostream>
#include <cstring>
#include <algorithm>
using namespace std;

const int N = 1000010;
int primes[N];
bool st[N];
int cnt;
int n;

void get_primes(int n)
{
    for(int i=2;i<=n;i++)
    {
        if(!st[i]) primes[cnt++]=i;
        for(int j=0;primes[j]<=n/i;j++)
        {
            st[i*primes[j]]=true;
            if(i%primes[j]==0) break;
        }
    }
}

int main()
{
    ios::sync_with_stdio(0),cin.tie(0),cout.tie(0);
    get_primes(1000000);
    int n;
    while(cin>>n,n)
    {
        int flag=0;
        for(int i=0;i<cnt;i++)
        {
            if(!st[n-primes[i]])
            {
                printf("%d = %d + %d\n",n,primes[i],n-primes[i]);
                flag=1;
                break;
            }

        }
        if(!flag) cout<<"Goldbach's conjectrue is wrong."<<endl;
    }
}

求组合数 I(数据量a<=2000,b<=2000,递推的方法)

• 性质1：$C_{a}^{0}=1$
• 性质2：$C_{a}^{b}= C_{a-1}^{b-1}+C_{a-1}^{b}$

性质2的证明(有点像01背包)

• $C_{a}^{b}$表示为在a个数里挑b个数的方案数
• 那么从第i个数去考虑,分为两种情况选第i个数和不选第i个数
• 情况1选第i个数：那么剩下b-1个数就在a-1个数里面选,即$C_{a-1}^{b-1}$
• 情况2不选第i个数：那么剩下b个数就在a-1个数里面选,即$C_{a-1}^{b}$
• 情况1和情况2的方案数相加即为$C_{a}^{b}$的方案数

#include <iostream>
#include <cstring>
#include <algorithm>
using namespace std;

const int N = 2010,mod=1e9+7;
int c[N][N];
int n;

void init()//预处理所有的组合数,时间复杂度在4*10^6
{
    for(int i=0;i<N;i++)//坑点：不能写成i<=N,不然会数组越界
        for(int j=0;j<=i;j++)
        {
            if(!j) c[i][j]=1;//这就是定义c[i][0]=1;
            else c[i][j]=(c[i-1][j-1]+c[i-1][j])%mod;//有点像01背包思想
        }
}

int main()
{
    ios::sync_with_stdio(0),cin.tie(0),cout.tie(0);
    cin>>n;
    init();
    while (n -- )
    {
        int a,b;
        cin>>a>>b;
        cout<<c[a][b]%mod<<endl;
    }

}

求组合数 I(数据量a<=2000,b<=2000,递推的方法)

• 性质1：$C_{a}^{0}=1$
• 性质2：$C_{a}^{b}= C_{a-1}^{b-1}+C_{a-1}^{b}$

性质2的证明(有点像01背包)

• $C_{a}^{b}$表示为在a个数里挑b个数的方案数
• 那么从第i个数去考虑,分为两种情况选第i个数和不选第i个数
• 情况1选第i个数：那么剩下b-1个数就在a-1个数里面选,即$C_{a-1}^{b-1}$
• 情况2不选第i个数：那么剩下b个数就在a-1个数里面选,即$C_{a-1}^{b}$
• 情况1和情况2的方案数相加即为$C_{a}^{b}$的方案数

#include <iostream>
#include <cstring>
#include <algorithm>
using namespace std;

const int N = 2010,mod=1e9+7;
int c[N][N];
int n;

void init()//预处理所有的组合数,时间复杂度在4*10^6
{
    for(int i=0;i<N;i++)//坑点：不能写成i<=N,不然会数组越界
        for(int j=0;j<=i;j++)
        {
            if(!j) c[i][j]=1;//这就是定义c[i][0]=1;
            else c[i][j]=(c[i-1][j-1]+c[i-1][j])%mod;//有点像01背包思想
        }
}

int main()
{
    ios::sync_with_stdio(0),cin.tie(0),cout.tie(0);
    cin>>n;
    init();
    while (n -- )
    {
        int a,b;
        cin>>a>>b;
        cout<<c[a][b]%mod<<endl;
    }

}

推导过程

• $ax \equiv b(mod\;m)$就可以推出$ax =my+b$
• 变形一下$ax-my =b$,不就是系数为a和-m的扩欧算法吗,即exgcd(a,-m,x,y)
• 怎么判断有没有解呢,根据裴蜀定理ax+my所能构造的最小正整数就是gcd(a,m),即$ax+my=gcd(a,m)$
• 那么要有解的情况下,b必须是gcd(a,m)的倍数,即$b\%gcd(a,m)=0$
• 由于exgcd(a,m,x,y)求得的x只是ax+my=gcd(a,m)的解
• 要求得ax+my=b的解,直接将x乘以b/gcd(a,m)的倍数即可,即$x^{\prime}=x\times \frac{b}{gcd(a,m)}$
• 最后一定要转换为longlong,防止爆int,并且要转换到int范围里，取余m即可

#include <iostream>
#include <cstring>
#include <algorithm>
using namespace std;

typedef long long LL;

int n;

int exgcd(int a,int b,int &x,int &y)
{
    if(!b)
    {
        x=1,y=0;
        return a;
    }
    int g=exgcd(b,a%b,x,y);
    int t=x;
    x=y;
    y=t-a/b*y;
    return g;
}

int main()
{
    cin>>n;
    while (n -- )
    {
        int a,b,m,x,y;
        cin>>a>>b>>m;
        int d=exgcd(a,-m,x,y);//这里求出来的x只是ax+my=gcd(a,m)的x,要求得ax+my=b;方程两边都扩大b/gcd(a,m)倍
        if(b%d!=0) cout<<"impossible"<<endl;
        else cout<<(LL)x*(b/d)%m<<endl;//放在m的解系中
    }

}

推导过程

• $ax \equiv b(mod\;m)$就可以推出$ax =my+b$
• 变形一下$ax-my =b$,不就是系数为a和-m的扩欧算法吗,即exgcd(a,-m,x,y)
• 怎么判断有没有解呢,根据裴蜀定理ax+my所能构造的最小正整数就是gcd(a,m),即$ax+my=gcd(a,m)$
• 那么要有解的情况下,b必须是gcd(a,m)的倍数,即$b\%gcd(a,m)=0$
• 由于exgcd(a,m,x,y)求得的x只是ax+my=gcd(a,m)的解
• 要求得ax+my=b的解,直接将x乘以b/gcd(a,m)的倍数即可,即$x^{\prime}=x\times \frac{b}{gcd(a,m)}$
• 最后一定要转换为longlong,防止爆int,并且要转换到int范围里，取余m即可

#include <iostream>
#include <cstring>
#include <algorithm>
using namespace std;

typedef long long LL;

int n;

int exgcd(int a,int b,int &x,int &y)
{
    if(!b)
    {
        x=1,y=0;
        return a;
    }
    int g=exgcd(b,a%b,x,y);
    int t=x;
    x=y;
    y=t-a/b*y;
    return g;
}

int main()
{
    cin>>n;
    while (n -- )
    {
        int a,b,m,x,y;
        cin>>a>>b>>m;
        int d=exgcd(a,-m,x,y);//这里求出来的x只是ax+my=gcd(a,m)的x,要求得ax+my=b;方程两边都扩大b/gcd(a,m)倍
        if(b%d!=0) cout<<"impossible"<<endl;
        else cout<<(LL)x*(b/d)%m<<endl;//放在m的解系中
    }

}