扩展中国剩余定理 exCRT 学习笔记

前言

由于 ({mathrm{CRT}}subseteq{mathrm{exCRT}}),而且 CRT 又太抽象了,所以直接学 exCRT 了。

摘自 huyufeifei 博客

lbwtxdy

这么抽象的东西我怎么可能会写

前置技能

  • gcd/lcm
  • exgcd
  • 快速乘

参考资料

用途

用于求一个关于 (x​) 的同余方程组

[left{ egin{matrix} xequiv a_1pmod{b_1}\ xequiv a_2pmod{b_2}\ cdots\ xequiv a_npmod{b_n} end{matrix} ight. ]

的解。

解法推导

考虑到如果 (x) 对原方程组成立,那么 (x) 对其中任意几个方程也都成立。那么如果要满足 (x) 对前 (i) 个方程都成立,一个必要条件是对前 (i-1) 个方程都成立。

采用一种类似“增量法”的思路来合并每个式子。

假定前 (i-1) 个式子已经合并完了,得到 (xequiv a_{i-1}pmod{b_{i-1}})。此时合并

[left{ egin{aligned} &xequiv a_{i-1}&pmod{b_{i-1}}\ &xequiv a_i&pmod{b_i} end{aligned} ight. ]

它的本质是存在 (k_1,k_2​),满足(上下颠倒了一下,方便下面推导)

[left{ egin{aligned} &x+k_1b_i=a_i\ &x+k_2b_{i-1}=a_{i-1} end{aligned} ight. ]

两式相减,得到

[k_1b_i-k_2b_{i-1}=a_i-a_{i-1} ]

此时,仅有 (k_1,k_2) 是变量,剩下的都已知。但是我们用 exgcd 解方程时,要求等式右边是 (gcd(b_i,b_{i-1}))。此时我们列出一个新方程

[k_1'b_i+k_2'b_{i-1}=gcd(b_i,b_{i-1}) ]

解出 (k_1'​) 的一个特解 (k_0​)

那么为了满足新方程和原方程都成立,那么

[frac {k_1}{k_1'}=frac{k_2}{k_2'}=frac{a_i-a_{i-1}}{gcd(b_i,b_{i-1})} ]

因此 (k_1) 有特解

[k_{1_0}=k_0 imes frac{a_i-a_{i-1}}{gcd(b_i,b_{i-1})} ]

(k_1) 的通解是

[k_{1_0}+t imes frac{b_{i-1}}{gcd(b_i,b_{i-1})}, tin mathbb Z ]

不定方程的通解

对于 (ax+by=c​) 这个方程,假定我们已经解出来一组解为

[left{ egin{matrix} x=x_0\ y=y_0 end{matrix} ight. ]

我们带入后对左侧式子变形,得到 (ax_0+S+by_0-S=c​)

我们要从 (ax_0+S​) 中提一个 (a​) 出来,从 (by_0-S​) 中提一个 (b​) 出来,得到 (a(x_0+frac Sa)+b(y_0-frac Sb)=c​)

此时最小的满足 (x_0+frac Sa, y_0-frac Sb​) 都是整数的 (S​) 就是 (a,b​) 的最小公倍数了。

那么每个相邻的 (x​) 之间相差 (frac Sa=frac{operatorname{lcm}(a,b)}a=frac b{gcd(a,b)}​),即 (x​) 的通解为 (x_0+t imes frac b{gcd(a,b)}, tin mathbb Z​)

后面的式子仅与 (a,b​) 有关。

又因为

[x+k_1b_i=a_i ]

我们带入 (k_1) 的通解,得

[x+left(k_{1_0}+t imes frac{b_{i-1}}{gcd(b_i,b_{i-1})} ight) imes b_i=a_i, tin mathbb Z ]

化简

[x+b_ik_{1_0}+t imes frac{b_ib_{i-1}}{gcd(b_i,b_{i-1})}=a_i, tin mathbb Z\ t imesoperatorname{lcm}(b_i,b_{i-1})=a_i-x-b_ik_{1_0}, tin mathbb Z ]

由于 (tin mathbb Z),所以

[a_i-x-b_ik_{1_0}|operatorname{lcm}(b_i,b_{i-1})\ a_i-x-b_ik_{1_0}equiv 0pmod{operatorname{lcm}(b_i,b_{i-1})}\ xequiv a_i-b_ik_{1_0}pmod{operatorname{lcm}(b_i,b_{i-1})} ]

(equiv) 符号右边的都是已知量,这时我们就把两个方程结合到一块了。

我们依次合并 ((2)=[(1),(2)],(3)=[(2),(3)],cdots,(n)=[(n-1),(n)]),最终得到第 (n) 个式子

[xequiv a_npmod{operatorname{lcm}(b_1,b_2,cdots,b_m)} ]

(a_n) 是合并 (n-1) 个式子后的新 (a_n),上述过程中每一步都会更新 (a_i​)

此时 (a_n) 就是解了,取模取正数就是最小正整数解。

(上下两部分实质相同)

解法简述

对于

[left{ egin{aligned} &xequiv a_{i-1}&pmod{b_{i-1}}\ &xequiv a_i&pmod{b_i} end{aligned} ight. ]

[left{ egin{aligned} &x+k_1b_i=a_i\ &x+k_2b_{i-1}=a_{i-1} end{aligned} ight. ]

解出

[k_1'b_i-k_2'b_{i-1}=gcd(b_i,b_{i-1}) ]

的特解 (k_1'​),两边同乘 (frac{a_i-a_{i-1}}{gcd(b_i,b_{i-1})}​),得到

[k_1b_i-k_2b_{i-1}=a_i-a_{i-1} ]

那么 (k_ 1') 也被乘了 (frac{a_i-a_{i-1}}{gcd(b_i,b_{i-1})}),因此 (k_1) 有特解

[k_{1_0}=k_1' imes frac{a_i-a_{i-1}}{gcd(b_i,b_{i-1})} ]

(k_1) 的通解是

[k_1=k_{1_0}+t imesfrac{b_{i-1}}{gcd(b_i,b_{i-1})} ]

代回去就是

[x+k_{1_0}b_i+t imes operatorname{lcm}(b_i,b_{i-1})=a_i ]

得到

[xequiv a_i-k_{1_0}b_ipmod{operatorname{lcm}(b_i,b_{i-1})} ]

作为第 (i) 个方程即可。

注意事项

(k_1,k_{1_0},a_i) 的计算中是有可能爆 long long 的,而这三个计算恰好又都是模意义下,所以使用快速乘。

一般情况下不考虑 (operatorname{lcm}(b_1,b_2,cdots,b_n))long long 的情况。

代码

#include<cstdio>
#define ll long long
ll exgcd(ll a,ll b,ll &x,ll &y)
{
    if(!b)
    {
        x=1,y=0;
        return a;
    }
    ll g=exgcd(b,a%b,y,x);
    y-=a/b*x;
    return g;
}
ll qmul(ll x,ll y,ll p)
{
    ll ans=0,f=1;
    if(x<0)
    {
        x=-x;
        f=-f;
    }
    if(y<0)
    {
        y=-y;
        f=-f;
    }
    while(y)
    {
        if(y&1)
            ans=(ans+x)%p;
        x=(x+x)%p;
        y>>=1;
    }
    return ans*f;
}
ll a[100100],b[100100];
int main()
{
    int n;
    ll x,y;
    scanf("%d",&n);
    for(int i=1;i<=n;++i)
        scanf("%lld%lld",&b[i],&a[i]);
    for(int i=2;i<=n;++i)
    {
        ll g=exgcd(b[i],b[i-1],x,y);
        ll t=b[i-1]/g;//k1 的最小波动 Δ
        x=qmul(x,(a[i]-a[i-1])/g,t);
        x=(x%t+t)%t;
        a[i]-=qmul(x,b[i],b[i]/g*b[i-1]);
        b[i]=b[i]/g*b[i-1];
        a[i]=(a[i]%b[i]+b[i])%b[i];
    }
    printf("%lld
",(a[n]%b[n]+b[n])%b[n]);
    return 0;
}
原文地址:https://www.cnblogs.com/wjyyy/p/excrt.html