数学知识(上)

质数

AcWing 866. 试除法判定质数

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
#include <iostream>
#include <algorithm>

using namespace std;

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;
}

int main()
{
    int n;
    cin >> n;
    
    while (n -- )
    {
        int x;
        cin >> x;
        if (is_prime(x)) puts("Yes");
        else puts("No");
    }
    
    return 0;
}

AcWing 867. 分解质因数

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
#include <bits/stdc++.h>

using namespace std;

void divide(int x)
{
    for (int i = 2; i <= x / i; i ++ )
        	if (x % i == 0)
            {
                int s = 0;
                while (x % i == 0) x /= i,s ++ ;
                cout << i << ' ' << s << endl;
            }
   	if (x > 1) cout << x << ' ' << 1 << endl;
    cout << endl;
}

int main()
{
    int n;
    cin >> n;
    while (n -- )
    {
        int x;
        cin >> x;
        divide(x);
    }
    return 0;
}

AcWing 868. 筛质数—朴素

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
#include <bits/stdc++.h>

using namespace std;

const int N = 1000010;

int primes[N], cnt;
bool st[N];

void get_primes(int n)
{
    for (int i = 2; i <= n; i ++ )
    {
        if (st[i]) continue;
        primes[cnt ++ ] = i;
        for (int j = i + i; j <= n; j += i)
            st[j] = true;
    }
}

int main()
{
    int n;
    cin >> n;
    
    get_primes(n);
    cout << cnt << endl;
    
    return 0;
}

AcWing 868. 筛质数—线性筛法

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
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
#include <bits/stdc++>

using namespace std;

const int N = 1000010;

int primes[N], cnt;  // primes[] 保存质数， cnt 是质数个数
bool st[N];   // st[i] = true  表示 i 不是质数

void get_primes(int n)
{
    for (int i = 2; i <= n; i ++ )
    {
        if (!st[i]) primes[cnt ++ ] = i;  // 如果 i 未被标记， i 是质数
        // i是质数就筛掉i与primes中所有质数的乘积；i是非质数筛掉i与primes中所有<=最小质因子的乘积
        for (int j = 0; primes[j] <= n / i; j ++ )
        {
            st[primes[j] * i] = true;  // 标记 primes[j] * i 为合数
            if (i % primes[j] == 0) break;  // i 含有 primes[j] 这个质因子， 终止
            // 如果i % primes[j] == 0   
            //	那么primes[j]一定是primes[j] * i的最小质因子
            // 如果i % primes[j] != 0   
            // 	那么primes[j]一定小于i的所有质因子 primes[j]也一定是primes[j] * i的最小质因子
           
            /*
            	举个例子i=12， 首先 st[24] = true；
				如果没有i%primes[j]的判断，那么会执行st[36] = true;
				当i等于18的时候 st[i * primes[0]] = true;  这就是标记两次
				思想就是，每次都用最小质数来标记，这也就不会出现重复标记的情况
			*/
        }
    }
}

int main()
{
    int n;
    cin >> n;
    
    get_primes(n);
    
    cout << cnt << endl;
    
    return 0;
}

/*
 第一次外层循环：i = 2

- primes[] 为空，st[] 所有元素为 `false`。

标记：st[2 * 2] = st[4] = true;`

- 标记 4 为合数。

检查：if (2 % primes[0] == 0) break;

- 2 % 2 == 0，跳出内层循环，表示 2 是最小质因子，后续无需检查更大的质数。

此时：

- primes[] = {2}
- st[] = {false, false, false, true, false, false, false, false, false, false, false}
- 合数 4 被标记为合数，跳出内层循环。

------

2. 第二次外层循环：i = 3

标记：st[3 * 2] = st[6] = true;

- 标记 6 为合数。

检查：if (3 % primes[0] == 0) break;

- 3 % 2 != 0，继续执行下一个质数。

标记：st[3 * 3] = st[9] = true;

- 标记 9 为合数。

检查：if (3 % primes[1] == 0) break;

- 3 % 3 == 0，跳出内层循环，表示 3 是 9 的最小质因子。

此时：

- primes[] = {2, 3}
- st[] = {false, false, false, true, false, false, true, false, false, true, false}
- 合数 6 和 9 被标记，跳出内层循环。

------

3. 第三次外层循环：i = 4

- 4 已经被标记为合数了，st[4] == true，跳过。

------

4. 第四次外层循环：i = 5

标记：st[5 * 2] = st[10] = true;

- 标记 10 为合数。

检查：if (5 % primes[0] == 0) break;

- 5 % 2 != 0，继续执行下一个质数。

标记：st[5 * 3] = st[15]，但 15 超出了范围。

此时：

- primes[] = {2, 3, 5}
- st[] = {false, false, false, true, false, false, true, false, false, true, true}
- 合数 10 被标记。
*/

快速幂

AcWing 875. 快速幂

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
#include <bits/stdc++.h>

using namespace std;

typedef long long LL;

LL qmi(int a, int b, int p)
{
    LL res = 1;
    while (b)
    {
        if (b & 1) res = (LL)res * a % p;  // 这里有个定理
        b >>= 1;
        a = (LL)a * a % p;  // 反复平方
    }
    return res;
}

int main()
{
    int n;
    scanf("%d", &n);
    while (n -- )
    {
        int a, b, p;
        scanf("%d%d%d", &a, &b, &p);
        printf("%lld", qmi(a, b, p));
    }
    return 0;
}

AcWing 876. 快速幂求逆元

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 / b = a * x(mod p)     a / b = a * b ^ -1    p是质数
// b * x = 1 (mod p)   求b * b ^ -1 = 1(mod p)
// 那么  x叫做  b 的 乘法逆元
// 费马定理 b^(p - 1) = 1(mod p)   b * b^(p - 2) = 1(mod p);这时候直接求b^(p - 2)即可
// 那么 b ^ -1 = b^(p - 2)     b * b^(p - 2) = 1(mod p) --->b ^ -1 = b^(p - 2)(mod p)

#include <bits/stdc++.h>

using namespace std;

typedef long long LL;

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

int main()
{
    int n;
    scanf("%d", &n);
    while (n -- )
    {
        int a, p;
        scanf("%d%d", &a, &p);
        if (a % p == 0) puts("impossible");
        else printf("%lld\n", qmi(a, p - 2, p));
    }
    return 0;
}

扩展欧几里得算法

AcWing 877. 扩展欧几里得算法

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 <bits/stdc++.h>

using namespace std;

// 根据费马定理，任意正整数a, b都存在x, y使得ax + by = gcd(a, b)
int exgcd(int a,int b, int &x, int &y)
{
    if (!b)  // 如果b = 0, 则gcd(a, b) = 1 * a + 0 * b
    {
        x = 1, y = 0;
        return a;
    }
    int d = exgcd(b, a % b, y, x);
    y -= a / b * x;
    return d;
}

int main()
{
    int n;
    scanf("%d", &n);
    
    while (n -- )
    {
        int a, b;
        scanf("%d%d", &a, &b);
        int x, y;
        exgcd(a, b, x, y);
        printf("%d %d\n", x, y);
    }
    return 0;
}

证明

AcWing 878. 线性同余方程

1
2
3
4
5
6
7
8
a * x ≡ b(mod m)   这个公式的意思(a * x) % m = b

a * x = m * y + b   --->  a * x - m * y = b
令 y' = -y    ----> a * x + m * y' = b
这里就是扩展欧几里得算法一样了 ，两边同时乘以 d / b 倍得到
a * x * d / b + m * y' * d / b = d
令 x0 = x * d / b， y0
后面x ≡ x0 * (b / d) (mod m)
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
#include <bits/stdc++.h>

using namespace std;

typedef long long LL;

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 main()
{
    int n;
    scanf("%d", &n);
    
    while (n -- )
    {
        int a, b, m;
        scanf("%d%d%d", &a, &b, &m);
        
        int x, y;
        int d = exgcd(a, m, x, y);
        if (b % d) puts("impossible");
        else printf("%d\n", (LL)b / d * x % m);
    }
    return 0;
}

AcWing 204. 表达整数的奇怪方式 （中国剩余定理）

x ≡ m1 (mod a1) 这个就是说x和m1对a1取模，结果是一样的，所以可以用a1和m1来表示x

—> x = k1 * a1 + m1

我们要求解如下形式的线性同余方程组

1
2
3
4
5
6
7
8
9
10
11
12
13
x ≡ m₁ (mod a₁)     x = k1 * a1 + m1
x ≡ m₂ (mod a₂)  	x = k2 * a2 + m2
...  
x ≡ mₙ (mod aₙ)		 x = kn * an + mn

k1 * a1 + m1 = k2 * a2 + m2
k1 * a1 - k2 * a2 = m2 - m1    求k1  k2   使用扩展欧几里得算法就能求了 
如果有解  等价于(a1, a2) | m2 - m1        x = k1 * a1 + m1

===>  k1 + k * a2 / d
	  k2 + k * a2 / d   这就是方程的所有解
	  x = k1 * a1 + m1 = (k1 + k * a2 / d) * a1 + m1
	    = a1*k1+m1+k[a1,a2] = x0 + k * a  ==>   x mod a ≡ x0  ===> x0 mod a

举个例子来看

1
2
3
4
31 mod 8 = 7
31 mod 11 = 9
验证  31 ÷ 8 = 3 余 7 ⇒ 31 ≡ 7 (mod 8)
	31 ÷ 11 = 2 余 9 ⇒ 31 ≡ 9 (mod 11)
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
#include <bits/stdc++.h>

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;
    
    // 读入方程个数 n，再读入第一个方程的模 a1 和余数 m1。
    bool has_answer = true;
    LL m1, a1;
    cin >> a1 >> m1;
    // 接下来逐个合并其余的 n-1 个方程。
    for (int i = 0; i < n - 1; i ++ )
    {
        LL m2, a2;
        cin >> a2 >> m2;
        LL k1, k2;
        // a1 * k1 + a2 * k2 = d = gcd(a1, a2)
        LL d = exgcd(a1, a2, k1, k2);
        
        if ((m2 - m1) % d){
            has_answer = false;
            break;
        }
        
        // 这里求的是k1*a1-k2*a2=d 我们要求的是 m2 - m1 所以都乘以个(m2 - m1) / d即可
        k1 *= (m2 - m1) / d; // 计算出通解中的 k1，化为最小非负解。
        LL t = a2 / d;   // k1+k*a2/d是一组解
        k1 = (k1 % t + t) % t;    // 这里就是变成最小的整数解
        
        m1 = k1 * a1 + m1;
        // a1 应该更新为 a1 和 a2的最大公约数
        a1 = abs(a1 / d * a2);  // 取绝对值更安全
    }
    
    // 得到的可能是负数，那么加上a1  再模就为正数了
    if (has_answer) cout << (m1 % a1 + a1) % a1 << endl;
    else puts("-1");
    
    return 0;
}