参考资料
题意简述
给定
1≤a,b,c≤109,求解不定方程:
ax+by=c
- 若无整数解,输出
-1;
- 若有正整数解,输出正整数解的数量、x 的最小值、y 的最小值、x 的最大值、y 的最大值;
- 若无正整数解,输出 x 的最小正数值和 y 的最小正数值。
基础知识
一次不定方程
一次不定方程(Linear Diophantine Equation)是形如以下形式的方程:
a1x1+a2x2+⋯+anxn=b
其中
ai 和
b 都是整数,我们会研究
xi 的整数解。
特别地,当
n=2 时为二元一次不定方程:
a1x1+a2x2=b
我们可以调整一下变量名:
ax+by=c
裴蜀定理
裴蜀定理(Bézout's Lemma)给出了二元一次不定方程 有整数解 的充要条件。
对于整数
a,b(不全为
0),存在整数
x,y 使得
ax+by=gcd(a,b)。
换句话说,
ax+by=c 有整数解当且仅当
gcd(a,b)∣c。
欧几里得算法
欧几里得算法(辗转相除法)可以求两个整数的 最大公约数(Greatest Common Divisor,GCD)。
gcd(a,b)=gcd(b,amodb)
证明:设
a=kb+c。若
d∣a,b,则
d∣c;若
d∣b,c,则
d∣a。故
a,b 与
b,c 的公约数集合相同,因此
gcd(a,b)=gcd(b,c)=gcd(b,amodb)。
显然
amodb<b,因此可以递归求解。特别地,
gcd(a,0)=a。
gcd(a,b)={agcd(b,amodb)b=0otherwise
CPPll Gcd(ll a,ll b)
{
return !b?a:Gcd(b,a%b);
}
扩展欧几里得算法
扩展欧几里得算法(Extended Euclidean algorithm,EXGCD)可以求
ax+by=gcd(a,b) 的
一组整数解。
ax+by=gcd(a,b)⟺ay+b(x−⌊ba⌋y)=gcd(a,b)
设一组整数解为
(x0,y0):
ax0+by0=gcd(a,b)
bx1+(amodb)y1=gcd(b,amodb)
根据欧几里得算法:
gcd(a,b)=gcd(b,amodb)
所以:
ax0+by0=bx1+(amodb)y1
又因为:
amodb=a−⌊ba⌋×b
所以:
ax0+by0=bx1+(a−⌊ba⌋×b)y1
ax0+by0=ay1+bx1−⌊ba⌋×by1=ay1+b(x1−⌊ba⌋y1)
所以:
x0=y1,y0=x1−⌊ba⌋y1
将
x1,y1 不断代入递归求解直至
gcd 为
0 递归
x=1,y=0 回去求解。
解题思路
我们先用
扩展欧几里得算法 求出最大公约数
g=gcd(a,b) 和一组整数解
(x0,y0)。
ax0+by0=g
根据裴蜀定理,当
c 不是
g 的倍数时无解,输出
-1。
将方程两边同时乘以
gc:
gc⋅ax0+gc⋅by0=gc⋅d=c
此时我们令
x0←gc⋅x0,y0←gc⋅y0,就得到了原问题的一组解
(x0,y0)。
dx=gb,dy=ga
而所有正整数解分别为:
(xmin,ymax),(xmin+dx,ymax−dy),…,(xmax,ymin)
显然
xmin∈[1,dx],ymin∈[1,dy],可以通过取模得到:
xmin=(x0−1)moddx+1,ymin=(y0−1)moddy+1
当
x 最小时
y 最大,
y 最小时
x 最大,因此:
xmax=ac−bymin,ymax=bc−axmin
因此正整数的数量为:
n=dxxmax−xmin+1=dyymax−ymin+1
参考代码
CPP#include <bits/stdc++.h>
using namespace std;
using ll=long long;
tuple<ll,ll,ll> exgcd(ll a,ll b)
{
if(!b)return {a,1,0};
auto [g,x,y]=exgcd(b,a%b);
return {g,y,x-a/b*y};
}
int main()
{
ios::sync_with_stdio(false);
cin.tie(nullptr);
int T;
cin>>T;
while(T--)
{
ll a,b,c;
cin>>a>>b>>c;
auto [g,x,y]=exgcd(a,b);
if(c%g){cout<<-1<<'\n';continue;}
x*=c/g;y*=c/g;
ll dx=b/g,dy=a/g;
ll x_min=(x%dx+dx-1)%dx+1,y_min=(y%dy+dy-1)%dy+1,x_max=(c-b*y_min)/a,y_max=(c-a*x_min)/b;
if(y_max>0)cout<<(x_max-x_min)/dx+1<<' '<<x_min<<' '<<y_min<<' '<<x_max<<' '<<y_max<<'\n';
else cout<<x_min<<' '<<y_min<<'\n';
}
return 0;
}