小指针

V1

2022/03/18阅读:45主题:凝夜紫

欧几里德算法及其扩展

最大公约数

若自然数d,同时是自然数ab的约数,我们称dab的公约数。而所有公约数中最大的,我们称为ab的最大公约数,记为gcd(a,b)

求最大公约数的算法有很多,比如:分解质因数、九章算术里的更相减损术、短除法等等。但是应用最广泛的还是今天的主角--辗转相除法。

辗转相除法

古希腊数学家欧几里得在其著作《The Elements》中最早描述了这种算法,所以又被命名为欧几里得算法。

辗转相除法依赖除法的一个基本性质(下面称为性质):

假设d|a,同时d|b,那么一定有d|(ax+by)

其中|表示整除的意思,其中d|a,表示d可以整除a,既,ad的倍数,比如3可以整除9,表示为3|9

所以,这个性质的意思是,如果d能同时整除a和b,那么d也能整除a的倍数加上b的倍数

例如:4能整除8,也能整除12,那么4也能整除8*3+12*5(其中3和5可以取任意不为0整数)。

通过这个性质,我们就能得到辗转相除法的最终公式:gcd(a, b)=gcd(b, a mod b)

证明推导

  1. a<b的时候,有a%b==a,那么最终公式自然成立。

  2. a>=b的时候。首先看模数部分,先把模运算改成普通的四则运算: (其中 用来表示a除以b并取整)。
    举例:

  3. 在上面的公式中,因为 一定是一个整数,我们用x来表示这个整数,那么,公式变成

  4. 把这个公式带入到最终结果中得到:

    如何证明两面是相等的呢?

    1. 先看左面,假设dab的任意公约数,那么显然d可以整除ab,再根据上面的除法性质,我们知道d也能整除a-x*b(其中a的整数倍是1b的整数倍是-x)。又因为d是任意公约数,那么d也可以是最大公约数,所以gcd(a,b)=d,gcd(b,a-x*b)=d,推出gcd(a,b)==>gcd(b, a-x*b)
    2. 再看右面,假设dba-x*b的任意公约数,那么显然d可以整除b,那么根据上面的除法性质,我们知道d也能整除a-x*b+x*b(其中a-x*b的整数倍是1b的整数倍是x),也就是d能整除a。同步骤1的推论,可以得出gcd(b, a-x*b)==>gcd(a,b)
    3. 综合步骤1和步骤2,得出最终结论gcd(a, b)=gcd(b, a mod b)

Code

int gcd(int a, int b) {
    return b ? gcd(b, a % b) : a;
}

扩展欧几里德算法

裴蜀定理:关于任意的正整数a,b,一定有非零数对x,y,可以让ax+ay=gcd(a, b)

扩展欧几里德算法的作用就是求出其中的系数x,y

证明推导

  1. 在上一节,我们已经知道了欧几里德算法的最终形态:gcd(a, b)=gcd(b, a mod b)

  2. 根据裴蜀定理,现在有公式ax+by=gcd(a, b)。联立两个公式可以推出:ax+by=gcd(b, a mod b)

  3. 根据裴蜀定理,从gcd(b, a mod b)又可以推出

  4. 代入 到上面的公式,并且联合最初的 ,从而得到

  5. 我们看一下其中 的关系。将上一步的公式展开:

    由此,我们得到了 的关系:其中 原封不动的变成了 ,而 变成了 。该规则说明最终的 ,一定能通过子问题 构造出来。

  6. 那么最初的 是什么呢?我们注意到在欧几里德算法中,辗转相除到最后的子问题一定是gcd(a, 0),那么对于等式 ,如果让 一定成立。且由于b此时是0,所以 可以是任意整数,为了方便我们就让它也是0。

Code

输入示例:
2
4 6
8 18
输出示例:
-1 1
-2 1

在代码实现上,在辗转相除计算最大公约数的基础上,同时把系数x,y也求出来。其中有两点需要注意一下:

  1. x,y必须传递引用,否则在递归的过程中,x,y都不相同的。如exgcd(6,4)中的x,ygcd(4,2)中的x,ygcd(2,0)中初始化出来的x=1,y=0都是不相同的,而我们最后要得到初始化的那个x=1,y=0中的x,y,所以在传参的时候要传引用。
  2. 因为求公约数是自顶向上的递归,而求系数x,y的时候是从1,0求到x,y,因此必须使用自底向上的递归,也就是把求y的部分放到求公约数之后。
#include <iostream>
using namespace std;

int exgcd(int a, int b, int &x, int &y) {
    if (!b) {
        x = 1, y = 0// 初始化x,y
        return a;
    }
    // x和y也要换位置,因为a
    int d = exgcd(b, a % b, y, x); 
    y -= a / b * x; // 自底向上的递归求y。
    return 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 << endl;
    }
}

对于没办法传递引用的语言(“对,这里说的就是你,python”),我们需要手动翻转x和y的值。

def exgcd(a, b):    
    if b == 0:          
        return a, 10
    else:         
        d, x, y = exgcd(b, a % b)
        # 此处手动反转x,y
        x, y = y, x - (a // b) * y
        return d, x, y

分类:

数学

标签:

数学编程

作者介绍

小指针
V1