浅谈扩展欧几里得(Exgcd)

前言

假设我们已经会了欧几里得算法
而且,真真真真的是浅谈

基本形式

[ax+by=gcd(a,b) ]

[a, bin mathbb{N}^* ]

扩展欧几里得 (Exgcd) 则是求解以上方程的整数解

求特解

观察基本形式 (ax+by=gcd(a,b))
并考虑欧几里得算法的结束条件 (a=gcd, b=0)
则最终状态是 (x=1, y=0)
继续分析条件

[ecause gcd(a, b)=gcd(b, a \% b) ]

[ herefore bx'+(a\% b)y'=gcd(b, a\% b) ]

假设我们已经求解了

[bx'+(a\% b)y'=gcd(b, a\% b) ]

并且要求出他的上一状态

[ax+by=gcd(a,b) ]

我们知道

[a\% b=a-left lfloor frac{a}{b} ight floor b ]

所以

[egin{array}{lcl} egin{aligned} gcd(a, b) & = bx'+(a-left lfloor frac{a}{b} ight floor b)y' \ & = bx'+ay'-left lfloor frac{a}{b} ight floor by' \ & = ay'+(x'-left lfloor frac{a}{b} ight floor y')b \ end{aligned} \ ecause ax+by=gcd(a, b) \ herefore x=y', y=x'-left lfloor frac{a}{b} ight floor y' end{array} ]

最后递归求解即可

根据特解求通解

以下式为例,并假设有解,即 (gcd(a, b) | c)

[ax+by=c ]

我们已知一组 (ax+by=gcd(a,b)) 的特解,记为 (x_0, y_0)
则有

[egin{array}{rcl} ax_0+by_0&=&gcd(a,b)\ afrac{x_0c}{gcd(a,b)}+bfrac{y_0c}{gcd(a,b)}&=&c\ end{array} ]

于是我们找到了原方程的一组特解,记为 (x_1)(y_1)

[x_1=frac{x_0c}{gcd(a,b)}, y_1=frac{y_0c}{gcd(a,b)} ]

于是我们设任意 (din mathbb{Q}),那么一定有

[a(x_1+db)+b(y_1-da)=c ]

(x_1+db)(y_1-da) 必须为整数
又因为 (x_1, y_1) 为整数,所以我们只需要保证 (db, dain mathbb{Z})
令当 (d) 取到最小的可能的正值时的 (d_x=db, d_y=da),那么任意解中这两个变量与 (x_1, y_1) 的偏差一定分别为 (d_x)(d_y) 的倍数
那么显然最小 (d=frac{1}{gcd(a, b)}) ,所以 (d_x=frac{b}{gcd(a, b)}, d_y=frac{a}{gcd(a, b)})
于是通解形式为

[egin{array}{rcl} x&=&x_1+kd_x\ y&=&y_1-kd_y end{array} ]

其中,(k) 是任意整数

例题

P5656 【模板】二元一次不定方程 (exgcd)
以上我们已经求出通解形式,于是分类讨论即可(注意细节)

#include <algorithm>
#include <iostream>
#include <cstring>
#include <cstdio>
#include <cmath>
#include <queue>
#include <stack>
#include <list>
#include <map>
#include <set>
using namespace std;
#define fore(i,now) for(reg int i=head[now];i;i=e[i].nxt)
#define fora(i,a,b,c) for(reg int i=a;i<=b;i+=c)
#define forb(i,a,b,c) for(reg int i=a;i>=b;i-=c)
#define uLL unsigned LL
#define INF 0x3f3f3f3f
#define LL long long
#define reg register
#define R read()
inline LL read(){
    LL s=0, f=1; char c=getchar();
    while(c<'0'||c>'9') { if(c=='-') f=-f; c=getchar(); }
    while(c>='0'&&c<='9') s=(s<<3)+(s<<1)+(c^48), c=getchar();
    return s*f;
}
LL Exgcd(LL a, LL b, LL &x, LL &y){
    if(b==0) { x=1, y=0; return a; }
    LL g=Exgcd(b, a%b, x, y);
    LL t=x; x=y, y=t-(a/b)*y;
    return g;
}
int main(){
    int T=R; while(T-->0){
        LL a=R, b=R, c=R, x, y, k;
        LL g=Exgcd(a, b, x, y);
        /* 判断有无解 */
        if(c%g) { puts("-1"); continue; }
        /* 求 x1 及 y1 */
        x*=c/g, y*=c/g; 
        /* q, p 同 dx, dy */
        LL q=b/g, p=a/g; 
        /* 将 x 调整为最小正整数解 */
        if(x<=0) k=ceil((1.0-x)/q), x+=k*q, y-=k*p;
        else k=(x-1)/q, x-=k*q, y+=k*p;
        /* 判断有无正整数解 */
        if(y>0)
            /* 根据询问回答即可, (y-1)%p+1 是为了使 y 是正整数 */
            printf("%lld %lld %lld %lld %lld
", (y-1)/p+1, x, (y-1)%p+1, x+(y-1)/p*q, y);
        else
            printf("%lld %lld
", x, y+(LL)ceil((1.0-y)/p)*p);
    }
    return 0;
}

后记

参考 dengyaotriangle 的博客 - 题解 P5656 【【模板】二元一次不定方程(exgcd)】

原文地址:https://www.cnblogs.com/FUXyao/p/14754914.html