一、欧几里得算法:
欧几里得算法(又称辗转相除法)是一种求解两个整数最大公约数(Greatest Common Divisor, GCD)的经典算法。它的基本思想如下:
两个整数的最大公约数等于其中较小的那个数与两数相除所得的余数的最大公约数
可以通过递归来实现:
//求最大公约数——欧几里得算法,即辗转相除法
//gcd(a,b) <=> gcd(b,a%b)
int gcd(int a, int b) {
if (b == 0) return a; //边界
return gcd(b, a % b); //假设已经求出了a%b的最大公约数
}
欧几里得算法的应用:——求整数点:
已知a(x1,y1), b(x2,y2)都是整数点, c(x3,y3)在a,b连线上并满足c点也是整数点,求c点的可能坐标:
//求整数点:已知a(x1,y1), b(x2,y2)都是整数点, c(x3,y3)在a,b连线上并满足c点是整数点,求c点的坐标
vector<pair<int,int>> GetC(int x1, int y1, int x2, int y2) {
int x = abs(x2 - x1);
int y = abs(y2 - y1);
int gcd_xy = gcd(x, y);//求曼哈顿距离差最大公约数
vector<pair<int,int>> res;//存放结果
int dx = x / gcd_xy;//x方向的步长
int dy = y / gcd_xy;//y方向的步长
int cnt = 1;
int dir_x = ((x1 - x2) < 0) ? 1 : -1;//判断从a到b,x方向的方向
int dir_y = ((y1 - y2) < 0) ? 1 : -1;
while (cnt<gcd_xy) {
res.push_back({ x1 + dir_x * cnt * dx,y1 + dir_y * cnt * dy });
cnt++;
}
return res;
}
二、贝祖(裴蜀)等式:
1.定义:
对任何整数 a a a、 b b b和它们的最大公约数 d d d,关于未知数 x x x和 y y y的线性方程(称为裴蜀等式):
a x + b y = m ax + by = m ax+by=m
有整数解当且仅当 m m m是 d d d的倍数。
裴蜀等式有解时必然有无穷多个整数解,每组解 x , y x, y x,y称为裴蜀数,可用扩展欧几里得算法(Extended Euclidean algorithm)获得。
特别地,方程 a x + b y = 1 ax + by = 1 ax+by=1 有整数解当且仅当整数 a a a和 b b b互素
2.求解裴蜀公式:
扩展欧几里得算法求裴蜀公式的一个解:
设 a x + b y = d ax+by=d ax+by=d,假设 d = g c d ( a , b ) d=gcd(a,b) d=gcd(a,b)
- 扩展欧几里得算法就是在求 a , b a,b a,b 的最大公约数 d = gcd ( a , b ) d=\text{gcd}(a,b) d=gcd(a,b) 的同时,可求出贝祖等式 a x + b y = d ax + by = d ax+by=d 的一个解 ( x 0 , y 0 ) (x_0,y_0) (x0,y0)
(1)如何推导?(可以利用相邻两步之间的过程进行推导,提出系数a,b,对照x和y即可得出,不太好打出来,可自行推导或查找资料)
- x = y 1 x = y_1 x=y1
- y = x 1 − a b y 1 y = x_1 - \frac{a}{b}y_1 y=x1−bay1
下面的ext_gcd
函数即是求解扩展的欧几里得算法;
调用LinearEquation
函数求解贝祖等式的一个解,并乘以倍数m/d,保证解的正确性。
- 如果 a x + b y = m ax+by=m ax+by=m,m不是 g c d ( a , b ) gcd(a,b) gcd(a,b)的倍数,则该方程无解:
代码如下:
#define _CRT_SECURE_NO_WARNINGS
#include<iostream>
#include<vector>
using namespace std;
//扩展欧几里得算法:求解ax+by=m的整数解
int x, y;
int ext_gcd(int a, int b) {
if (b == 0) {
x = 1;
y = 0;
return a;
}
int res = ext_gcd(b, a % b); //假设已经求出了b,a%b的整数解
int temp = x;
x = y;
y= temp - a / b * y;
return res;
}
void LinearEquation(int a, int b, int m) {
int d = ext_gcd(a, b);
if (m % d != 0) {
printf("无解\n");
return;//无解
}
int n = m / d;
x *= n;
y *= n;
cout << "gcd: "<<d<<endl;
}
int main(void)
{
int a, b, m;
cin >> a >> b >> m;
LinearEquation(a, b, m); //ax+by=m的整数解
cout <<"x: " <<x << " " <<"y: " << y << endl;
return 0;
}
求裴蜀公式的通解:
(2)通解?
- x = x 0 + b g c d ∗ t x = x_0 + \frac{b}{gcd}*t x=x0+gcdb∗t − − > --> −−>所有的 x 对 b 模
- y = y 0 − a g c d ∗ t y = y_0 - \frac{a}{gcd}*t y=y0−gcda∗t − − > --> −−>所有的 y 对 a 模
或者
- x = x 0 − b g c d ∗ t x = x_0 - \frac{b}{gcd}*t x=x0−gcdb∗t − − > --> −−>所有的 x 对 b 模
- y = y 0 + a g c d ∗ t y = y_0 + \frac{a}{gcd}*t y=y0+gcda∗t − − > --> −−>所有的 y 对 a 模
即x如果加的话,y就要对应的减,从而保证二者的平衡,反之亦然。
由于贝祖公式一旦有解,则有无穷解,所以这里并未给出具体的求解代码,可根据需求求解相应的解。
(3)如果想得到 x > 0 x > 0 x>0 的第一个解?
- b / = d b/=d b/=d
- x = ( x 0 % b + b ) % b x = (x_0 \% b + b) \% b x=(x0%b+b)%b
利用基本的取模运算即可求解,x和y只需求出一个即可得出另一个。
3.裴蜀公式的应用:
让我们来看一看2016年第七届蓝桥杯B组C/C++决赛题——一步之遥
题目如下:
从昏迷中醒来,小明发现自己被关在X星球的废矿车里。
矿车停在平直的废弃的轨道上。
他的面前是两个按钮,分别写着“F”和“B”。
小明突然记起来,这两个按钮可以控制矿车在轨道上前进和后退。
按F,会前进97米。按B会后退127米。
透过昏暗的灯光,小明看到自己前方1米远正好有个监控探头。
他必须设法使得矿车正好停在摄像头的下方,才有机会争取同伴的援助。
或许,通过多次操作F和B可以办到。
矿车上的动力已经不太足,黄色的警示灯在默默闪烁...
每次进行 F 或 B 操作都会消耗一定的能量。
小明飞快地计算,至少要多少次操作,才能把矿车准确地停在前方1米远的地方。
请填写为了达成目标,最少需要操作的次数。
求解方法也很简单,就是解一个贝祖等式:
97
x
−
127
y
=
1
97x-127y=1
97x−127y=1
这里的正负号其实影响不大,看怎么理解了:
可以理解为:
- 97 97 97和 127 127 127求 g c d gcd gcd,得 x = 55 , y = − 42 x=55,y=-42 x=55,y=−42;即向前 55 ∗ 97 m 55*97m 55∗97m,再向前 − 42 ∗ 127 -42*127 −42∗127m;
- 97 97 97和 − 127 -127 −127求 g c d gcd gcd,得 x = 55 , y = 42 x=55,y=42 x=55,y=42;即向前 55 ∗ 97 55*97 55∗97m,再向后 42 ∗ 127 42*127 42∗127m;
代码如下,和扩展欧几里得算法基本一致:
#define _CRT_SECURE_NO_WARNINGS
#include<iostream>
#include<vector>
using namespace std;
//扩展欧几里得算法:求解ax+by=m的整数解
int x, y;
int ext_gcd(int a, int b) {
if (b == 0) {
x = 1;
y = 0;
return a;
}
int res = ext_gcd(b, a % b); //假设已经求出了b,a%b的整数解
int temp = x;
x = y;
y= temp - a / b * y;
return res;
}
bool LinearEquation(int a, int b, int m) {
int d = ext_gcd(a, b);
if (m % d != 0) {
printf("无解\n");
return false;//无解
}
int n = m / d;
x *= n;
y *= n;
cout << "gcd: "<<d<<endl;
return true;
}
int main(void)
{
int a, b, m;
cin >> a >> b >> m;
if (LinearEquation(a, b, m)) {
//ax+by=m的整数解
cout << "x: " << x << " " << "y: " << y << endl;
cout << "res: " << abs(x) + abs(y) << endl;
}
return 0;
}
当然,题目中求的是最少操作次数,上面并没有对其进行判断。但通过上面的(2):通解,可以验证这就是最少的操作次数,因为操作次数是x和y的绝对值的和,在当前情况下,无论是增加x,减少y;还是减少x,增加y,都会增大操作次数。
扩展:
所以,通过扩展欧几里得求得的一组解的绝对值和一定是最小的。
证明:
对整个递归过程进行讨论:
(1)递的过程——求出每层的a和b:(设层数c从1->n)
假设我们输入的a,b满足
a
>
=
b
a>=b
a>=b,若不满足就交换a,b。那么在求gcd的过程中a总是大于b,即
a
/
b
a/b
a/b总是大于1。
(2)归的过程——更新每层的x和y:(层数c从n->1)
因为递归求解gcd的最后一步必定是
d
x
+
0
y
=
d
,
d
=
g
c
d
(
a
,
b
)
dx + 0y = d,d=gcd(a,b)
dx+0y=d,d=gcd(a,b)
此时
x
=
1
,
y
=
0
x=1,y=0
x=1,y=0;
然后利用递推公式求解上一层的x和y。
- x = y 1 x = y_1 x=y1
- y = x 1 − a b y 1 y = x_1 - \frac{a}{b}y_1 y=x1−bay1
对每层来说
a
/
b
a/b
a/b大于1,记为q
注意:
- 不同层的a,b大小是不一定相同的;
- q只是一个大于1的记号,也不一定相同;
- 设c为层数;
- sa为sum of abs ;
迭代如下:
c
=
n
,
x
=
1
,
y
=
0
,
s
a
=
1
c=n,x=1,y=0,sa=1
c=n,x=1,y=0,sa=1
c
=
n
−
1
,
x
=
0
,
y
=
1
,
s
a
=
1
c=n-1,x=0,y=1,sa=1
c=n−1,x=0,y=1,sa=1
c
=
n
−
2
,
x
=
1
,
y
=
−
q
,
s
a
=
q
+
1
c=n-2,x=1,y=-q,sa=q+1
c=n−2,x=1,y=−q,sa=q+1
c
=
n
−
3
,
x
=
−
q
,
y
=
1
+
q
∗
q
,
s
a
=
q
2
+
q
+
1
c=n-3,x=-q,y=1+q*q,sa=q^2+q+1
c=n−3,x=−q,y=1+q∗q,sa=q2+q+1
c
=
n
−
4
,
x
=
1
+
q
∗
q
,
y
=
−
q
−
q
(
1
+
q
∗
q
)
,
s
a
=
q
3
+
q
2
+
2
q
+
1
c=n-4,x=1+q*q,y=-q-q(1+q*q),sa=q^3+q^2+2q+1
c=n−4,x=1+q∗q,y=−q−q(1+q∗q),sa=q3+q2+2q+1
……
通过归纳可以得到随着层数的降低,sa一定是逐渐增加的,并且每次的迭代都是最优的,所以当迭代到c=1时,得到的sa一定也是最优的。
证明不一定正确,欢迎大家前来交流。
当然,不放心的话也可以
(1)求出
x
>
0
x > 0
x>0 的第一个解,然后求出y,即
- b / = d b/=d b/=d
- x = ( x 0 % b + b ) % b x = (x_0 \% b + b) \% b x=(x0%b+b)%b
(2)求出
y
>
0
y > 0
y>0 的第一个解,然后求出x
(3)求出
y
<
0
y < 0
y<0 的第一个解,然后求出x
(4)求出
x
<
0
x < 0
x<0 的第一个解,然后求出y
通过比较上面4组x,y绝对值的大小得出最优的组合。
这样求解出来的x,y的abs和也一定是最优的