扩展欧几里德
问题:求解 s1 + v1*t = s2 + v2*t - k*m (v1<v2)
已知:s1, s2, v1, v2, m
求解该式子的算法我们称为扩展欧几里德算法。
该算法分为两个部分:
(1) 判定是否存在解
对于形如"Ax+By=C"的式子,其存在解的条件为C为A和B最大公约数的整数倍。
我们将A和B的最大公约数记为gcd(A,B)。因此其有解的条件是C=n*gcd(A,B)。
那么我们应该如何来求解gcd(A,B)呢?
一个朴素的算法是枚举1~min(A,B),最大的一个能同时被A,B整除的数即gcd(A,B)。显然这个算法是非常没有效率的。
为了求解gcd(A,B),欧几里德提出了一个辗转相除法:
首先要证明一个定理:gcd(A,B) = gcd(B, A mod B)
证明: 假设A = k*B+r,有r = A mod B。不妨设d为A和B的一个任意一个公约数,
则有A = pd, B = qd。 由r = A - k*B = pd - k*qd = (p - kq)*d,所以有d也为r的约数,
因此d是B和A mod B的公约数。 由于对任意一个A和B的公约数都满足这个性质,gcd(A,B)也满足,
因此有gcd(A,B)=gcd(B,A mod B)。
(2) 求解
在判定有解之后,我们需要在其基础上再求出一组(x,y)。由于A,B,C均是gcd(A,B)的整数倍,因此可以将它们都缩小gcd(A,B)倍。即A'=A/gcd(A,B),B'=B/gcd(A,B),C'=C/gcd(A,B)。
化简为A'x+B'y=C',gcd(A',B')=1,即A',B'互质。
此时,我们可以先求解出A'x+B'y=1的解(x',y'),再将其扩大C'倍,即为我们要求的最后解(x,y)=(C'x', C'y')。
那么接下来我们来研究如何求解A'x+B'y=1:
假设A>B>0,同时我们设:
A * x[1] + B * y[1] = gcd(A, B) B * x[2] + (A mod B) * y[2] = gcd(B, A mod B)
已知gcd(A,B)=gcd(B, A mod B),因此有:
A * x[1] + B * y[1] = B * x[2] + (A mod B) * y[2] => A * x[1] + B * y[1] = B * x[2] + (A - kB) * y[2] // A = kB + r => A * x[1] + B * y[1] = A * y[2] + B * x[2] - kB * y[2] => A * x[1] + B * y[1] = A * y[2] + B * (x[2] - ky[2]) => x[1] = y[2], y[1] = (x[2] - ky[2])
利用这个性质,我们可以递归的去求解(x,y)。
其终止条件为gcd(A, B)=B,此时对应的(x,y)=(0,1)
将这个过程写成伪代码为:
extend_gcd(A, B): If (A mod B == 0) Then Return (0, 1) End If (tempX, tempY) = extend_gcd(B, A mod B) x = tempY y = tempX - (A / B) * tempY Return (x, y)
由于A'B'是互质的,所以可以将A'x'+B'y'=1扩展为:
A'x'+B'y'+(u+(-u))A'B'=1 => (x' + uB')*A' + (y' - uA')*B' = 1 => X = x' + uB', Y = y' - uA'
可以求得最小的X为(x'+uB') mod B',(x'+uB'>0)
同时我们还需要将X扩大C'倍,因此最后解为:
x = (x'*C') mod B'
若x<0,则不断累加B',直到x>0为止。
solve(s1, s2, v1, v2, m): A = v1 - v2 B = m C = s2 - s1 If (A < 0) Then A = A + m // 相对距离变化 End If D = gcd(A, B) If (C mod D) Then Return -1 End If A = A / D B = B / D C = C / D (x, y) = extend_gcd(A, B) x = (x * C) mod B While (x < 0) x = x + B End While Return x
解:(v1 - v2)*t + k*m = s2 - s1
令 A = v1 - v2, B = m, C = s2 - s1;
当 A 小于 0 时, 相对距离变化, A += m
所以问题变成求解 A*x + B*y = C
令 D = gcd(A, B), 如果 C%D!= 0 则无解
A' = A/D, B' = B/D, C' = C/D
所以问题变成求解 A'*x + B'*y = C'
这时,gcd(A', B') = 1, 即 A', B' 互质
令 (x', y') = (x/C', y/C')
所以问题变成求解 A'*x' + B'*y' = 1
对于这个问题可以由递归解决。
因为 gcd(a, b) == gcd(b, a%b) == 1
所以 a*x[1] + b*y[1] = gcd(a, b)
b*x[2] + (a%b)*y[2] = gcd(b, a%b)
所以,x[1] = y[2], y[1] = x[2] - k*y[2] (k = a / b)
// 扩展欧几里德 求 Ax + By = C 的一组解 #include <iostream> #include <string> #include <algorithm> #define ll long long using namespace std; typedef struct node{ ll x, y; }Node; ll s1, s2, v1, v2, m; // 辗转相除法求最大公约数 ll gcd(ll a, ll b){ if(a%b==0) return b; else return gcd(b, a%b); } // 求解 Ax + By = 1 的一组解 A,B互质 // gcd(A, B) = gcd(B, A%B) = 1 // A*x[1] + B*y[1] = B*x[2] + (A%B)*y[2] // x[1] = y[2], y[1] = x[2] - k*y[2] Node extend_gcd(ll A, ll B){ Node t; t.x = 0; t.y = 1; if(A%B == 0){ return t; } Node tt = extend_gcd(B, A%B); t.x = tt.y; t.y = tt.x - (A/B)*tt.y; return t; } // 求 Ax + By = C 的一组解 ll solve(ll s1, ll s2, ll v1, ll v2, ll m){ ll A = v1 - v2; ll B = m; ll C = s2 - s1; if(A < 0){ A = A + m;// 相对距离变化 } ll D = gcd(A, B); if (C % D){ return -1; } A = A / D; B = B / D; C = C / D; Node tt = extend_gcd(A, B); tt.x = (tt.x * C)%B; while(tt.x < 0) tt.x += B; return tt.x; } // 求解 (v1 - v2)*t + k*m = s2 - s1 // A*t + k*B = C int main(){ cin >> s1 >> s2 >> v1 >> v2 >> m; cout << solve(s1, s2, v1, v2, m) << endl; }