Sanshuiii
0%

解不定方程

Posted on Edited on In

传送门：http://acm.hdu.edu.cn/showproblem.php?pid=1356

由题意，就是要解一个不定方程ax+by=d，要求（abs(x)+abs(y)）最小。

一.exgcd

先从exgcd说起。由裴蜀定理可知，ax+by=gcd(a,b)必定存在整数解(x,y)。那么有：

ax+by=gcd(a,b)=gcd(b,a%b)=bx'+(a%b)y'                （1.1）

其中等号1，3用到了裴蜀定理，等号2用到了gcd的性质。

由a%b=a-a/b(本文中所有的除法都是整数除法)那么我们把等式最右边变形，可以得到：

ax+by=bx'+ay'-(a/b)y'=ay'+b(x'-(a/b)*y')            （1.2）

于是得到一组解：

x=y'                                （1.3）
y=x'-(a/b)*y'                            （1.4）

所以可以用递归的方式求解。

最后考虑一下边界条件，当b==0时，从线性同余方程的角度看x=1，y=0显然是一组可行解。

1
2
3
4
5
6
7
8
9
10
11
int exgcd(int a,int b,int &x,int &y){
	if(!b){
		x=1;
		y=0;
		return a;
	}
	int d=exgcd(a,b,x,y);
	swap(x,y);//此时x=y',y=x'
	y-=a/b*x;//考虑x，y的实际情况，也就是 y=x'-(a/b)*y'
	return d;//顺便吧gcd也求出来
}

二.求通解

再考虑这道题，等号右边并不是gcd，而是一个d，因为题目告诉你一定有解，所以gcd一定是d的因子，我们把x和y都乘上一个rate=(d/gcd(a,b)),得到一组初始的

x0=x*rate                （2.1）
y0=y*rate                （2.2）

接下来就要去找对x0和y0加加减减，来求这个MIN(abs(x0+dx)+abs(y0-dy))。

我们先考虑

再考虑dx和dy的关系：

a*x0+b*y0=d                    （2.3）

对方程变形得到：

a*(x0+dx)+b*(y0-dy)=d                （2.4）

把里面的式子提出来得到：

a*dx-b*dy=0                    （2.5）

那么就是:

dx/dy=b/a                    （2.6）

考虑gcd(a,b)，可以得到：

dx/dy=(b/gcd)/(a/gcd)                （2.7）

此时(b/gcd)与(a/gcd)是互质的所以可以取到的最小的dx和dy就是：

dx=(b/gcd)                    （2.8）
dy=(a/gcd)                    （2.9）

那么方程的通解一定满足如下条件：

ansx=x0+k*dx                    （2.10）
ansy=y0-k*dy                    （2.11）

三.问题转化

知道了通解，那么哪一组解才是最优的呢？

因为dx！=dy，不妨设dx>dy，我们假设此时x满足

|x|>>|dx|(>>表示远大于)

那么我们让

x+=t*dx                    （3.1）
y-=t*dy                    （3.2）

在x和y均不变号的情况下，答案会变小t*(dx-dy)，是一组更优的解。所以最优解一定满足:

|x|<|dx|                （3.3）

那么问题就可以转化为变成了求方程关于x的最小正整数解。

但是问题在于，根据式(3.3)，去掉绝对值后，x可能是负的。但是问题在于，我们要求的方程

a*x+b*y=d                (3.4)

可以变形为

y=(d-a*x)/b                (3.5)

我们上面已经证明，通过上述方程构造出来的x，上面方程中的y一定是有解的，那么对于一个确定的x，我们总希望y最小，因为dx>dy，那么如果x取到负数，那么一定会大大增加y的值，所以很有可能会更糟。所以我们只要考虑x取最小正整数的情况就可以了。（我直言了，我不会证明，但是我的直觉告诉我这样做很对，oj也告诉我这样做很对）。

反之，如果dy>dx，那么考虑y的最小正整数解即可。

四.求最小正整数解

那这就非常容易了（说白了就是取模）：

ansx=（ansx%dx+dx）%dx

然后用式（3.5）：

ansy=(d-a*ansx)/b

就完事了。

五.贴代码

看清题意，多写两个if就完事了

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
#include<iostream>
using namespace std;
#define ABS(x) (((x)>0)?(x):-(x))

int x,y,x1,x2,y1,y2,dx,dy;
int exgcd(int a,int b,int &x,int &y){
	if(!b){x=1,y=0;return a;}
	int d=exgcd(b,a%b,x,y);
	int t=x;x=y,y=t-(a/b)*x;
	return d;
} 

void read(int &x)
{
	x=0;int f=1;char ch=getchar();
	while(ch<'0'||ch>'9'){if(ch=='-')f=-1;ch=getchar();}
	while(ch>='0'&&ch<='9'){x=x*10+ch-'0';ch=getchar();}
	x*=f;return;
}

int main(){
	int a,b,d;
	while(1){
		read(a),read(b),read(d);
		if(!a&&!b&&!d)break;
		int gcd=exgcd(a,b,x,y);
		int rate=d/gcd;
		dx=b/gcd,dy=a/gcd;
		x1=x2=x*rate;
		y1=y2=y*rate;
		x1=(x1%dx+dx)%dx,y1=ABS(d-a*x1)/b;
		y2=(y2%dy+dy)%dy,x2=ABS(d-b*y2)/a;
		if(x1+y1<x2+y2||(x1+y1==x2+y2&&a*x1+b*y1<a*x2+b*y2))cout<<x1<<" "<<y1<<endl;
		else cout<<x2<<" "<<y2<<endl;
	}
}