论求解线性方程

线性方程有很多种,最常见的是二元一次方程。

问题一:

  给定一个一元二次方程ax+by=c,求出一个整数解。

该方程有解的充要条件是gcd(a,b)|c,若不满足此条件,则方程一定无整数解。

可以根据扩展欧几里得求出ax+by=gcd(a,b)的一组特解,再将解乘c/gcd(a,b)即得到一组解。

扩展欧几里得算法:

  设z|x,z|y,根据整除的传递性,z|(x-y),进而可知,当z为gcd(x,y)时,z=gcd(y,x-y),可得gcd(x,y)=gcd(y,x-y)。

  进而可知gcd(x,y)=gcd(y,x-y)=gcd(y,x-2*y)=gcd(y,x-3*y)=……=gcd(y,x%y)。

  gcd(a,b)=gcd(b,a%b)

  p*a+q*b=gcd(a,b)=gcd(b,a%b)=p*b+q*a%b

  a%b=a-a/b*b

  整理可得p*a+q*b=q*a+(p-a/b*q)*b

  在求gcd的过程中递归即可。

模板:

 

 1 #include<iostream>
 2 #include<cstdio>
 3 using namespace std;
 4 long long a,b,c;
 5 long long ex_gcd(long long a,long long b,long long &x,long long &y)
 6 {
 7     if(b==0)
 8     {
 9         x=1;
10         y=0;
11         return a;
12     }
13     long long gcd=ex_gcd(b,a%b,x,y);
14     long long temp=x;
15     x=y;
16     y=temp-a/b*y;
17     return gcd;
18 }
19 int main()
20 {
21     cin>>a>>b;
22     long long x,y;
23     ex_gcd(a,b,x,y);
24     cout<<x<<" "<<y;
25     return 0;
26 }

 

问题二:

  给定同余方程ax≡b(mod p),求其最小正整数解。

该方程可转化为ax+py=b用扩展欧几里得求得特解后,将其转化为最小正整数解。

对于所有的x,其差值最小为p/gcd,设delta=p/gcd,则最小正整数解x0=(x%delta+delta)%delta。

问题三:

  给定方程ax+by=c,求其正整数解的数量。

将c转化为正数,这样比较好处理。 

由解析几何的知识可知,该方程在坐标系中等价于直线y=-a/b*x+c/b。

由此可以判断一些特殊解。

用扩展欧几里得可求出一组解,我们还知道x,y的最小增量dex,dey。

如果方程有有限个解,x有最大和最小值。

x和y的增长关系是相反的,即当x=x-dex时,y=y+dey。当x为最小整数解时,x=minx,当y为最小整数解时,x=maxx。

求得minx和maxx后,两者相减,再除以最小增量,加一即为解的个数。

Code:

 1 #include<iostream>
 2 #include<cstdio>
 3 #define LL long long
 4 using namespace std;
 5 int T,opa,opb;
 6 LL a,b,c;
 7 LL exgcd(LL a,LL b,LL &x,LL &y)
 8 {
 9     if(b==0){
10         x=1;y=0;
11         return a;
12     }
13     LL gcd=exgcd(b,a%b,x,y);
14     LL temp=x;
15     x=y;
16     y=temp-a/b*y;
17     return gcd;
18 }
19 bool judge(LL x,LL y)
20 {
21     if(x<0)    x=-x;
22     if(y<0)    y=-y;
23     LL xx,yy;
24     LL gcd=exgcd(x,y,xx,yy);
25     if(c%gcd==0)    return true;
26     else    return false;
27 }
28 int main()
29 {
30     scanf("%d",&T);
31     while(T--)
32     {
33         scanf("%lld%lld%lld",&a,&b,&c);
34         if(c<0){
35                 a=-a;b=-b;c=-c;
36         }
37         if(a==0&&b==0)    {
38             if(c==0)    cout<<"ZenMeZheMeDuo"<<endl;
39             else    cout<<0<<endl;
40             continue;
41         }    
42         else if(a==0){
43             if(c%b==0&&b>0)    cout<<"ZenMeZheMeDuo"<<endl;
44             else    cout<<0<<endl;
45             continue;
46         }
47         else if(b==0){
48             if(c%a==0&&a>0)    cout<<"ZenMeZheMeDuo"<<endl;
49             else    cout<<0<<endl;
50             continue;
51         }
52         if(b<0&&a<0){
53             cout<<0<<endl;
54             continue;
55         }
56         else if(a*b<0){
57             if(!judge(a,b))    cout<<0<<endl;
58             else    cout<<"ZenMeZheMeDuo"<<endl;
59             continue;
60         }
61         LL x,y;
62         LL gcd=exgcd(a,b,x,y);
63         if(c%gcd!=0){
64             cout<<0<<endl;
65             continue;
66         }
67         LL pow=c/gcd;
68         x=x*pow;y=y*pow;
69         LL dex=b/gcd,dey=a/gcd;
70         LL k=x/dex;
71         x=x-k*dex;y=y+k*dey;
72         if(x<=0){
73             x+=dex;y-=dey;
74         }
75         if(y<=0){
76             cout<<0<<endl;
77             continue;
78         }
79         LL minx=x;
80                 k=y/dey;
81         x=x+k*dex;y=y-k*dey;
82         if(y<=0){
83             x-=dex;y+=dey;
84         }
85         LL maxx=x;
86         LL ans=(maxx-minx)/dex+1;
87         if(ans>65535)    cout<<"ZenMeZheMeDuo"<<endl;
88         else if(ans<=0)    cout<<0<<endl;
89         else    cout<<ans<<endl;
90     }
91     return 0;
92 }                
View Code

 

原文地址:https://www.cnblogs.com/hz-Rockstar/p/11226311.html