gcd算法以及exgcd

本文介绍了欧几里得算法,即辗转相除法,用于求最大公约数,还可求最小公倍数。同时阐述了扩展算法,包括解ax+by=gcd(a,b)、ax+by=c及解集的方法,最后通过练习展示了扩展gcd在特定问题中的应用。

摘要生成于 C知道 ,由 DeepSeek-R1 满血版支持, 前往体验 >

1.欧几里得算法

欧几里得是求最大公约数的经典算法,又称辗转相除法。
gcd函数就是用来求(a,b)的最大公约数的。
gcd函数的基本性质:
gcd(a,b)=gcd(b,a)=gcd(-a,b)=gcd(|a|,|b|)
gcd(a,b)=gcd(b,a mod b)

证明:a可以表示成a = kb + r,则r = a mod b
	 假设d是a,b的一个公约数,
	 则有 d | a , d | b ,而r = a - kb,
	 因此 d | r 因此d是(b,a mod b)的公约数
	 假设d 是(b,a mod b)的公约数,
	 则 d | b , d | r ,但是a = kb +r
	 因此d也是(a,b)的公约数
	 因此(a,b)和(b,a mod b)的公约数是一样的,其最大公约数也必然相等,得证

所以通过递归,最终可以使得其中一个参数为零,则另一个不为零的数则为最大公约数
板子如下:

typedef long long ll;
ll gcd(ll a,ll b){
	return b==0?a:gcd(b,a%b);
}

当然也可以求最小公倍数
最小公倍数lcm=a * b/gcd(a,b)

2.扩展算法
1.解ax+by=gcd(a,b)

对于不完全为 0 的非负整数 a,b,gcd(a,b)表示 a,b 的最大公约数,必然
存在整数对 x,y ,使得 gcd(a,b)=ax+by。

求解 x,y的方法的理解:
	设 a>b。
1,显然当 b=0,gcd(a,b)=a。此时 x=1,y=0;
2,a>b>0 时:
	设 ax1+ by1= gcd(a,b);
	bx2+ (a mod b)y2= gcd(b,a mod b);
根据朴素的欧几里德原理有 gcd(a,b) = gcd(b,a mod b);
	则:ax1+ by1= bx2+ (a mod b)y2;
	即:ax1+ by1= bx2+ (a - [a / b] * b)y2=ay2+ bx2- [a / b] * by2;
说明: 
	a-[a/b]*b即为mod运算。[a/b]代表取小于a/b的最大整数。
也就是ax1+ by1 == ay2+ b(x2- [a / b] *y2);
根据恒等定理得:
		x1=y2;
		y1=x2- [a / b] *y2;
这样我们就得到了求解 x1,y1 的方法:x1,y1 的值基于 x2,y2.
上面的思想是以递归定义的,因为 gcd 不断的递归求解一定会有个时候 b=0,所以递归可以结束。

板子如下:

int gcd(int a,int b,int &x,int &y){
    if (b==0){
        x=1,y=0;
        return a;
    }
    int q=gcd(b,a%b,y,x);
    y-=a/b*x;
    return q;
}

可以这样思考:
对于a’=b,b’=a%b 而言,我们求得 x, y使得 a’x+b’y=Gcd(a’,b’)
由于b’=a % b=a - a / b * b
(这里的 / 是程序设计语言中的除法)
那么可以得到:
ax+by=gcd(a,b) =>
bx+(a - a / b * b)y = gcd(a’, b’) = gcd(a, b) =>
ay +b(x - a / b * y) = Gcd(a, b)
因此对于a和b而言,他们的相对应的p,q分别是 y和 (x-a/b*y)。
板子如下:

nt gcd(int a,int b,int &x,int &y){
    if (b==0){
        x=1,y=0;
        return a;
    }
    int q=gcd(b,a%b,x,y);
    int t=x;x=y;
    y=t-a/b*y;
    return q;
}
2.解ax+by=c

使用扩展欧几里德算法解决不定方程的办法
对于不定整数方程pa+qb=c,若 c mod Gcd(a, b)=0,则该方程存在整数解,
否则不存在整数解。
通过上面的程序即可求出一组x,y使得gcd(a,b)=ax+by
找其所有解的方法如下:
在找到x * a+y * b = Gcd(a, b)的一组解x0,y0后,
可以得到x * a+y * b = c的一组解:
x1 = x0 * (c/Gcd(a,b)),
y1 = y0 * (c/Gcd(a,b)),

3.解集

假设我们得出ax+by=c一组解x与y,如果x每减去一个t(t为常数),
那么y至少必须加上一个a’* t/b’
(其中a’=a/gcd(a,b),b’=b/gcd(a,b)),
才能使等式两边成立,即
  a(x-t)+b(y+a’ * t/b’)=ax+by=c
由于我们需要变换后的x和y仍然是整数,由于t已经是整数,即x-t也一定是整数,
现在只需考虑a’ * t/b’是否为整数,由于gcd(a’,b’)=1,
即a’一定不是b’的倍数,故必须令t=t * b’(即c必须是b’的倍数),才使得a’ * t/b’为整数,即解集为:
  x=x-t * b/gcd(a,b),y=y+t * a/gcd(a,b),t为任意整数
特别的,当gcd(a,b)=1,解集为:x=x-t * b,y=y+t * a

练习

1.扩展gcd-时间复杂性

题目内容:
计算循环语句的执行频次 for(i=A; i!=B ; i+=C) x+=1;
其中A,B,C,i都是k位无符号整数。

输入描述

A B C k, 其中0<k<32

输出描述

输出执行频次数,如果是无穷,则输出“forever”

输入样例

3 7 2 16

输出样例

2

k位无符号整数则数据范围为0~2k(具体含义可以搜索无符号整数溢出)
所以设y,题目意思可以理解为
A+x * C - B=y * 2k => x * C+y * 2k=B-A
于是可以用exgcd求解
令a=C,b=2k,c=B-A
首先判断B-A 是否是gcd(a,b)倍数,如果不是则无法求出一组x。
求出一组x0,y0后,则通解为
x=x0 * a/gcd(a,b)*t
y=y0 - b/gcd(a,b)*t

令n=a/gcd(a,b),m= b/gcd(a,b)
则可以得到
①x0 * a+y0 * b=c
②a(x0+n * t)+b(y0-m * t)=c
联立上式
可得
a * t * n - b * t * m = 0
a * n - b * m=0
也就是a * n=b * m,要想n,m最小就是使得a * n,b * m都为a,b的最小公倍数:
a * n = (a * b)/gcd(a,b) => n = b/gcd(a,b);
同理:m = a /gcd(a,b);
则:
x0 = (x % n + n)%n
y0 = (y % m + m)%m

#include<iostream>
using namespace std;
//A + xC - B = y * 2^k
int exgcd(int a,int b,int &x,int &y){
	if(b==0){
		x=1;y=0;
		return a;
	}
		int d=exgcd(b,a%b,x,y);
		int t=x;
		x=y;
		y=t-a/b*y;
	return d;
}
int main(){
	int A,B,C,k;
	cin >> A >> B >> C >> k;
	//特判A=B=C=K=0 
	if(A==0&&A==B&&B==C&&C==k) return !(cout << "forever");
	int b=1<<k,a=C,c=B-A,x=10,y=10;
	int d=exgcd(a,b,x,y);
	if(c%d){
		return !(cout << "forever");
	} 
	x=(x*(c/d))%b;
	x=(x%(b/d)+(b/d))%(b/d);//防止出现负数解
	cout << x;
	return 0;
} 
评论 8
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值