求解模线性方程C语言,求解模线性方程

求解模线性方程组

求解线性不定方程组

ax + by = c

先求出一组解， 然后考虑如何表示通解， 设d = gcd(a, b), 假设c不是d的倍数， 则左边是d的倍数而右边不是， 则方程无解， 所以方程有解当且仅当d | c.

设c = c’ * d, 我们先考虑方程 ax + by = d, 这样由扩展gcd便可求出一组解 (x’, y’), 则(c’x’, c’y’)就是原方程的一组解，然后考虑通解：

假设有两组解(x1, y2) , (x2, y2), 有 ax1 + by1 == ax2 + by2 = c, 移项得： a(x1 - x2) == b(y2 - y1), 消去d后有 a’(x1 - x2) == b’(y2 - y1),

此时a’ 和 b’ 互素， 所以(x1 - x2)一定是b’的倍数， 而(y2 - y1)一定是a’的倍数， 由此可得到通解：给一组特解(x, y), 通解为(x - kb’, y + ka’).

求解模线性方程

ax = b(mod n)

其实方程等价于 ax - ny = b, 标准模线性方程，但是得考虑剩余系。

算法导论上有两个定理：

定理一：设d = gcd(a, n), 假定对整数x’, y’, 有d = ax’ + ny’， 如果d | b, 则方程ax = b(mod)有一个解的值为x0, 满足：、

x0 = x'(b / d)(mod n)

定理二：假设方程ax = b(mod n)有解， x0是方程的任意一个解， 则方程对模n恰有d个不同的解， 分别为： xi = x0 + i * (n / d), 其中 i = 1,2,3……d - 1

有了这两个定理， 解方程就不难了。

代码实现

扩展欧几里得算法

ll exgcd(ll a,ll b,ll &x,ll &y)

{

if(b==0)

{

x=1;

y=0;

return a;

}

ll r=exgcd(b,a%b,x,y);

ll t=x;

x=y;

y=t-a/b*y;

return r;

}

求解模线性方程求的的只是一组解给一组特解(x, y), 通解为(x - kb’, y + ka’).改代码返回满足x>0&&y>0的最小解。

bool linear_equation(ll a,ll b,ll c,ll &x,ll &y)

{

ll d=exgcd(a,b,x,y);

if(c%d)

return false;

ll k=c/d;

x*=k; y*=k; //求得的只是其中一组解

ll b1 = b/d;

ll a1 = a/d;

ll i = 0;

//cout<

if(y<0){

while (y<=0){

y+=a1;

x-=b1;

}

}

while (y-a1>=0){

y-=a1;

x+=b1;

}

if(y>=0&&x>=0){

return true;

} else{

return false;

}

}

中国剩余定理求解模线性方程组

互质版本

中国剩余定理是中国古代求解一次同余方程组的方法，是数论中的一个重要定理。

设m1，m2，m3，…，mk是两两互素的正整数，即gcd(mi,mj)=1，i!=j，i,j=1,2,3,…,k.

则同余方程组：

x = a1 (mod n1)

x = a2 (mod n2)

…

x = ak (mod nk)

模[n1,n2,…nk]有唯一解，即在[n1,n2,…,nk]的意义下，存在唯一的x，满足：

x = ai mod [n1,n2,…,nk], i=1,2,3,…,k。

解可以写为这种形式：

x = sigma(ai* mi*mi’) mod(N)

其中N=n1*n2*…*nk，mi=N/ni，mi’为mi在模ni乘法下的逆元。

非互质版本采用方程合并的思想

//

// Created by 王若璇 on 16/7/11.

//

//中国剩余定理,非互质版本

#include

#include

#include

#include

using namespace std;

typedef long long ll;

ll gcd(ll a,ll b){

return b==0?a:gcd(b,a%b);

}

ll exgcd(ll a,ll b,ll &x,ll &y)

{

if(b==0)

{

x=1;

y=0;

return a;

}

ll r=exgcd(b,a%b,x,y);

ll t=x;

x=y;

y=t-a/b*y;

return r;

}

ll cal(ll a,ll m){

ll x,y;

ll r = exgcd(a,m,x,y);

if(1%r!=0){

return -1;

}

//x*=1/r;

//m = abs(m);

ll ans = (x%m+m)%m;

//if(ans<=0) ans+=m;

return ans;

}

bool merge(ll a1,ll r1,ll a2,ll r2,ll &aa,ll &rr){

ll d = gcd(a1,a2);

ll c = r2-r1;

if(c%d!=0){

return false;

}

c = (c%a2+a2)%a2;

c/=d;

a1/=d;

a2/=d;

c*=cal(a1,a2);

c = c%a2;

c*=a1*d;

c+=r1;

aa = a1*a2*d;

rr = c;

return true;

}

int main(){

ios::sync_with_stdio(false);

cin.tie(NULL);

int n;

while (cin>>n){

ll k,m;

cin>>k>>m;

bool flag = true;

ll a,r,aa,rr;

for(int i = 1;i

cin>>a>>r;

if(flag&&merge(k,m,a,r,aa,rr)){

k = aa;

m = rr;

} else{

flag = false;

}

}

if(flag){

cout<

} else{

cout<

}

}

return 0;

}