[自用] 数论和组合计数类数学相关（定理&证明&板子）

0 写在前面

　　本文受 NaVi_Awson 的启发，有些地方相似，一些地方甚至直接引用，特此说明(感谢dalao)。

1 数论

　　1.0 gcd

　　　　1.0.0 gcd

　　　　　　$gcd(a,b) = gcd(b,a\;mod\;b)$

　　　　　证明：设 $c\mid a$，$c\mid b$，则 $c\mid (b-a)$。

　　　　　　　　设 $c\nmid a$，则 $c$ 不是 $a,b-a$ 的公因子。

　　　　　　　　设 $c\mid a$，$c\nmid b$，则 $c$ 不是 $a,b-a$ 的公因子。

 

1 int gcd(int a,int b){
2     if(!b) return a;
3     return gcd(b,a%b);
4 }

 　　　　1.0.1 exgcd

　　　　　　　对于一组不定方程 $a\times b+b\times y = c$。若 $gcd(a,b)\mid c$ 则这组方程有无数组解，否则无解。

　　　　　　　求出一组 $a\times b+b\times y = gcd(a,b) $ 的解需要 $exgcd$，是在 $gcd$ 递归函数过程中实现的。

　　　　　　　设 $x_0,y_0$ 是这一层的解，$x_1,y_1$ 是上一层的解。

　　　　　　　首先边界条件：当$b=0$时显然$x_0=1,y_0=0$。

　　　　　　　其他情况下有：$x_0=y_1,y_0=x_1- {\lfloor {\frac{a}{b}} \rfloor} \times y_1$。

　　　　   证明：

　　　　　　　　　　　　　　　　　　　　　　　$$gcd(a,b)=x_0 \times a + y_0 \times b$$

　　　　　　　　　　　$$gcd(b,a\;mod\;b)=x_1 \times b + y_1 \times {(a\;mod\;b)} = x_1 \times b + y_1 \times {(a- {\lfloor {\frac{a}{b}} \rfloor)} \times b}$$

　　　　　　　因为 　　　　　　　　　　　　　　　　$$gcd(a,b)=gcd(b,a\;mod\;b)$$

　　　　　　　得 　　　　　　　　　　　　　　　　　$$x_0 \times a + y_0 \times b =x_1 \times b + y_1 \times {(a-{\lfloor {\frac{a}{b}} \rfloor)} \times b}$$

　　　　　　　合并得　　　　　　 　　　　　　　　　$$x_0 \times a+ y_0 \times b=y_1 \times a+ (x_1-{\lfloor {\frac{a}{b}} \rfloor} \times y_1) \times b$$

　　　　　　　根据恒等定理 　　　　　　　　　　　　$$x_0=y_1,y_0=x_1-{\lfloor {\frac{a}{b}} \rfloor} \times y_1$$

 1 int exgcd(int a, int b, int &x, int &y) {
 2     if (!b) {
 3         x=1; y=0;
 4         return a;
 5     }
 6     int c=exgcd(b,a%b,x,y);
 7     int t=x;
 8     x=y;
 9     y=t-a/b*y;
10     return c;
11 }

　　　　upd：对于方程 $a\times x+b\times y=c$ 解出两个特解 $x_0,y_0$ 之后，方程的通解是 $x=x_0+k\frac{b}{gcd(a,b)},y=y_0+k\frac{a}{gcd(a,b)}$。

　  证明：对于两组解 $x_0,y_0$ 和 $x_1,y_1$。显然满足 $ax_0+by_0=c,ax_1+by_1=c$。我们想要求通解，实际上就是想求一个 $\delta x$ 和 $\delta y$，这样通解就是 $x_0+k\delta x$ 和 $y_0+k\delta y$。

  　　　　推一下式子： $$a(x_1-x_0)+b(y_1-y_0)=0$$

$$a\delta x+b\delta y=0$$

$$\frac{a}{gcd(a,b)}\delta x+\frac{b}{gcd(a,b)}\delta y=0$$

$$\delta x=\frac{b}{gcd(a,b)},\delta y=\frac{a}{gcd(a,b)}$$

　  证毕。

 

　　1.1 逆元

对于正整数 $a$ 和 $m$，如果有 $ a \times x \equiv 1 \quad (mod\;m) $，那么把这个同余方程 $x$ 中的最小正整数解叫做 $a$ 模 $m$ 的逆元。

逆元一般用扩展欧几里得算法来求得，如果 $m$ 为素数，那么还可以根据费马小定理得到逆元为 $a^{m-2} \equiv m$。

推导过程如下：

$$ a^{m-1} \equiv 1 \quad (mod\;m) \Rightarrow a \times a^{m-2} \equiv 1 \quad (mod\;m) \Rightarrow a^{m-2} \equiv \frac{1}{a} \quad (mod\;m)$$

另外，逆元还有线性算法：

首先，$1^{-1} \equiv 1\;(mod\;p)$

然后，我们设 $p=k \times i+r , r<i , 1<i<p$ ，再将这个式子放到 $mod\;p$ 意义下，就有： $$ k \times i + r \equiv 0 \quad (mod\;p) $$

在等号两边同乘上 $i^{-1},r^{-1}$，就会得到：

$$ k\times r^{-1} +i^{-1} \equiv 0 \quad (mod\;p) $$

$$ i^{-1} \equiv -k \times r^{-1} \quad (mod\;p) $$

$$ i^{-1} \equiv - {\lfloor {\frac{p}{i}} \rfloor} \times (p\;mod\;i)^{-1} \quad (mod\;p) $$

递推代码就是:

1 inv[i]=(p-p/i)*inv[p%i];
2 (inv[i]+=mod)%=mod;

　　　　  因为有可能是负数，所以还要加上一个 $mod$

　    1.2 中国剩余定理

　　　　1.2.0 中国剩余定理

　　在《孙子算经》中有这样一个问题：“今有物不知其数，三三数之剩二（除以 3 余 2 ），五五数之剩三（除以 5 余3），七七数之剩二（除以 7 余 2 ），问物几何？”这个问题被称为“孙子问题”，该问题的一般解法国际上称为“中国剩余定理”。具体解法分三步：

　　1.找出三个数：从 3 和 5 的公倍数中找出被 7 除余 1 的最小数 15，从 3 和 7 的公倍数中找出被 5 除余 1 的最小数 21，最后从 5 和 7 的公倍数中找出除 3 余 1 的最小数 70 。

　　2.用 15 乘以 2（ 2 为最终结果除以 7 的余数），用 21 乘以 3（3 为最终结果除以 5 的余数），同理，用 70 除以 2（2 为最终结果除以 3 的余数），然后把三个乘积相加 $15 \times 2 + 21 \times 3 + 70 \times 2$ 得到和 233 .

　　3.用 233 除以 3,5,7, 三个数的最小公倍数 105，得到余数 23，即 $233 \% 105 =23$。这个余数 23 就是符合条件的最小数。

　　就这么简单。我们在感叹神奇的同时不禁想知道古人是如何想到这个方法的，有什么基本的数学依据吗？

　　我们将“孙子问题”拆分成几个简单的小问题，从零开始，试图揣测古人是如何推导出这个解法的。

　　首先，我们假设 $n_1$ 是满足除以 3 余 2 的一个数，比如 2,5,8 等等，也就是满足 $3 \times k +2\;(k \geq 0)$ 的一个任意数。同样，我们假设 $n_2$ 是满足除以 5 余 3 的一个数，$n_3$ 是满足除以 7 余 2 的一个数。

　　有了前面的假设，我们先从 $n_1$这个角度出发，已知 $n_1$ 满足除以 3 余 2，能不能使得 $n_1+n_2$ 的和仍然满足除以 3 余 2？进而使得  $n_1+n_2+n_3$ 的和仍然满足除以 3 余 2？

　　这就牵涉到一个最基本数学定理，如果有 $a\%b=c$，则有 $(a + k\times b)\%b = c\;(k \in N^*)$，换句话说，如果一个除法运算的余数为 $c$，那么被除数余 $k$ 倍的除数相加（或相减）的和（差）再与除数相除，余数不变。这是很好证明的。

　　以此定理为依据，如果 $n_2$ 是 3 的倍数， $n_1+n_2$ 就依然满足除以 3 余 2。同理，如果 $n_3$ 也是 3 的倍数，那么 $n_1+n_2+n_3$ 的和就满足除以 3 余 2。这是从 $n_1$ 的角度考虑的。再从 $n_2,n_3$ 的角度出发，我们可推导出以下三点：

　　1.为使 $n_1+n_2+n_3$ 的和满足除以 3 余 2，$n_2$ 和 $n_3$ 必须是 3 的倍数。

　　2.为使 $n_1+n_2+n_3$ 的和满足除以 5 余 3，$n_1$ 和 $n_3$ 必须是 5 的倍数。

　　3.为使 $n_1+n_2+n_3$ 的和满足除以 7 余 2，$n_1$ 和 $n_2$ 必须是 7 的倍数。

　　因此，为使 $n_1+n_2+n_3$ 的和作为“孙子问题”的一个最终解，需满足：

　　1. $n_1$ 除以 3 余 2，且是 5 和 7 的公倍数。

　　2. $n_2$ 除以 5 余 3，且是 3 和 7 的公倍数。

　　3. $n_3$ 除以 7 余 2，且是 3 和 5 的公倍数。

　　所以，孙子问题解法的本质是从 5 和 7 的公倍数中找一个除以 3 余 2 的数 $n_1$，从 3 和 7 的公倍数中找一个除以 5 余 3 的数 $n_2$，从 3 和 5 的公倍数中找一个除以 7 余 2 的数 $n_3$，再将三个数相加得到解。在求 $n_1,n_2,n_3$ 时又使用了一个小技巧，以 $n_1$ 为例，并非从 5 和 7 的公倍数中直接找一个除以 3 余 2 的数，而是先找一个除以 3 余 1 的数，再乘以 2。也就是先求出 5 和 7 的公倍数模 3 下的逆元，再用逆元去乘余数。

　　这里又有一个数学公式，如果 $a\%b=c$，那么 $$(a \times k) \%b =(a\%b + a\%b + a\%b + ... + a\%b)\%b = (c + c + c +...+ c)\%b = (k \times c) \%b\;(k>0)$$.也就是说，如果一个除法的余数为 $c$，那么被除数的 $k$ 倍与除数相除的余数为 $(k \times c) \%b$。展开式中已证明。

　　最后，我们还要清楚一点，$n_1+n_2+n_3$ 只是问题的一个解，并不是最小的解。如何得到最小解？我们只需要从中最大限度的减掉 3,5,7 的公倍数 105 即可。道理就是前面讲过的定理“如果 $a\%b=c$，则有 $(a - k \times b) \%b =c $”.所以 $(n_1+n_2+n_3)\%105$ 就是最终的最小解。　

　　　　　　这样一来就得到了中国剩余定理的公式：

设正整数$m_1,m_2,...,m_k$两两互素，则同余方程组

$$x \equiv a_1 \quad (mod\;m_1)$$

$$x \equiv a_2 \quad (mod\;m_2)$$

$$x \equiv a_3 \quad (mod\;m_3)$$

$$.$$

$$.$$

$$.$$

$$x \equiv a_k \quad (mod\;m_k)$$

有整数解。并且在模$M=m_1 \times m_2 \times ... \times m_k$下的解是唯一的，解为

$$x \equiv (a_1M_1M_1^{-1} + a_2M_2M_2^{-1} + ... + a_kM_kM_k^{-1}) \quad (mod\;M)$$

其中$M_i = M/m_i$，而$M_i^{-1}$为$M_i$模$m_i$的逆元。

 

 1 int exgcd(int a,int b,int &x,int &y){
 2     if(!b){
 3         x=1;y=0;
 4         return a;
 5     }
 6     int c=exgcd(b,a%b,x,y);
 7     int t=x;
 8     x=y;
 9     y=t-a/b*y;
10     return c;
11 }
12 
13 int inv(int a,int b){
14     int x,y;
15     exgcd(a,b,x,y);
16     return (x%b+b)%b;
17 }
18 
19 int CRT(){
20     int M=1,ans=0;
21     for(int i=1;i<=n;i++) M*=m[i];
22     for(int i=1;i<=n;i++)
23         (ans+=(M/m[i])*inv(M/m[i],m[i])*a[i])%=M;
24     return ans;
25 }

　　　　1.2.1 扩展中国剩余定理

普通的中国剩余定理要求所有的 $m_i$ 互质，那么如果不互质呢，怎么求解同余方程组？

这种情况就采用两两合并的思想，假设要合并如下两个方程:

$$x = r_1 + m_1 \times x_1$$

$$x = r_2 + m_2 \times x_2$$

那么得到：

$$r_1 + m_1 \times x_1 = a_2 + m_2 \times x_2 \Rightarrow m_1 \times x_1 + m_2 \times x_2 = r_2 - r_1$$

我们要求出一个最小的 $x$ 使它满足：

$$x = r_1 + m_1 \times x_1 = r_2 + m_2 \times x_2$$

那么 $x_1$ 和 $x_2$ 就要尽可能的小，于是我们用扩展欧几里得算法求出 $x_1$ 的最小正整数解，将它代回 $a_1 + m_1 \times x_1$，得到 $x$ 的一个特解 $x'$，当然也是最小正整数解。

所以 $x$ 的通解一定是 $x'$ 加上 $lcm(m_1,m_2) \times k$（$k$ 是一个常数），这样才能保证 $x$ 模 $m_1$ 和 $m_2$ 的余数是 $r_1$ 和 $r_2$。由此，我们把这个 $x'$ 当做新的方程的余数，把 $lcm(m_1,m_2)$ 当做新的方程的模数。

合并完成：

$$x \equiv x' \quad (mod\;lcm(m_1,m_2))$$

poj 2891代码


#include<cstdio>
#define int long long
const int N=1e5+5;

int n,a[N],r[N];int exgcd(int a,int b,int &x,int &y){
    if(!b){x=1,y=0;return a;}
    int d=exgcd(b,a%b,y,x);y-=a/b*x;return d;
}

int CRT(){
    int A=a[1],R=r[1];
    for(int i=2;i<=n;i++){
        int x,y,d=exgcd(A,a[i],x,y),z=A/d*a[i];
        if((r[i]-R)%d) return -1;
        x*=(r[i]-R)/d;x=(x%(a[i]/d)+a[i]/d)%(a[i]/d);//求出通解之后算出最小整数解
        (R+=x*A)%=z;A=z;
    } return R;
}

signed main(){
    while(~scanf("%lld",&n)){
        for(int i=1;i<=n;i++) scanf("%lld%lld",&a[i],&r[i]);
        printf("%lld\n",CRT());
    } return 0;
}

 

 

 

 

1.3 高斯消元


　　假设有一个线性方程组是长这样的：

$$\begin{Bmatrix}3x&+&2y&+&z&=&10\\5x&+&y&+&6z&=&25\\2x&+&3y&+&4z&=&20\end{Bmatrix}$$

　　$emmm$ 这就是一个很简单的三元一次方程，让我们想想常规方法该怎么做（先不谈 $ code$ ）

　　初中老师说过：我们可以加减消元或者代入消元,但是我们需要在程序里实现的时候，需要一种有规律可循的算法。所以我们选择加减消元，但用代入消元回带。

　　整体思路就是我们可以先在某一个式子里，用这个式子的 $x$ 消去其他式子里的 $x$ ,然后在剩下的两个式子里再选择一个式子里的 $y$ ，用这个 $y$ 消去最后剩下的式子里的$ y$ 。那么现在最后一个方程里就只有一个未知数 $z$ 了。倘若 $z$ 的系数是 $1$ ，那么我们就可以直接得出答案来了（别觉得这句话是废话）。

　　比如刚才这个方程，我们用第二个式子去消1、3式里的 $x$ ：

$$\begin{Bmatrix} 0 \times x&+& \frac{7}{5}y&+&(-\frac{13}{5})&=&-5\\5x&+&y&+&6z&=&25\\0\times x&+&\frac{13}{5}&+&\frac{8}{5}z&=&10\end{Bmatrix}$$

　　整理之后再用第三个式子里的 $y$ 消去第一个式子里的 $y$ (注意，由于第二个式子作为消元用式，所以接下来的运算不再考虑二式):

$$\begin{Bmatrix}0\times y&+&(-\frac{225}{65}z)&=&-\frac{135}{13}\\ \frac{13}{5}y&+&\frac{8}{5}z&=&10\end{Bmatrix}$$

　　那么我们发现在 $1$ 式中只剩下一个未知数了，那么就可解得： $$z=3$$
　　带回三式里解出
$$y=2$$
　　再将 $x$、$y$ 带回最早被消掉的二式里，解得
$$x=1$$
　　好像这个方法再数学逻辑上讲是特别麻烦的，但是却是一个通用性强的方法

　　那么放在程序实现上来讲，我们可以先用一个 $n\times (n+1)$ 的矩阵来记录每一个式子的系数以及结果。譬如上式就可以用矩阵表示成这样：
$$\begin{Bmatrix} 3&2&1&\mid&10\\5&1&6&\mid&25\\2&3&4&\mid&20\end{Bmatrix}$$
　　左边记录系数，右边记录每个式子的结果。
 
　　那么首先我们需要用第一列中（所有的 $x$ 中）系数最大的来消其他两个式子。而这里为了方便起见，我们将这个选中的系数置为 $1$ ，方便上例中地不断带回原式的操作（这样在回带的时候就可以不考虑原本的系数了）。
 
代码：

 

 

#include<cmath>
#include<cstdio>
#include<algorithm>
#define db double
#define eps 1e-8

int n;
db ans[105];
db a[105][105];

signed main(){
    scanf("%d",&n);
    for(int i=1;i<=n;i++){
        for(int j=1;j<=n+1;j++)
            scanf("%lf",&a[i][j]);
    }
    for(int i=1;i<=n;i++){
        int r=0;
        for(int j=i;j<=n;j++){
            if(fabs(a[j][i]-a[r][i])) r=j; 
        }
        if(!r){
            printf("No Solution");
            return 0;
        }
        std::swap(a[i],a[r]);
        for(int j=i+1;j<=n;j++){
            db tmp=a[j][i]/a[i][i];
            for(int k=i;k<=n+1;k++)
                a[j][k]-=tmp*a[i][k];
        }
    }
    ans[n]=a[n][n+1]/a[n][n];
    for(int i=n-1;i;i--){
        for(int j=n;j>i;j--)
            a[i][n+1]-=ans[j]*a[i][j];
        ans[i]=a[i][n+1]/a[i][i];
    }
    for(int i=1;i<=n;i++) printf("%.2lf\n",ans[i]);
    return 0;
}

 

　　　　　upd on 2018.7.8

　　　　辗转相减高斯消元

　　　　　　刚做这题学到了个新姿势。

　　　　　　高斯消元的过程中如果要求模一个非质数该怎么做呢？

　　　　　　考虑更相减损术。对于元素 a,b，想把 a 变成 0，那么就要不断地令 a=a%b,swap(a,b)，直到有一个数变成 0 为止。

luogu4111 code

#include<cstdio>
#include<cctype>
#include<cstring>
#include<iostream>
const int mod=1e9;
#define int long long
#define min(A,B) ((A)<(B)?(A):(B))
#define max(A,B) ((A)>(B)?(A):(B))
#define swap(A,B) ((A)^=(B)^=(A)^=(B))

int n,m,cnt;
int mp[12][12];
int a[105][105],b[105][105];
int dx[]={-1,0,1,0},dy[]={0,1,0,-1};

int getint(){
    int x=0,f=0;char ch=getchar();
    while(!isdigit(ch)) f|=ch=='-',ch=getchar();
    while(isdigit(ch)) x=(x<<1)+(x<<3)+(ch^48),ch=getchar();
    return x;
}

signed main(){
    n=getint(),m=getint();
    for(int i=1;i<=n;i++){
        for(int j=1;j<=m;j++){
            char ch;std::cin>>ch;
            if(ch=='*') continue;
            mp[i][j]=++cnt;
        }
    }
    for(int i=1;i<=n;i++){
        for(int j=1;j<=m;j++){
            if(!mp[i][j]) continue;
            for(int k=0;k<4;k++){
                int nx=i+dx[k];
                int ny=j+dy[k];
                if(nx<1 or nx>n or ny<1 or ny>m) continue;
                if(!mp[nx][ny]) continue;
                a[mp[i][j]][mp[i][j]]++;
                a[mp[i][j]][mp[nx][ny]]--;
            }
        }
    }
    int ans=1;
    for(int i=1;i<cnt;i++){
        for(int j=i+1;j<cnt;j++){
            while(a[j][i]){
                int t=a[i][i]/a[j][i];
                for(int k=i;k<cnt;k++){
                    a[i][k]-=t*a[j][k];
                    a[i][k]=(a[i][k]+mod)%mod;
                    swap(a[j][k],a[i][k]);
                }
                ans*=-1;
            }
        }
        ans=ans*a[i][i]%mod;
        ans=(ans+mod)%mod;
    }
    printf("%lld\n",ans);
    return 0;
}

 

 

　   1.4 Baby Step，Giant Step

 　　　　　1.4.0 BSGS

　　　　　  问题：给定整数 $a,b,p$，其中 $gcd(a,p)=1$，求一个最小的非负整数 $x$，使 $a^x\equiv b\;(mod\;p)$

　　　　　由欧拉定理易得 $x\in [0,\phi(p)-1]$。设 $x=x_1+x_2$，那么 $a^{x_1}\equiv b\times (a^{x_2})^{-1}$。

　　　　　建立一个 $Hash$ 表，存下 $a^q=n$ 的情况下 $n$ 对应的最小正整数 $q$。

　　　　　我们枚举 $x_2$，同时算出 $m=b\times (a^{x_2})^{-1}$，然后查表，看是否有整数 $p$ 对应 $a^p=m$，如果有，那么可以用 $p+x_2$ 更新答案。

　　　　　据均值不等式，枚举 $x_2=k\times \sqrt{\phi(p)}$ 时最优。

　　　　　upd on 2018.7.4  改了一个小错误感谢sx小伙伴提醒

　　　　1.4.1 扩展BSGS

　　　　　要求 $gcd(a,p)=1$ 太特殊了，对于更一般的情况，如果 $a,p$ 不互质，那该怎么办呢？

　　　　　由上面互质的做法得到启发，我们是不是可以把它变成互质的，然后继续做呢？

　　　　　具体做法如下：

　　　　1.当 $gcd(a,p)\ne 1$ 时， $a,b,p$ 同时除以 $gcd(a,p)$，直到  $b\nmid gcd(a,p)$ 或者 $gcd(a,p)=1$，并记录 $num$ 为除了多少次。

　　　　2.如果 $b\nmid gcd(a,p)$ 并且 $b!=1$ ，那么同余方程显然无解。

　　　　3.将式子转化为 $a_0^{x-num}\equiv b_0\times (a_1)^{-1}\;(mod\;p_0)$

　　　　4.观察到这已经是 $BSGS$ 可以解决的范畴了。但是还有一点要注意的是直接 $BSGS$ 的话可能会漏掉 $x=[0,1...num-1]$ 的解，所以这里暴力跑一下快速幂　　判断即可。

　　　　　　　SPOJ3105 Code

　　　　　　　

#include<map>
#include<cmath>
#include<cstdio>
#include<cstring>
#define int long long
#define min(A,B) ((A)<(B)?(A):(B))

int a,b,p;
std::map<int,int> mp;

int gcd(int x,int y){
    if(!y) return x;
    return gcd(y,x%y);
}

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

int inv(int a,int b){
    int x,y;
    exgcd(a,b,x,y);
    return (x%b+b)%b;
}

int ksm(int a,int b,int c){
    int ans=1;
    while(b){
        if(b&1) ans=ans*a%c;
        a=a*a%c;
        b>>=1;
    }
    return ans;
}

int phi(int n){
    int ans=n,sq=sqrt(n);
    for(int i=2;i<=sq;i++){
        if(n%i==0){
            ans=ans/i*(i-1);
            while(n%i==0)
                n/=i;
        }
    }
    if(n>1)
        ans=ans/n*(n-1);
    return ans;
}

int solve(){
    int aa=a,bb=b,pp=p,a2=1,t;
    int num=0;
    while((t=gcd(a,p))!=1 and b%t==0){
        a2*=a/t; b/=t; p/=t; 
        num++; a2%=p;
    }
    b=b*inv(a2,p)%p;
    if(b==1)
        return num;
    int d=gcd(a,p);
    if(b%d)
        return -1;
    int now=phi(p);
    int block=sqrt(now);
    int k=ksm(a,block,p);
    int last=1;
    for(int i=0;i<=num;last=last*a%pp,i++){
        if(last==bb)
            return i;
    }
    last=1;
    for(int i=0;i<=block;last=last*a%p,i++){
        int q=last;
        if(!mp.count(q))
            mp[q]=i;
    }
    last=1;
    for(int i=0;i<=block;i++,last=last*k%p){
        int q=b*inv(last,p)%p;
        if(mp.count(q))
            return i*block+mp[q]+num;
    }
    return -1;
}

signed main(){
    while(1){
        scanf("%lld%lld%lld",&a,&p,&b);
        if(!a)
            return 0;
        b%=p;
        int k=solve();
        printf(k==-1?"No Solution\n":"%lld\n",k);
        mp.clear();
    }
}

 

 

 

2 组合数学

　　2.1 Catalan数

　　　　　2.1.0 带限制条件的路径总数

　首先我们来看一个问题：

　在一个平面直角坐标系中，只能往右或往上走一个单位长度，问有多少种不同的路径可以从左下角 (1, 1) 走到右上角 (n, n)，并且要求路径不能经过直线 y = x 上方的点，下图中的路径都是合法的（图片来源 Wikipedia）

　如果没有限制条件，那么从左下角走到右上角一共有 2n 步，有 n 步是往右，另外 n 步是往上，那么路径方案数就是 2n 步中选择 n 步往右，一共有 {2n} \choose {n} （即 C_{2n}^n）种方案

　那么我们考虑一下这里面有多少种方案是不合法的

　首先对于每一种不合法的方案，它的路径一定与 y = x + 1 有交。我们找到它与 y = x + 1 的第一个交点，然后将这个点后面部分的路径关于 y = x + 1 做一个对称。由于原来路径到达 (n, n)，新的对称之后的路径就会到达 (n - 1, n + 1)。这样我们把每一种不合法方案都对应到了一条从 (1, 1) 到 (n - 1, n + 1) 的路径，现在再来看是否每一条这样的路径都能对应到一种不合法方案，如果是，那么这就建立了一个一一映射的关系，也就是它们的方案总数相同。这是肯定的，因为每一条这样的路径必定与 y = x + 1 有交，那么对称回去，就得到一条不合法方案

　由于从 (1, 1) 到 (n - 1, n + 1) 的路径有 {n - 1 + n + 1} \choose {n - 1} 条，那么合法的　　方案就是

C_n = {{2n} \choose {n}} - {{2n} \choose {n - 1}} = \frac{1}{n + 1} {{2n} \choose {n}}

　我们把这个方案数记为 C_n，这就是著名的 Catalan 数

　我们来看看它的前几项（n 从 0 开始）

　1, 1, 2, 5, 14, 42, 132, 429, 1430, 4862, 16796, 58786, 208012, 742900, 2674440, 9694845, 35357670, 129644790, 477638700, 1767263190

 

　2.1.1 括号序列计数

　再来看一个问题：有多少种不同的长度为 n 的括号序列？

　首先一个括号序列是指 (), ()(), (())() 这样的由括号组成的序列，并且没有左右括号无法匹配的情况

　我们可以将长度为 2n 的括号序列映射成刚刚所说的路径：首先对于左括号，那么就向右走一个单位长度，对于右括号，那么就向上走一个单位长度，由于括号序列合法，那么每次向右走的次数不会少于向上的次数，也就是这条路径不会在 y = x 之上。再考虑每一条这样的路径，也能够对应到一种合法的括号序列，因此，长度为 2n 的括号序列的方案数就是 C_n

　2.1.2 出栈顺序

　现在来考虑你有 n 个元素（元素之间是没有区别的）和一个栈，每次可以将一个元素入栈，或者将栈顶元素弹出，问有多少种可能的操作序列，这可以将问题对应成括号序列，入栈为左括号，出栈为右括号，因此方案数也是 C_n

   2.1.3 排队问题

　现在有 2n 个人，他们身高互不相同，他们要成两排，每一排有 n 个人，并且满足每一排必须是从矮到高，且后一排的人要比前一排对应的人要高，问有多少种排法

　我们考虑先把这些人从矮到高排成一排，那么现在来分配哪个人在前，哪个人在后，例如有 6个人，身高是 1, 2, 3, 4, 5, 6

　那么我们用 1 表示这个人应该在后排，0 表示这个人应该在前排，例如说 100110 表示两排分别是 2, 3, 6 和 1, 4, 5 这是不合法的

　那么合法方案应该是怎么样的呢？后排要比前排对应的人高，也就是说 0 的出现次数在每一个地方都不应该小于 1，这恰恰又是一个括号序列，因此，方案仍然是 Catalan 数

　2.1.4 二叉树计数

　现在你需要统计有多少种不同的 n 个结点的二叉树

 

　图上的是 3 个结点的二叉树，一共有 5 种方案

　朴素的想法是由于二叉树是递归定义的，可以考虑使用递推方法

　我们可以用 f_n 表示有 n 个结点的二叉树的方案数（空树算一种，即 $f_0=1$ ），那么枚举子树大小可以得到方程

f_n = \sum_{i = 0}^{n - 1} f_if_{n - i - 1}

　如果直接计算，你需要 \mathcal O(n^2) 的时间

　现在我们换一个角度来想，对这棵二叉树进行遍历，并且考虑一个括号序列，当第一次遇到这个结点的时候，在括号序列末尾添加一个左括号，在从左子树回到这个结点的时候，在括号序列中添加一个右括号，这样，将每一种不同的二叉树都对应到了一种不同的括号序列，同样对于每一种不同的括号序列都可以找到对应的一种不同的二叉树，因此，有 n 个结点的二叉树的数量也是 C_n。

2.2 $Lucas$ 定理 (摘自Navi-Awson)

　　2.2.0 $Lucas$

　　　　$C_n^m\equiv C^{\lfloor{\frac{m}{p}}\rfloor}_{\lfloor{\frac{n}{p}}\rfloor}\times C_{n\;mod\;p}^{m\;mod\;p}\quad (mod\;p)$，$p$ 是质数。

　　证明：需要用到多项式同余的那套理论，然而我并不会。

　　2.2.1 扩展 $Lucas$

　　　　一般的 $Lucas$ 是在模数 $P$ 是质数的条件下适用的。我们来考虑 $P$ 不是质数的条件。

　　　　我们对 $P$ 进行唯一分解，记 $P=p_1^{k_1}p_2^{k_2}...p_q^{k_q}$，由于形同 $p_i^{k_i}$ 的部分是互质的，显然我们可以用 $CRT$ 合并。

　　　　列出方程组：

$$ans\equiv c_1\quad (mod\;p_1^{k_1})$$

$$ans\equiv c_2\quad (mod\;p_2^{k_2})$$

$$...$$

$$ans\equiv c_q\quad (mod\;p_q^{k_q})$$

　　　　对于每个 $c_i$，表示 $C_n^m$ 在 $mod\;p_i^{k_i}$ 下的结果。由解的唯一性，我们可以证明这个 $ans$ 就是我们要求的。

　　　　根据 $C_n^m=\frac{n!}{m!(n-m)!}$ 我们只要求出 $n!\;mod\;p_i^{k_i},m!\;mod\;p_i^{k_i},(n-m)!\;mod\;p_i^{k_i}$，再用逆元的那套理论就可以求 $c_i$ 了。

　　　　考虑如何求 $n!\;mod\;p_i^{k_i}$ 。容易发现 $n!= (\Pi_{j=1,p_i\nmid j}^{n}j )\times  (p_i^{\lfloor{\frac{n}{p_i}}\rfloor} )\times  (\lfloor{\frac{n}{p_i}}\rfloor !)$

　　　　上述式子分为三个部分，第一个部分显然在模 $p_i^{k_i}$ 下，是以 $p_i^{k_i}$ 为周期的。可以周期内找循环节算，周期外的暴力算；第二部分可以直接算，第三部分可以递归求解。

　　　　另外注意的是求组合逆元的是，存在阶乘中的某一个数可能还有 $p_i$ 这个质因子，不能直接算。直接把 $p_i$ 全部提出来，最后求完逆元再补回去。求 $n!$ 内质因子 $p$ 的个数可以用 $\sum _{i=1}^{\infty}\lfloor{\frac{n}{p^i}}\rfloor$ 来求。

　　　　代码：

int ksm(int a,int b,int x){
    int ans=1;a%=x;
    while(b){
        if(b&1) ans=(ans*a)%x;
        a=(a*a)%x;
        b>>=1;
    }
    return ans;
}

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

int inv(int a,int b){
    int x,y;
    exgcd(a,b,x,y);
    x=(x%b+b)%b;
    return x;
}

int solve(int a,int cur){
    if(!a) return 1;
    int ans=1;
    for(int i=2;i<c[cur];i++){
        if(i%p[cur])
            ans=(ans*i)%c[cur];
    }
    ans=ksm(ans,a/c[cur],c[cur]);
    for(int i=2;i<=a%c[cur];i++){
        if(i%p[cur])
            ans=(ans*i)%c[cur];
    }
    return ans*solve(a/p[cur],cur)%c[cur];
}

int C(int cur,int a,int b){
    if(a<b) return 0;
    int x=solve(a,cur);
    int y=solve(b,cur);
    int z=solve(a-b,cur);
    int k=0;
    for(int i=a;i;i/=p[cur])
        k+=i/p[cur];
    for(int i=b;i;i/=p[cur])
        k-=i/p[cur];
    for(int i=a-b;i;i/=p[cur])
        k-=i/p[cur];
    int ans=(x*inv(y,c[cur])*inv(z,c[cur])*ksm(p[cur],k,c[cur]))%c[cur];
    return ans;
}

int CRT(){
    int ans=0;
    for(int i=1;i<=cnt;i++)
        ans=(ans+inv(mod/b[i],b[i])*(mod/b[i])*a[i])%mod;
    return ans;
}

 

posted @ 2018-04-22 16:51  YoungNeal  阅读(1020)  评论(6)  编辑  收藏  举报