十三、数论基础

模运算

模运算是大数运算中的常用操作。如果一个数太大,无法直接输出,或者不需要直接输出,可以把它取模后缩小数值再输出。定义取模运算为a除以m的余数,记为:
a    m o d    m   =   a   %   m a\ \ mod\ \ m\ =\ a\ \%\ m a  mod  m = a % m
取模的结果满足 0 ≤ a % m ≤ m − 1 0≤a\%m≤m-1 0a%mm1,题目用给定的 m m m限制计算结果的范围。例如 m = 10 m=10 m=10,就是取计算结果的个位数。取模操作满足以下性质:
加: ( a + b ) %   m = ( ( a %   m ) + ( b %   m ) ) %   m 加:(a+b)\%\ m=((a\%\ m)+(b\%\ m))\%\ m 加:(a+b)% m=((a% m)+(b% m))% m
减: ( a − b ) %   m = ( ( a %   m ) − ( b %   m ) ) %   m 减:(a-b)\%\ m=((a\%\ m)-(b\%\ m))\%\ m 减:(ab)% m=((a% m)(b% m))% m
乘: ( a   ×   b ) %   m = ( ( a %   m )   ×   ( b %   m ) ) %   m 乘:(a\ ×\ b)\%\ m=((a\%\ m)\ ×\ (b\%\ m))\%\ m 乘:(a × b)% m=((a% m) × (b% m))% m
然而,对除法取模进行下面的类似操作是错误的:
( a / b ) %   m = ( ( a %   m ) / ( b %   m ) ) % m (a/b)\%\ m=((a\%\ m)/(b\%\ m))\%m (a/b)% m=((a% m)/(b% m))%m
如: ( 100 / 50 ) % 20 = 2 , ( 100 % 20 ) / ( 50 % 20 ) % 20 = 0 (100/50)\%20=2,(100\%20)/(50\%20)\%20=0 (100/50)%20=2(100%20)/(50%20)%20=0,两者不相等。除法的取模需要用到逆元。

快速幂

快速幂以及扩展的矩阵快速幂,由于应用场景比较常见,也是竞赛中常见的题型。幂运算 a n a^n an n n n a a a相乘。快速幂就是高效地算出 a n a^n an。当 n n n很大时,例如: n = 1 0 9 n=10^9 n=109,计算 a n a^n an这样大的数,一是数字太大,二是计算时间很长。
下面先考虑如何缩短计算时间,如果用暴力的方法直接算 a n a^n an,即逐个做乘法,复杂度是 O ( n ) O(n) O(n),即使能算出来,也会超时。
读者很容易想到快速幂的办法:先算 a 2 a² a2,然后继续算平方 ( a 2 ) 2 (a²)² (a2)2,一直算到 n n n次幂。这
是分治法的思想,复杂度为 O ( l o g 2 n ) O(log₂n) O(log2n)。下面是代码,请读者自己理解:

int fastPow(int a, int n) {
    if (n == 1) return a;
    int temp = fastPow(a, n / 2);   //分治
    if (n % 2 == 1) {               //奇数个a,此处也可以写为if(n &1)
        return temp * temp * a;
    } 
    else {                          //偶数个a
        return temp * temp;
    }
}

程序中的递归,层数只有 l o g 2 n log_2n log2n,不用担心溢出的问题。

矩阵快速幂

给定一个 m × m m×m m×m的矩阵 A A A,求它的 n n n次幂 A n A^n An,这也是常见的计算。同样有矩阵快速幂的算法,原理是把矩阵当作变量来操作,程序和上面的很相似。
首先需要定义矩阵的结构体,并且定义矩阵相乘的操作。注意矩阵相乘也需要取模:

const int MAXN = 2;                     //根据题目要求定义矩阵的阶,本例中是2
const int MOD = 1e4;                    //根据题目要求定义模
struct Matrix {                         //定义矩阵
    int m[MAXN][MAXN];
    Matrix() {                          //初始化矩阵
        memset(m, 0, sizeof(m));
    }
};

Matrix Multi(Matrix a, Matrix b) {      //矩阵的乘法
    Matrix res;
    for (int i = 0; i < MAXN; i++)
        for (int j = 0; j < MAXN; j++)
            for (int k = 0; k < MAXN; k++)
                res.m[i][j] = (res.m[i][j] + a.m[i][k] * b.m[k][j]) % MOD;
    return res;
}

下面是矩阵快速幂的程序代码,和前面单变量的快速幂的代码非常相似:

Matrix fastm(Matrix a, int n) {
    Matrix res;
    for (int i = 0; i < MAXN; i++) res.m[i][i] = 1;
			//初始化为单位矩阵,相当于前面程序中的"int res =1;"
    while (n) {
        if (n & 1)res = Multi(res, a);
        a = Multi(a, a);
        n /= 2;
    }
    return res;
}

最大公约数,最小公倍数

最大公约数GCD和最小公倍数LCM是竞赛中常见的知识点,虽然这两个知识点很容易理解,但往往会与其他知识点结合起来出综合题。

int GCD(int a, int b){
    if(a % b == 0)return b;
    return GCD(b, a % b);
}

int LCM(int a, int b){
    return a * b / GCD(a, b);
}

扩展欧几里得算法与二元一次方程的整数解

给出整数 a 、 b 、 n a、b、n abn,问方程 a x + b y = n ax+by=n ax+by=n什么时候有整数解?如何求所有的整数解?
有解的充分必要条件是: g c d ( a , b ) gcd(a,b) gcd(a,b)可以整除 n n n。简单证明如下:

令:a = gcd(a, b)a'、b = gcd(a, b)b';
有:ax + by = gcd(a, b) (a'x + b'y) = n;
如果x、y、a'、b'都是整数,
那么:n必须是gcd(a,b)的倍数才有整数解。

例如 4 x + 6 y = 8 4x+6y=8 4x+6y=8 2 x + 3 y = 4 2x+3y=4 2x+3y=4有整数解, 4 x + 6 y = 7 4x+6y=7 4x+6y=7则没有整数解。如果确定有解,一种解题方法是先找到一个解 ( x 0 , y 0 ) (x_0,y_0) (x0,y0),那么通解公式如下:
{ x = x 0 + b t y = y 0 − a t     t 为任意整数 \begin{cases} x = x_0 + bt\\ y = y_0 -at\ \ \ t为任意整数 \end{cases} {x=x0+bty=y0at   t为任意整数
所以,问题转化为如何求 ( x 0 , y 0 ) (x_0,y_0) (x0,y0)。利用扩展欧几里得算法可以求出这个特解。

扩展欧几里得算法

当方程符合 a x + b y = g c d ( a , b ) ax+by=gcd(a,b) ax+by=gcd(a,b)时,可以用扩展欧几里得算法求 ( x 0 , y 0 ) (x_0,y_0) (x0,y0)。程序如下:

void ex_gcd(int a, int b, int &x, int &y) { //x, y变化后形成特解(x0, y0)
    if (b == 0) {
        x = 1, y = 0;
        return;
    }
    ex_gcd(b, a % b, x, y);
    int tmp = x;
    x = y;
    y = tmp - (a / b) * y;
}

有时候为了简化描述,在 a x + b y = g c d ( a , b ) ax+by=gcd(a,b) ax+by=gcd(a,b)两边除以 g c d ( a , b ) gcd(a,b) gcd(a,b),得到 c x + d y = 1 cx+dy=1 cx+dy=1
其中: c = a / g c d ( c , b ) , d = b / g c d ( a , b ) c=a/gcd(c,b),d=b/gcd(a,b) c=a/gcd(c,b),d=b/gcd(a,b)
很明显, c 、 d c、d cd是互质的。 c x + d y = 1 cx+dy=1 cx+dy=1的通解如下:
{ x = x 0 + d t y = y 0 − c t     t 为任意整数 \begin{cases} x = x_0 + dt\\ y = y_0 -ct\ \ \ t为任意整数 \end{cases} {x=x0+dty=y0ct   t为任意整数

求任意方程ax+by=n的一个整数解

用扩展欧几里得算法求解 a x + b y = g c d ( a , b ) ax+by=gcd(a,b) ax+by=gcd(a,b)后,利用它可以进一步解任意方程 a x + b y = n ax+by=n ax+by=n,得到一个整数解。其步骤如下:

  1. 判断方程 a x + b y = n ax+by=n ax+by=n是否有整数解,有解的条件是 g c d ( a , b ) gcd(a,b) gcd(a,b)可以整除 n n n
  2. 用扩展欧几里得算法求 a x + b y = g c d ( a , b ) ax+by=gcd(a,b) ax+by=gcd(a,b)的一个解 ( x 0 , y 0 ) (x_0,y_0) (x0,y0)
  3. a x 0 + b y 0 = g c d ( a , b ) ax_0+by_0=gcd(a,b) ax0+by0=gcd(a,b)两边同时乘以 n / g c d ( a , b ) n/gcd(a,b) n/gcd(a,b)得: a x 0 n / g c d ( a , b ) + b y 0 n / g c d ( a , b ) = n ax_0n/gcd(a,b)+by_0n/gcd(a,b)=n ax0n/gcd(a,b)+by0n/gcd(a,b)=n
  4. 对照 a x + b y = n ax+by=n ax+by=n,得到它的一个解 ( x ‘ , y ‘ ) (x`,y`) (x,y)是:
    { x ‘ = x 0 ∗ n / g c d ( a , b ) y ‘ = y 0 ∗ n / g c d ( a , b ) \begin{cases} x` = x_0*n/gcd(a,b)\\ y` = y_0*n/gcd(a,b) \end{cases} {x=x0n/gcd(a,b)y=y0n/gcd(a,b)
应用场合

扩展欧几里得算法是一个很有用的工具,在竞赛题目中常用于以下场合:

  1. 求解不定方程
  2. 求解模的逆元
  3. 求解同余方程
    虽然用扩展欧几里得算法可以算 a x + b y = g c d ( a , b ) ax+by=gcd(a,b) ax+by=gcd(a,b)的通解,不过一般没有这个需求,而是用于求某些特殊的解,例如求解逆元,逆元是除法取模操作常用的工具。

同余与逆元

同余在数论中非常有用,它用类似处理等式的方式来处理整除关系,非常简便。

同余

两个整数 a 、 b a、b ab和一个正整数 m m m,如果 a a a除以 m m m所得的余数和 b b b除以 m m m所得的余数相等。
即: a % m = b % m a\%m=b\%m a%m=b%m,称 a a a b b b对于 m m m同余, m m m称为同余的模。
同余的概念也可以这样理解: m ∣ ( a − b ) m|(a-b) m(ab)
即: a − b a-b ab m m m的整数倍。例如 6 ∣ ( 23 − 5 ) 6|(23-5) 6∣(235) 23 23 23 5 5 5对模 6 6 6同余。
同余的符号记为 a ≡ b ( m o d   m ) a\equiv b(mod\ m) ab(mod m)

一元线性同余方程

a x ≡ b ( m o d   m ) ax\equiv b(mod\ m) axb(mod m),即 a x ax ax除以 m , b m,b m,b除以 m m m的余数相同,这里 a 、 b 、 m a、b、m abm都是整数,求解 x x x的值。
方程也可以这样理解: a x − b ax-b axb m m m的整数倍。设 y y y是倍数,那么 a x − b = m y ax-b=my axb=my
移项得到: a x − m y = b ax-my=b axmy=b。因为 y y y可以是负数,改写为: a x + m y = b ax+my=b ax+my=b
这就是在扩展欧几里得算法中的二元一次不定方程
当且仅当 g c d ( a , m ) gcd(a,m) gcd(a,m)能整除 b b b时有整数解。例如 15 x + 6 y = 9 15x+6y=9 15x+6y=9,有整数解 x = 1 , y = − 1 x=1,y=-1 x=1,y=1。当 g c d ( a , m ) = b gcd(a,m)=b gcd(a,m)=b时,可以直接用扩展欧几里得算法求解 a x + m y = b ax+my=b ax+my=b
如果不满足 g c d ( a , m ) = b gcd(a,m)=b gcd(a,m)=b,还能用扩展欧几里得算法求解 a x + m y = b ax+my=b ax+my=b吗?答案是肯定的,但是需要结合逆元。

逆元

给出 a a a m m m,求解方程 a x ≡ 1 ( m o d   m ) ax\equiv 1(mod\ m) ax1(mod m),即 a x ax ax除以 m m m余数是 1 1 1
根据前面的讨论,有解的条件是 g c d ( a , m ) = 1 gcd(a,m)=1 gcd(a,m)=1,即 a 、 m a、m am互质。该问题等价于求解$ax+my=1,可以用上一节的扩展欧几里得算法求解。例如 8 x ≡ 1 ( m o d   31 ) 8x\equiv 1(mod\ 31) 8x1(mod 31),等价于求解 8 x + 31 y = 1 8x+31y=1 8x+31y=1,用扩展欧几里得算法求得一个特解是 x = 4 , y = − 1 x=4,y=-1 x=4,y=1
方程 a x ≡ 1 ( m o d   m ) ax\equiv 1(mod\ m) ax1(mod m)的一个解 x x x,称 x x x a a a m m m。注意,这样的 x x x有很多,把它们都称为
求逆元的代码如下:

	int inverse_mod(int a, int m){  
		int x, y;  
		ex_gcd(a, m, x, y);  
		return(m + x % m ) % m;  
}

试除法判断质数

问题:输入一个很大的数 n n n,判断它是不是素数。
素数定义:一个数 n n n,如果不能被 [ 2 , n − 1 ] [2,n-1] [2,n1]内的所有数整除, n n n就是素数。当然,并不需
要把 [ 2 , n − 1 ] [2,n-1] [2,n1]内的数都试一遍,这个范围可以缩小到 [ 2 , n ] [2,\sqrt n] [2,n ]
给定 n n n,如果它不能整除 [ 2 , n ] [2,\sqrt n] [2,n ]内的所有数,它就是素数。证明如下:

n = a × b n=a×b n=a×b,有 m i n ( a , b ) ≤ n min(a,b)≤\sqrt n min(a,b)n ,令 a ≤ b a≤b ab。只要检查 [ 2 , n ] [2,\sqrt n] [2,n ]内的数,如果 n n n不是素数,就能找到一个 a a a。如果不存在这个 a a a,那么 [ n , n − 1 ] [\sqrt n,n-1] [n ,n1]内也不存在 b b b

以上判断的范围可以再缩小一点: [ 2 , n ] [2,\sqrt n] [2,n ]内所有的素数。其原理很简单,读者在学过下文提到的埃式筛法之后更容易理解。
用试除法判断素数,复杂度是 O ( n ) O(\sqrt n) O(n ),对于 n ≤ 1 0 12 n≤10^{12} n1012的数是没有问题的。下面是试除法判断素数的代码:

bool is_pri(long long n){  
	if(n <= 1)return false;  
	for(int i = 2; i * i <= n; i++){  
		if(n % i == 0)return false;  
	}  
	return ture;  
}

相关推荐

  1. 数论基础

    2024-06-19 09:36:04       13 阅读
  2. 算法基础

    2024-06-19 09:36:04       36 阅读
  3. 【算法集训】基础数据结构:、哈希表

    2024-06-19 09:36:04       42 阅读
  4. 贪心算法基础题(第天)

    2024-06-19 09:36:04       16 阅读
  5. GO基础进阶篇 ()、泛型

    2024-06-19 09:36:04       35 阅读
  6. 网络工程师:新兴科技基础知识面试题(

    2024-06-19 09:36:04       37 阅读

最近更新

  1. python之参数注解介绍

    2024-06-19 09:36:04       0 阅读
  2. 学习opencv

    2024-06-19 09:36:04       0 阅读
  3. DP学习——中介者模式

    2024-06-19 09:36:04       0 阅读
  4. 交换机(Switches)和桥(Bridges)的区别

    2024-06-19 09:36:04       0 阅读
  5. 测试面试宝典(二十一)—— get和post的区别

    2024-06-19 09:36:04       0 阅读
  6. ESP8266AT指令查看有哪些指令可用(3)

    2024-06-19 09:36:04       0 阅读

热门阅读

  1. Ruby 数据库访问 - DBI 教程

    2024-06-19 09:36:04       17 阅读
  2. 安卓交叉编译——ndk

    2024-06-19 09:36:04       16 阅读
  3. Swarm 集群管理

    2024-06-19 09:36:04       14 阅读
  4. PostgreSQL源码分析——创建用户

    2024-06-19 09:36:04       14 阅读
  5. Linux 上的 TTY 是什么?

    2024-06-19 09:36:04       16 阅读
  6. USB 端点停止

    2024-06-19 09:36:04       14 阅读
  7. 通信基础知识

    2024-06-19 09:36:04       14 阅读
  8. 如何给vue开发的网站做seo?

    2024-06-19 09:36:04       16 阅读