算法·数学知识

数论一定要算时间复杂度以防超时

质数

质数的判定——试除法

时间复杂度一定为

O(x)O(\sqrt x) 

1
2
3
4
5
6
7
8
bool is_prime(int x)
{
    if (x < 2) return false;
    for (int i = 2; i <= x / i; i ++ )//防溢出
        if (x % i == 0)
            return false;
    return true;
}

试除法分解质因数

从小到大枚举所有约数

底数是指指数，指数是指每个底数出现的次数

n中最多只包含一个大于sqrt(n)的质因子

时间复杂度

O(x)−O(logn)O(\sqrt x)-O(logn) 

1
2
3
4
5
6
7
8
9
10
11
12
void divide(int x)
{
    for (int i = 2; i <= x / i; i ++ )
        if (x % i == 0)//i一定是质数
        {
            int s = 0;
            while (x % i == 0) x /= i, s ++ ;
            cout << i << ' ' << s << endl;
        }
    if (x > 1) cout << x << ' ' << 1 << endl;
    cout << endl;
}

朴素筛法求素数（埃筛）

时间复杂度是nlogn

1
2
3
4
5
6
7
8
9
10
11
12
13
14
int primes[N], cnt;     // primes[]存储所有素数
bool st[N];         // st[x]存储x是否被筛掉

void get_primes(int n)
{
    for (int i = 2; i <= n; i ++ )
    {
        if (st[i]) continue;//如果被筛过，继续找质数
        //没有被筛过，是质数
        primes[cnt ++ ] = i;//把质数i加入质数数组里面并且指向下个位置
        for (int j = i + i; j <= n; j += i)//把i的倍数删掉
            st[j] = true;
    }
}

线性筛法求素数

质数定理：1-n中有n/lnn个质数

时间复杂度是O(nloglogn)

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
int primes[N], cnt;     // primes[]存储所有素数
bool st[N];         // st[x]存储x是否被筛掉

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[primes[j] * i] = true;//只用筛质数的倍数，当一个数不是质数的时候就不需要筛掉他的所有倍数（因为质数会筛掉），这样可以避免重复。每次把质数和i的乘积筛掉
            if (i % primes[j] == 0) break;//意味着primes[j]一定是i的最小质因子，因为是从小到大枚举的质因子，pj也已i当时pj*i的最小质因子
            //如果i%pj！=0，说明pj一定小于i的所有质因子。pj也一定是pj*i的最小质因子
        }
    }
}

朴素筛法和线性筛法的最大区别在于，朴素筛法可能出现重复筛的过程（一个数有多个质因数），因为是从1-n找出质数然后在进行筛除，而线性筛法是将所有质数先进行，然后乘以倍数，找出质数在进行下一轮质数排查。一个是横向一个是纵向

约数

int范围内约数最多的大概是1500个

试除法求约数

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
vector<int> get_divisors(int x)
{
    vector<int> res;
    for (int i = 1; i <= x / i; i ++ )
        if (x % i == 0)
        {
            res.push_back(i);
            if (i != x / i) res.push_back(x / i);//防止重复输入
        }
    sort(res.begin(), res.end());
    return res;
}
int main(){
    int n;
    cin >> n;
    while(n--){
        int x;
        cin>>x;
        auto res= get_divisors(x);
        for(auto t :res)cout<<t<<' ';
        cout<<endl;
    }
    return 0;
}

约数个数和约数之和

1
2
3
如果 N = p1^c1 * p2^c2 * ... *pk^ck
约数个数： (c1 + 1) * (c2 + 1) * ... * (ck + 1)
约数之和： (p1^0 + p1^1 + ... + p1^c1) * ... * (pk^0 + pk^1 + ... + pk^ck)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
//约数个数
//先把乘积因式分解求出来
#include <iostream>
#include <algorithm>
#include <unordered_map>
using namespace std;
typedef long long LL;
const int mod=1e9+7;
int main(){
    int n;
    cin>> n;
    unordered_map<int,int>primes;
    while(n--)
    {
        int x;
        cin >> x;
        for(int i =2;i<=x/i;i++)
            while(x%i==0){
                x/=i;
                primes[i]++;//此处primes[i]对应的i是质数即底数，i是first元素，primes[i]的值对应的是质数出现的次数即指数，对应的是second元素，即first对应的是p，second对应的是c
            }
        if(x>1)primes[x]++;
    }
    LL res =1;
    for(auto prime:primes) res= res*(prime.second+1)%mod;
    cout<<res<<endl;
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
//求约数之和
//约数个数
//先把乘积因式分解求出来
#include <iostream>
#include <algorithm>
#include <unordered_map>
using namespace std;
typedef long long LL;
const int mod=1e9+7;
int main(){
    int n;
    cin>> n;
    unordered_map<int,int>primes;
    while(n--)
    {
        int x;
        cin >> x;
        for(int i =2;i<=x/i;i++)
            while(x%i==0){
                x/=i;
                primes[i]++;
            }
        if(x>1)primes[x]++;
    }
    LL res =1;
    for(auto prime:primes) {
         int p = prime.first,a=prime.second;
        LL t = 1;
        while(a--) t=(t*p+1)%mod;
        res=res*t%mod;
    }
    cout<<res<<endl;
    return 0;
}

欧几里得算法（辗转相除法）

d能整除a，d能整除b，d能整除a+b

gcd(a,b)=gcd(b,a mod b)

1
2
3
4
5
6
int gcd(int a, int b)
{
    return b ? gcd(b, a % b) : a;
    //三目运算符 表达式1？表达式2:表达式3; 1条件，2if 3else
    //b不是0，返回b和a%b，是的话返回a
}

欧拉函数

互质是公约数只有1的两个整数，叫做互质整数

1-N中与N互质的数的个数被称为欧拉函数,记为Φ(N)。若在算术基本定理中：

N=p1α1∗p2α2∗...∗pmαmN=p1^{α^1}*p2^{α^2}*...*pm^{α^m}

则：

Φ(N)=N∗p1−1p1∗p2−1p2∗...∗pm−1pmΦ(N)=N*\frac {p1-1}{p1}*\frac {p2-1}{p2}*...*\frac {pm-1}{pm}

Φ(N)=N∗(1−1p1)∗(1−1p2)∗...∗(1−1pm)Φ(N)=N*(1-\frac{1}{p1})*(1-\frac{1}{p2})*...*(1-\frac{1}{pm})

  • 从1-N中去掉p1,p2…pk的所有倍数(都含有公约数质数本身)
  • 加上所有pi*pj的倍数（所有两个质数的组合）
  • 减去所有pi *pj *pk的倍数
  • 加上pi * pj * pk * pl的倍数
  • 依此类推

欧拉函数

时间复杂度是

O(n)O(\sqrt n) 

1
2
3
4
5
6
7
8
9
10
11
12
13
int phi(int x)
{
    int res = x;
    for (int i = 2; i <= x / i; i ++ )
        if (x % i == 0)
        {
            res = res / i * (i - 1);
            while (x % i == 0) x /= i;//分解约数
        }
    if (x > 1) res = res / x * (x - 1);//最后一个带2/x的约数

    return res;//欧拉函数
}

筛法求欧拉函数

求1-n中每一个数的欧拉函数

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
int primes[N], cnt;     // primes[]存储所有素数
int euler[N];           // 存储每个数的欧拉函数
bool st[N];         // st[x]存储x是否被筛掉


void get_eulers(int n)
{
    euler[1] = 1;
    for (int i = 2; i <= n; i ++ )
    {
        if (!st[i])
        {
            primes[cnt ++ ] = i;//找出质数
            euler[i] = i - 1;//质数的欧拉函数就是本身-1，因为他和他前面的数没有公因数，互质
        }
        for (int j = 0; primes[j] <= n / i; j ++ )
        {
            int t = primes[j] * i;//*i * prime[j]的最小质因子就是prime[j]，则prime[j] * i就不是质数了，把它true掉
            st[t] = true;//筛掉
            if (i % primes[j] == 0)//这里求i * prime[j]的欧拉函数
            {/*这里求i * prime[j]的欧拉函数，根据欧拉公式推导如下:
1.euler[i]  = i * (1 - 1 / p1) * (1 - 1 / p2) * (1 - 1 / p3) * ... * (1 - 1 / pk)
2.euler[i * prime[j]] = (i * prime[j]) * (1 - 1 / p1) * (1 - 1 / p2) * (1 - 1 / p3) * ... * (1 - 1 / pk)
                由此可得euler[i * prime[j]] = prime[j] * euler[i]*/
                euler[t] = euler[i] * primes[j];
                break;
            }
            /*此时prime[j]不是i的质因子，则推导如下：
1.euler[     i      ]  =      i        * (1 - 1 / p1) * (1 - 1 / p2) * (1 - 1 / p3) * ... * (1 - 1 / pk)
2.euler[i * prime[j]] = (i * prime[j]) * (1 - 1 / p1) * (1 - 1 / p2) * (1 - 1 / p3) * ... * (1 - 1 / pk) * (1 - 1 / prime[j])
            由此可得euler[i * prime[j]] = euler[i] * (prime[j] - 1);*/
            euler[t] = euler[i] * (primes[j] - 1);
        }
    }
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
#include <iostream>
using namespace std;
const int N =1000010;
int primes[N],cnt;
typedef long long LL;
bool st[N];
int phi[N];
LL get_eluers(int n){
    phi[1]=1;
    for(int i = 2;i<=n;i++){
        if(!st[i]){//是质数 
            primes[cnt++]=i;
            phi[i]=i-1;
        }
        for(int j = 0;primes[j]<=n/i;j++){
            st[primes[j]*i=true;
            if(i%primes[j]==0){
                phi[primes[j]*i]=phi[i]*primes[j];
                break;
            }
            phi[primes[j]*i]=phi[i]*(primes[j]-1);
        }
    }
    LL res=0;
    for(int i = 1;i<=n;i++) res+=phi[i];
    return res;
}
int main(){
    int n;
    cin >> n;
    cout<<get_eluers(n);<<endl;
    return 0;
}

欧拉定理

  • ≡同余，两个整数除以同一个整数，若得相同余数，则二整数同余。

若a与n互质，则有

aφ(n)≡1(modn)a^{φ(n)}≡1(mod n)

当n是质数的时候，有：

an−1≡1(modn)a^{n-1}≡1(mod n)

快速幂

1
2
3
4
5
6
7
8
9
10
11
12
求 m^k mod p，时间复杂度 O(logk)。
int qmi(int m, int k, int p)
{
    int res = 1 % p, t = m;
    while (k)
    {
        if (k&1) res = res * t % p;
        t = t * t % p;
        k >>= 1;
    }
    return res;
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
#include <iostream>
#include <algorithm>
using namespace std;
typedef long long LL;
int qmi(int a,int k,int p){//返回的就是a^k%p的结果
    int res = 1;
    while(k){
        if(k&1) res = (LL)res*a%p;
        k>>=1;
        a=(LL)a*a %p;
    }
    return res;
}
int main(){
     ios::sync_with_stdio(false);
     cin.tie(0);
     itn n;
     cin >> n;
     while(n--){
     int a,k,p;
     cin >>a >> k>>p;
     cout<<qmi(a,k,p);
     }
     return 0;
}

逆元

若整数b,m互质，并且b|a，则存在一个整数x，使得

ab≡a∗x(modm)\frac {a}{b}≡a*x(mod m)

则x叫做b的模m乘法逆元，记为b^-1(mod m)

ab≡a∗b−1(modm)\frac ab≡a*b^{-1}(modm)

b存在乘法逆元的充要条件是b与模数m互质，当模数m为质数时，b^(m-2)即为b的乘法逆元

b如果是m的倍数一定无解，即b%m==0无解。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
#include <iostream>
using namespace std;
typedef long long LL;
int qmi(int a,int k,int p){//返回的就是a^k%p的结果
    int res = 1;
    while(k){
        if(k&1) res = (LL)res*a%p;
        k>>=1;
        a=(LL)a*a %p;
    }
    return res;
}
int main(){
    int n;
    cin >> n;
    while(n--){
        int a,k,p;
        cin >> a>>p;
        int res = qmi(a,p-2,p);
        if(a%p) cout<<res;
        else cout<<"impossible";dddd
    }
}

扩展欧几里得算法

裴蜀定理

有一对正整数a,b，那么存在非零整数x,y，使得ax+by=gcd(a,b)(gcd是最大公约数)，且为a和b能凑出来的最小的正整数，若一个数能整除gcd(a,b)则有解，否则无解

  • 整数中的裴蜀定理：

    对于任意两个整数a、b，假设d是他们的最大公约数，ax+by=m，当且仅当m是d的倍数时有整数解。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
// 求x, y，使得ax + by = gcd(a, b)
int exgcd(int a, int b, int &x, int &y)
{
    if (!b)
    {
        x = 1; y = 0;/*如果b=0，则gcd(a, b) = 1 * a + 0 * b*/
        return a;
    }
    int d = exgcd(b, a % b, y, x);   /*此时我们假设已经求出来了x, y，ax + by = d，根据辗转相除法（欧几里得算法）分析：
    此时为gcd(a, b)， 则下一循环为gcd(b, a % b)，此时我们求它的x1, y1可得
    b * x1 + (a % b) * y1 
    = b * x1 + (a - (a / b) * b) * y1 
    = b * x1 + a * y1 – (a / b) * b * y1 
    = a * y1 + b * (x1 – a / b * y1) 
    我们发现: x = y1, y = x1 - a / b * y1因此我们根据这个规律可以推导出x和y*/ 
  
    y -= (a/b) * x;
    return d;//d是最大公约数
}
int main(){
    int n;
    cin >> n;
    while(n--)
    {
        int a,b,x,y;
        cin>> a >> b;
        exgcd(a,b,x,y);
        cout<<x<y;
        
    }
    return 0;
}

求余数方程

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
// 求x, y，使得ax + by = gcd(a, b)
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;//d是最大公约数
}
int main(){
    int n;
    cin >> n;
    while(n--)
    {
        int a,b,m;
        cin>> a >> b>>m;
        int x,y;
        int d =exgcd(a,m,x,y);
        if(b%d) puts("impossible");
        else cout<<(LL)x*(b/d)%m;
    }
    return 0;
}

中国剩余定理

假设整数m1,m2,…,mn两两互质，则对任意的整数a1,a2,…,an,方程组(S)有解

M=m1∗m2∗m3∗...∗mkM=m_1*m_2*m_3*...*m_k

Mi=MmiM_i= \frac {M}{m_i}

Mi−1表示Mi mod miM_i^{-1}表示M_i ~mod~m_i

x=a1∗M1∗M1−1+a2∗M2∗M2−1+...+aK∗MK∗MK−1x = a_1*M_1*M_1^{-1}+a_2*M_2*M_2^{-1}+...+a_K*M_K*M_K^{-1}

也被称为中国剩余定理

  • 求逆即Mi^-1相当于解ax≡1（mod m），用扩展欧几里得算法可以求

  • a mod b = c 等价于 (a + nb) mod b = c;
a mod b = c 等价于 2a mod 2b = 2c;

    1
    2
    3
    4
    5
    6
    7
    8
    9
    10
    11
    12
    13
    14
    15
    16
    17
    18
    19
    20
    21
    22
    23
    24
    25
    26
    27
    28
    29
    30
    31
    32
    33
    34
    35
    36
    37
    38
    39
    40
    41
    42
    43
    44
    45
    46
    47
    48
    49
    50
    51
    52
    53
    54
    55
    56
    57
    58
    59
    60
    61

    ```c++
    //表达整数的奇怪方式
    #include <iostream>
    #include <algorithm>

    using namespace std;

    typedef long long LL;


    LL exgcd(LL a, LL b, LL &x, LL &y)
    {
        if (!b)
        {
            x = 1, y = 0;
            return a;
        }

        LL d = exgcd(b, a % b, y, x);
        y -= a / b * x;
        return d;
    }


    int main()
    {
        int n;
        cin >> n;

        LL x = 0, m1, a1;
        cin >> m1 >> a1;
        for (int i = 0; i < n - 1; i ++ )
        {
            LL m2, a2;
            cin >> m2 >> a2;
            LL k1, k2;
            LL d = exgcd(m1, m2, k1, k2);
            if ((a2 - a1) % d)
            {
                x = -1;
                break;
            }

            k1 *= (a2 - a1) / d;
            k1 = (k1 % (m2/d) + m2/d) % (m2/d);

            x = k1 * m1 + a1;

            LL m = abs(m1 / d * m2);
            a1 = k1 * m1 + a1;
            m1 = m;
        }

        if (x != -1) x = (a1 % m1 + m1) % m1;

        cout << x << endl;

        return 0;
    }

高斯消元

时间复杂度O(n^3)

  • 把某一行乘一个非零的数
    • 相当于等式两边同时乘以一个非零的数
  • 交换某两行
    • 相当于把方程组内的某两个方程交换一下位置
  • 把某行的若干倍加到另一行去
    • 方程相加
    • 初等行列变换
  • 完美阶梯型——唯一解
    • 0=非零——无解
    • 0=0——无穷解
  • 算法步骤
    • 枚举每一列
      • 找到这一列绝对值最大的一行
      • 将该行换到最上面去
      • 将该行第一个数变成1（同时除以常数）
      • 将下列所有行的第c列消成0
      • 重复上述过程
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
// a[N][N]是增广矩阵
int gauss()
{
    int c, r;
    for (c = 0, r = 0; c < n; c ++ )
    {
        int t = r;
        for (int i = r; i < n; i ++ )   // 找到绝对值最大的行
            if (fabs(a[i][c]) > fabs(a[t][c]))
                t = i;

        if (fabs(a[t][c]) < eps) continue;

        for (int i = c; i <= n; i ++ ) swap(a[t][i], a[r][i]);      // 将绝对值最大的行换到最顶端
        for (int i = n; i >= c; i -- ) a[r][i] /= a[r][c];      // 将当前行的首位变成1
        for (int i = r + 1; i < n; i ++ )       // 用当前行将下面所有的列消成0
            if (fabs(a[i][c]) > eps)
                for (int j = n; j >= c; j -- )
                    a[i][j] -= a[r][j] * a[i][c];

        r ++ ;
    }

    if (r < n)
    {
        for (int i = r; i < n; i ++ )
            if (fabs(a[i][n]) > eps)
                return 2; // 无解
        return 1; // 有无穷多组解
    }

    for (int i = n - 1; i >= 0; i -- )
        for (int j = i + 1; j < n; j ++ )
            a[i][n] -= a[i][j] * a[j][n];

    return 0; // 有唯一解
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
#include <iostream>
#include <algorithm>
#include <iomanip>
#include <cmath>
using namespace std;
int n ;
double a[N][N];
const int N=110;
const double eps=1e-6;
int gauss(){
    int c,r;//c是列，r是行
    for(c=0,r=0;c<n;c++){//枚举每一列
        int t = r;//记录行数
        for(int i =r;i<n;i++)//从该行开始枚举到最后一行
            	if(fabs(a[i][c])>fabs(a[t][c]))//对比绝对值，t为当前记录的最大绝对值的行，如果该行大于绝对值，更新记录的行序
                    t=i;
        if(fabs(a[t][c])<eps)//当前列是不是0 
            continue
        for(int i = c;i<=n;i++) swap(a[t][i],a[r][i]);//把绝对值最大的行的每一个数字换到最上面去
        for(int i = n;i>=c;i--) a[r][i]/=a[r][c];//枚举最大绝对值的行的每一个数，将该行第一个数变成1，同时同行改变其他数字
        for(int i = r+1;i<n;i++)
            if(fabs(a[i][c])>eps)//当前列不为0
                for(int j = n;j>=c;j--)
                    a[i][j] -=a[r][j]*a[i][c];//将下列所有行的第c列消成0
        r++;
    }
    if(r<n){
        for(int i = r;i<n;i++)
            if(fabs(a[i][n])>eps)
                return 2;//无解
        return 1;//有无穷多组解
    }
    for (int i = n - 1; i >= 0; i -- )//枚举行，从尾开始枚举
        for (int j = i + 1; j < n; j ++ )
            a[i][n] -= a[i][j] * a[j][n];//后面所有系数消成0
    return 0;//有唯一解
}
int main(){
cin >> n;
for(int i =0;i<n;i++)
	for(int j = 0;j<n;j++)
		cin >> a[i][j];
    int t =gauss();
    if(t == 0){
        for(int i = 0;i<n;i++) cout<<a[i][n];
    }
    else if(t==1)puts("Infinite group solutions");
    else puts("No solution");
    return 0;
}

组合计数

  • 有很多种方式，根据数据范围进行选择

Cnm=AnmAmm=n(n−1)(n−2)...(n−m+1)m!=n!m!(n−m)!C_n^m=\frac {A_n^m}{A_m^m}=\frac {n(n-1)(n-2)...(n-m+1)}{m!}=\frac {n!}{m!(n-m)!}

递推法求组合数

Cab=Ca−1b+Ca−1b−1C_a^b=C_{a-1}^b+C_{a-1}^{b-1}

10万组 1<=b<=a<=2000

时间复杂度n²

1
2
3
4
5
// c[a][b] 表示从a个苹果中选b个的方案数
for (int i = 0; i < N; i ++ )
    for (int j = 0; j <= i; j ++ )
        if (!j) c[i][j] = 1;
        else c[i][j] = (c[i - 1][j] + c[i - 1][j - 1]) % mod;

通过预处理逆元的方式求组合数

1万组 1<=b<=a<=10^5

Cab=fact[a]∗infact[b−a]∗infact[b]C_a^b=fact[a]*infact[b-a]*infact[b]

核心思想：预处理出来阶乘

时间复杂度nlogn

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
首先预处理出所有阶乘取模的余数fact[N]，以及所有阶乘取模的逆元infact[N]
如果取模的数是质数，可以用费马小定理求逆元
int qmi(int a, int k, int p)    // 快速幂模板
{
    int res = 1;
    while (k)
    {
        if (k & 1) res = (LL)res * a % p;
        a = (LL)a * a % p;
        k >>= 1;
    }
    return res;
}

// 预处理阶乘的余数和阶乘逆元的余数
fact[0] = infact[0] = 1;
for (int i = 1; i < N; i ++ )
{
    fact[i] = (LL)fact[i - 1] * i % mod;
    infact[i] = (LL)infact[i - 1] * qmi(i, mod - 2, mod) % mod;
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
#include <iostream>
#include <algorithm>
using namespace std;
typedef long long LL;
const int N =100010,mod = 1e9+7;
int fact[N],infact[N];
int qmi(int a, int k, int p)    // 快速幂模板
{
    int res = 1;
    while (k)
    {
        if (k & 1) res = (LL)res * a % p;
        a = (LL)a * a % p;
        k >>= 1;
    }
    return res;
}
int main(){
   fact[0]=infact[0]=1;
    for(int i = 1;i<N;i++){
        fact[i]=(LL)fact[i-1]*i%mod;
        infact[i]=(LL)infact[i-1]*qmi(i,mod-2,mod)%mod;
    }
    int n;
    cin >> n;
    while(n--){
        int a,b;
        cin >> a >> b;
        cout<<(LL)fact[a]*infact[b]%mod*infact[a-b]%mod;
    }
    return 0;
}

Lucas定理

时间复杂度：

plogN*Logp

20组 1<=b<=a<=10^18 1<=p<=10^5

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
若p是质数，则对于任意整数 1 <= m <= n，有：
    C(n, m) = C(n % p, m % p) * C(n / p, m / p) (mod p)

int qmi(int a, int k, int p)  // 快速幂模板
{
    int res = 1 % p;
    while (k)
    {
        if (k & 1) res = (LL)res * a % p;
        a = (LL)a * a % p;
        k >>= 1;
    }
    return res;
}

int C(int a, int b, int p)  // 通过定理求组合数C(a, b)
{
    if (a < b) return 0;

    LL x = 1, y = 1;  // x是分子，y是分母
    for (int i = a, j = 1; j <= b; i --, j ++ )
    {
        x = (LL)x * i % p;
        y = (LL) y * j % p;
    }

    return x * (LL)qmi(y, p - 2, p) % p;
}

int lucas(LL a, LL b, int p)
{
    if (a < p && b < p) return C(a, b, p);
    return (LL)C(a % p, b % p, p) * lucas(a / p, b / p, p) % p;
}
int main(){
    int n;
    cin>> n;
    while(n--){
        LL a,b;
        cin >> a >> b>>p;
        cout<<lucas(a,b)<<endl;
    }
    return 0;
}

分解质因数求组合数

大答案需要高精度

  • 分解质因数
  • 高精度乘法
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
当我们需要求出组合数的真实值，而非对某个数的余数时，分解质因数的方式比较好用：
    1. 筛法求出范围内的所有质数
    2. 通过 C(a, b) = a! / b! / (a - b)! 这个公式求出每个质因子的次数。 n! 中p的次数是 n / p + n / p^2 + n / p^3 + ...
    3. 用高精度乘法将所有质因子相乘
#include <iostream>
#include <algorithm>
    using namespace std;
const int N =5010;
int primes[N], cnt;     // 存储所有质数
int sum[N];     // 存储每个质数的次数
bool st[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[primes[j] * i] = true;
            if (i % primes[j] == 0) break;
        }
    }
}


int get(int n, int p)       // 求n！中的质数的指数
{
    int res = 0;
    while (n)
    {
        res += n / p;
        n /= p;
    }
    return res;
}


vector<int> mul(vector<int> a, int b)       // 高精度乘低精度模板
{
    vector<int> c;
    int t = 0;
    for (int i = 0; i < a.size(); i ++ )
    {
        t += a[i] * b;
        c.push_back(t % 10);
        t /= 10;
    }

    while (t)
    {
        c.push_back(t % 10);
        t /= 10;
    }

    return c;
}

get_primes(a);  // 预处理范围内的所有质数

for (int i = 0; i < cnt; i ++ )     // 求每个质因数的次数
{
    int p = primes[i];
    sum[i] = get(a, p) - get(b, p) - get(a - b, p);
}

int main(){
    int a,b;
    cin>>a>>b;
    get_primes(a);
    for(int i = 0;i<cnt;i++)//遍历质数
    {
        int p =primes[i];
        sum[i]=get(a)-get(b)-get(a-b);//求出质数的指数
    }
    vector<int> res;
    res.push_back(1);
    for(int i = 0;i<cnt;i++)// 用高精度乘法将所有质因子相乘，先遍历质数
        	for(int j =0;j<sum[i];j++)//遍历该质数的指数
                res = mul(res,primes[i]);
    for(int i = res.size()-1;i>=0;i--) cout<<res[i];//高精度乘法的输出，倒着存正着输出
    puts("");
    return 0;
}

卡特兰数

给定n个0和n个1，它们按照某种顺序排成长度为2n的序列，满足任意前缀中0的个数都不少于1的个数的序列的数量为： Cat(n) = C(2n, n) / (n + 1)

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
#include <iostream>
#include <algorithm>
using namespace std;
const int mod =1e9+7;
typedef long long LL;
int qmi(int a ,int b,int p){
   int res= 1;
   while(k){
   if(k&1) res = (LL) res*a %p;
       a= (LL)a* a %p;
       k>>=1;
   }
    return res;
    
}
int main(){
    int n;
    cin>>n;
    int a = 2*n,b=n;
    int res=1;
    for(int i =a;i>a-b;i--)res = (LL)res*i%mod;
    for(int i = 1;i<=b;i++) res= (LL) res*qmi(i,mod-2,mod)%mod;
    res = (LL)res*qmi(n+1,mod-2,mod)%mod;
    cout<<res<<endl;
    return 0;
}

容斥原理

时间复杂度2^n

从n个数中选任意多个数的方案数

Cn0+Cn1+...+Cnn=2nC^0_n+C_n^1+...+C_n^n=2^n

一般枚举所有的选法，用位运算来枚举

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
#include <iostream>
#include <algorithm>
using namespace std;
const int N =20;
int n,m;
int p[N];
typedef long long LL;
int main(){
	cin >> n>>m;
	for(int i =0;i<m;i++) cin>>p[i];
	int res = 0;    // 利用位运算来统计所有集合被选的状态 将i用二进制表示 其中1为当前位集合被选 0为没有选 

	for(int i = 1;i<1<<m;i++){// 这里的 1 << m 是等价于2 的 m次方 ， 有2的m次方个选择个数；
		int t =1,cnt=0;    //t表示所有容斥数的乘积，cnt用来表示有多少个集合 
		for(int j =0;j<m;j++)// m代表的是把i进行2进制移位操作，需要把i化为2进制 m位;
			if(i>>j&1){// 移位判断为0还是1;
				cnt ++;
				if((LL)t*p[j]>n){ // 更新t值;
					t=-1;
					break;
				}
                t*=p[j];
			}
			if(t!=-1)//说明刚好整除
			{
				if(cnt%2) res+=n/t;// 判断s是奇数还是偶数， 奇数加偶数减
				else res-=n/t;
			}
		
		
	}
	cout<<res<<endl;
	return 0;
}

简单博弈论

先手必胜状态：可以走到某一个必败状态

先手必败状态：走不到任何一个必败状态

NIM游戏

给定N堆物品，第i堆物品有Ai个。两名玩家轮流行动，每次可以任选一堆，取走任意多个物品，可把一堆取光，但不能不取。取走最后一件物品者获胜。两人都采取最优策略，问先手是否必胜。

我们把这种游戏称为NIM博弈。把游戏过程中面临的状态称为局面。整局游戏第一个行动的称为先手，第二个行动的称为后手。若在某一局面下无论采取何种行动，都会输掉游戏，则称该局面必败。
所谓采取最优策略是指，若在某一局面下存在某种行动，使得行动后对面面临必败局面，则优先采取该行动。同时，这样的局面被称为必胜。我们讨论的博弈问题一般都只考虑理想情况，即两人均无失误，都采取最优策略行动时游戏的结果。
NIM博弈不存在平局，只有先手必胜和先手必败两种情况。

定理： NIM博弈先手必胜，当且仅当 A1 ^ A2 ^ … ^ An != 0

NIM博弈先手必败：A1 ^ A2 ^ … ^ An = 0

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
#include <iostream>
#include <algorithm>
using namespace std;
int main(){
	int n;
	int res=0;
	cin >>n;
	while(n--){
		int x;
		cin >> x;
		res^=x;
	}
	if(res)puts("Yes");
	else puts("No");
	return 0;
}

游戏ICG

若一个游戏满足：

由两名玩家交替行动；
在游戏进程的任意时刻，可以执行的合法行动与轮到哪名玩家无关；
不能行动的玩家判负；
则称该游戏为一个公平组合游戏。
NIM博弈属于公平组合游戏，但城建的棋类游戏，比如围棋，就不是公平组合游戏。因为围棋交战双方分别只能落黑子和白子，胜负判定也比较复杂，不满足条件2和条件3。

有向图游戏

给定一个有向无环图，图中有一个唯一的起点，在起点上放有一枚棋子。两名玩家交替地把这枚棋子沿有向边进行移动，每次可以移动一步，无法移动者判负。该游戏被称为有向图游戏。
任何一个公平组合游戏都可以转化为有向图游戏。具体方法是，把每个局面看成图中的一个节点，并且从每个局面向沿着合法行动能够到达的下一个局面连有向边。

Mex运算

设S表示一个非负整数集合。定义mex(S)为求出不属于集合S的最小非负整数的运算，即：
mex(S) = min{x}, x属于自然数，且x不属于S

SG函数

在有向图游戏中，对于每个节点x，设从x出发共有k条有向边，分别到达节点y1, y2, …, yk，定义SG(x)为x的后继节点y1, y2, …, yk 的SG函数值构成的集合再执行mex(S)运算的结果，即：
SG(x) = mex({SG(y1), SG(y2), …, SG(yk)})
特别地，整个有向图游戏G的SG函数值被定义为有向图游戏起点s的SG函数值，即SG(G) = SG(s)。

SG（终点）=0

必败态是0

任何一种非零状态一定会有某种方式到0，任何一种0的方式是走不到0的

n个图的局面：

如果有多个独立的局面的话，可以分别求每个sg的值，然后异或起来，就是当前整个游戏sg的值

定理

有向图游戏的某个局面必胜，当且仅当该局面对应节点的SG函数值大于0。
有向图游戏的某个局面必败，当且仅当该局面对应节点的SG函数值等于0。

有向图游戏的和

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
//集合-nim
#include <iostream>
#include <algorithm>
#include <cstring>
#include <unordered_set>
using namespace std;
const int N=110,M=10010;
int n,m;
int s[N],f[N];//s存储操作数量，f存储sg函数
int sg(int x){
	if(f[x]!=-1) return f[x];//说明已经找出该点的sg函数，不重复计算
	
	unordered_set<int> S;//用来存储所有能够到的局面
	for(int i = 0;i<m;i++){
		int sum=s[i];//从中取出s[i]个石子
		if(x>=sum)S.insert(sg(x-sum));//如果取走石子的个数没有超过剩余石子的个数，说明可以取走，往记录图中加入该局面
		
	}
    //mex操作
	for(int i = 0;;i++)
		if(!S.count(i))//如果该局面没有出现过
			return f[x]=i;//说明该点的sg就是它本身
	
	
}
int main(){
	cin >>m;
	for(int i =0;i<m;i++)
		cin >>s[i];//读入每次可取操作数
		cin>>n;
		memset(f,-1,sizeof f);//初始化
		int res =0;
		for(int i = 0;i<n;i++){
			int x;
			cin>>x;//读入石子数量
			res^=sg(x);//有向图往下一个节点走
		}
		if(res) puts("Yes");
		else puts("No");
		return 0;
}
算法
#算法
算法·数学知识
http://example.com/2024/03/27/算法·数学知识/
作者
Kugeln
发布于
2024年3月27日
许可协议