中国剩余定理
发布日期:2021-06-29 05:37:46 浏览次数:4 分类:技术文章

本文共 3205 字,大约阅读时间需要 10 分钟。

参考自:

中国剩余定理(CRT)的表述如下

 

设正整数两两互素,则同余方程组

 

                             

 

有整数解。并且在模下的解是唯一的,解为

 

                               

 

其中,而的逆元。

代码:

int CRT(int a[],int m[],int n)  {      int M = 1;      int ans = 0;      for(int i=1; i<=n; i++)          M *= m[i];      for(int i=1; i<=n; i++)      {          int x, y;          int Mi = M / m[i];          extend_Euclid(Mi, m[i], x, y);          ans = (ans + Mi * x * a[i]) % M;      }      if(ans < 0) ans += M;      return ans;  }
习题:
代码:

#include 
#include
#include
using namespace std; int a[4], m[4]; void extend_Euclid(int a, int b, int &x, int &y) { if(b == 0) { x = 1; y = 0; return; } extend_Euclid(b, a % b, x, y); int tmp = x; x = y; y = tmp - (a / b) * y; } int CRT(int a[],int m[],int n) { int M = 1; int ans = 0; for(int i=1; i<=n; i++) M *= m[i]; for(int i=1; i<=n; i++) { int x, y; int Mi = M / m[i]; extend_Euclid(Mi, m[i], x, y); ans = (ans + Mi * x * a[i]) % M; } if(ans < 0) ans += M; return ans; } int main() { int p, e, i, d, t = 1; while(cin>>p>>e>>i>>d) { if(p == -1 && e == -1 && i == -1 && d == -1) break; a[1] = p; a[2] = e; a[3] = i; m[1] = 23; m[2] = 28; m[3] = 33; int ans = CRT(a, m, 3); if(ans <= d) ans += 21252; cout<<"Case "<
<<": the next triple peak occurs in "<

普通的中国剩余定理要求所有的互素,那么如果不互素呢,怎么求解同余方程组?

 

这种情况就采用两两合并的思想,假设要合并如下两个方程

 

      

 

那么得到

 

       

 

在利用扩展欧几里得解出的最小正整数解,再带入

 

       

 

得到后合并为一个方程的结果为

 

       

即新方程为X = A mod M,其中A就是上式的x,M为lcm(m1,m2)。

 

这样一直合并下去,最终可以求得同余方程组的解。

题目:

#include 
#include
#include
using namespace std; typedef long long LL; const int N = 1005; LL a[N], m[N]; LL gcd(LL a,LL b) { return b? gcd(b, a % b) : a; } void extend_Euclid(LL a, LL b, LL &x, LL &y) { if(b == 0) { x = 1; y = 0; return; } extend_Euclid(b, a % b, x, y); LL tmp = x; x = y; y = tmp - (a / b) * y; } LL Inv(LL a, LL b) { LL d = gcd(a, b); if(d != 1) return -1; LL x, y; extend_Euclid(a, b, x, y); return (x % b + b) % b; } bool merge(LL a1, LL m1, LL a2, LL m2, LL &a3, LL &m3) { LL d = gcd(m1, m2); LL c = a2 - a1; if(c % d) return false; c = (c % m2 + m2) % m2; m1 /= d; m2 /= d; c /= d; c *= Inv(m1, m2); c %= m2; c *= m1 * d; c += a1; m3 = m1 * m2 * d; a3 = (c % m3 + m3) % m3; return true; } LL CRT(LL a[], LL m[], int n) { LL a1 = a[1]; LL m1 = m[1]; for(int i=2; i<=n; i++) { LL a2 = a[i]; LL m2 = m[i]; LL m3, a3; if(!merge(a1, m1, a2, m2, a3, m3)) return -1; a1 = a3; m1 = m3; } return (a1 % m1 + m1) % m1; } int main() { int n; while(scanf("%d",&n)!=EOF) { for(int i=1; i<=n; i++) scanf("%I64d%I64d",&m[i], &a[i]); LL ans = CRT(a, m, n); printf("%I64d\n",ans); } return 0; }

转载地址:https://blog.csdn.net/zhj_fly/article/details/76234869 如侵犯您的版权,请留言回复原文章的地址,我们会给您删除此文章,给您带来不便请您谅解!

上一篇:组合数取模
下一篇:HDU5019Revenge of GCD

发表评论

最新留言

哈哈,博客排版真的漂亮呢~
[***.90.31.176]2024年04月19日 13时30分51秒