前往小程序,Get更优阅读体验!
立即前往
首页
学习
活动
专区
工具
TVP
发布
社区首页 >专栏 >数学--数论--中国剩余定理+扩展中国剩余定理(孙子定理)

数学--数论--中国剩余定理+扩展中国剩余定理(孙子定理)

作者头像
风骨散人Chiam
发布2020-10-28 11:57:30
6380
发布2020-10-28 11:57:30
举报
文章被收录于专栏:CSDN旧文

中国剩余定理

问题

求解同余方程组

在这里插入图片描述
在这里插入图片描述

其中 m 1 , m 2 , m 3 . . . m k m_1,m_2,m_3...m_k m1​,m2​,m3​...mk​为两两互质的整数 求x的最小非负整数解

定理 令 M = ∏ i = 1 k m i M=\prod_{i=1}^km_i M=∏i=1k​mi​,即M是所有 m i m_i mi​的最小公倍数 t i t_i ti​为同余方程 M t / m i ≡ 1 ( m o d   m i ) M_t/m_i ≡ 1(mod\ m_i ) Mt​/mi​≡1(mod mi​)的最小非负整数解 则有一个解为 x = ∑ i = 1 k a i M m i t i x=\sum_{i=1}^ka_i\frac{M}{m_i}t_i x=∑i=1k​ai​mi​M​ti​ 通解为 x + i ∗ M ( i ∈ Z ) x+i*M(i\in Z) x+i∗M(i∈Z) 特别的,最小非负整数解为(x%M+M)%M 证明

在这里插入图片描述
在这里插入图片描述

代码实现 同余式的求解可以用扩展欧几里得

代码语言:javascript
复制
void exgcd(int a,int b,int &x,int &y)
{
    if(b==0){ x=1; y=0; return;}
    exgcd(b,a%b,x,y);
    int tp=x;
    x=y; y=tp-a/b*y;
}

int china()
{
    int ans=0,lcm=1,x,y;
    for(int i=1;i<=k;++i) lcm*=b[i];
    for(int i=1;i<=k;++i)
    {
        int tp=lcm/b[i];
        exgcd(tp,b[i],x,y);
        x=(x%b[i]+b[i])%b[i];//x要为最小非负整数解
        ans=(ans+tp*x*a[i])%lcm;
    }
    return (ans+lcm)%lcm;
}

扩展中国剩余定理

求解同余方程组

在这里插入图片描述
在这里插入图片描述

其中 m 1 , m 2 , m 3 . . . m k m_1,m_2,m_3...m_k m1​,m2​,m3​...mk​不一定为两两互质的整数 求x的最小非负整数解 求解: 求解 假设已经求出前k-1个方程组成的同余方程组的一个解为x 且有 M = ∏ i − 1 k − 1 m i M=\prod_{i-1}^{k-1}m_i M=∏i−1k−1​mi​ 则前k-1个方程的方程组通解为x+i∗M(i∈Z)x+i*M(i\in Z)x+i∗M(i∈Z)

那么对于加入第k个方程后的方程组 我们就是要求一个正整数t,使得 x + t ∗ M ≡ a k ( m o d   m k ) x+t∗M≡a_k(mod \ m_k) x+t∗M≡ak​(mod mk​)

转化一下上述式子得 t ∗ M ≡ a k − x ( m o d   m k ) t∗M≡a_k−x(mod\ m_k ) t∗M≡ak​−x(mod mk​)

对于这个式子我们已经可以通过扩展欧几里得求解t 若该同余式无解,则整个方程组无解 若有,则前k个同余式组成的方程组的一个解解为 x k = x + t ∗ M x_k=x+t*M xk​=x+t∗M

所以整个算法的思路就是求解k次扩展欧几里得

代码语言:javascript
复制
#include<iostream>
using namespace std;
#define LL long long
LL mi[1100],ai[1100];//mi为要模的数,ai为余数。
LL gcd(LL a, LL b)
{
    return b == 0 ? a : gcd(b, a%b);
}
void exgcd(LL a, LL b, LL &d, LL &x, LL &y)
{
    if(!b)
    {
        d = a, x = 1, y = 0;
    }
    else
    {
        exgcd(b, a%b, d, y, x);
        y -= x * (a / b);
    }
}
LL CRT(LL l, LL r, LL *mi, LL *ai)
{
    LL lcm = 1;
    for(LL i = l; i <= r; i++)
        lcm = lcm / gcd(lcm, mi[i]) * mi[i];
    for(LL i = l+1; i <= r; i++)
    {
        LL A = mi[l], B = mi[i], d, x, y, c = ai[i] - ai[l];
        exgcd(A, B, d, x, y);
        if(c % d)
            return -1;
        LL mod = mi[i] / d;
        LL k = ((x * c / d) % mod + mod) % mod;
        ai[l] = mi[l] * k + ai[l];
        mi[l] = mi[l] * mi[i] / d;
    }
    /*if(ai[l] == 0)
        return lcm;*/ //保证结果为正整数
    return ai[l];
}
int main()
{
    LL t,n,i,aa;
    while(cin>>n>>aa)
    {
		if(n==0) break; 
        for(i=1; i<=n; i++){
            cin>>mi[i];
			ai[i]=mi[i]-aa;
		}
        cout<<CRT(1ll,n,mi,ai)<<endl;
    }
    return 0;
}
本文参与 腾讯云自媒体同步曝光计划,分享自作者个人站点/博客。
原始发表:2019/12/06 ,如有侵权请联系 cloudcommunity@tencent.com 删除

本文分享自 作者个人站点/博客 前往查看

如有侵权,请联系 cloudcommunity@tencent.com 删除。

本文参与 腾讯云自媒体同步曝光计划  ,欢迎热爱写作的你一起参与!

评论
登录后参与评论
0 条评论
热度
最新
推荐阅读
目录
  • 中国剩余定理
    • 问题
    • 扩展中国剩余定理
    领券
    问题归档专栏文章快讯文章归档关键词归档开发者手册归档开发者手册 Section 归档