晚上做到的一个扩欧的水题(?)
wa了好几发感觉自己药丸…重新推了一遍公式才发现自己打错了orz
借此复习一下扩欧吧…orz
题目:http://www.lydsy.com/JudgeOnline/problem.php?id=1407
看到这里的人应该都知道扩欧是干嘛的吧…如果不知道也没关系下面有…
欧几里得辗转相除:$gcd(a,b)=gcd(b,a mod b)$,证明可以看下面。
扩展欧几里得:
考虑方程$ax+by=gcd(a,b)$,如果$b=0$则直接有$x=1,y=0$,否则通过上面的式子可以得到$$ax+by=gcd(a,b)=gcd(b,amod b)=bx+(amod b)y=bx+(a-b*lfloor frac{b}{a} floor)y$$然后根据$a,b$合并同类项:$$ax+by=ay+b*(x-lfloor frac{a}{b} floor y)$$ 看,他出现了递归的形式!所以我们就可以递归的求了(返回的是$gcd(a,b)$):
inline 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 tx,ty;tx=y;ty=x-a/b*y; x=tx;y=ty; return d; }
好了刚刚说过下面要证明辗转相除法,假如我们通过辗转相除法求出了$d=gcd(a,b)=ax+by$,我们怀疑有一个更大的公约数$d'$因为某些原因没有被找到,但是这是不存在的~因为如果有这样的一个$d'$,作为$a,b$的$gcd$一定整除$ax+by$,一定整除$d$,所以一定$leq d$。此外也容易验证求出来的$d$一定整除$a,b$,所以就是真的最大公约数啦~
可以证明复杂度是$O(log n)$的。
不过我们一般不是直接用上面的方程,实际上更多的是讨论$ax+by=c$这个不定方程的解。
把方程两边同除掉$gcd(a,b)$:$$frac{a}{gcd(a,b)}x+frac{b}{gcd(a,b)}y=frac{c}{gcd(a,b)}$$ 这里的其他项都是整数,所以要有解的充要条件就是$gcd(a,b)|c$了。如果它有解的话就先求出$ax+by=gcd(a,b)$的一组解$(x_0,y_0)$,将解代到方程里之后再将两边都乘上$k=frac{c}{gcd(a,b)}$就得到了$kax_0+kby_0=c$,从而原方程的一组特解就是$x=kx_0,y=ky_0$了。
如果我们要求出其他的解(通解)怎么办?也就是我们想要找一组最小的$Delta x,Delta y$使得$a(x+Delta x)+b(y+Delta y)=c$,结合$ax+by=c$我们可以得到$a Delta x+b Delta y =0$
哎呀这个要怎么弄呀,构造!令$d=gcd(a,b)$,根据$ab-ba=0$我们给他同除掉$d$就可以得到一个有用的等式:$afrac{b}{d}-bfrac{a}{d}=0$,很显然$Delta x=frac{b}{gcd(a,b)},Delta y=-frac{a}{gcd(a,b)}$,而且这里的$Delta x,Delta y$一定是满足条件的最小值了,因为不能再除了。
好了现在回到我们的题目,$n$非常小,又保证了答案的范围,我们可以直接枚举答案$m$然后判断是否对于所有$i ot = j$都满足同余方程$c_i+x*p_i equiv c_j+x*p_j pmod m$无解或解比$min(l_i,l_j)$还小,这个方程可以直接转换成$x(p_i-p_j)-my = c_j-c_i$来做
#include<cstdio> #include<cstring> #include<algorithm> using namespace std; const int N=17; int n,m; int c[N],p[N],l[N]; inline int gcd(int a,int b) { return !b?a:gcd(b,a%b); } inline 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 tx,ty;tx=y;ty=x-a/b*y; x=tx;y=ty; return d; } inline bool check() { for(register int i=1;i<=n;i++) { for(register int j=i+1;j<=n;j++) { int ta,tb,tc,t,x,y; ta=p[i]-p[j];tb=m;tc=c[j]-c[i]; t=gcd(ta,tb); if(tc%t)continue; ta/=t;tb/=t;tc/=t;exgcd(ta,tb,x,y); tb=abs(tb);x=((x*tc)%tb+tb)%tb;if(!x)x+=tb; if(x<=min(l[i],l[j]))return 0; } } return 1; } int main() { scanf("%d",&n); for(register int i=1;i<=n;i++)scanf("%d%d%d",&c[i],&p[i],&l[i]),m=max(m,c[i]); while(1) { if(check())break; m++; } printf("%d",m); return 0; }