与数论的厮守04:扩展中国剩余定理

  根据数论的尿性,一般的定理解决不了的问题怎么办?那就拓展一下呗。

  我们先看中国剩余定理,它解决的是一堆同余方程:nΞa1(mod b1),nΞa2(mod b2)...nΞak(mod bk),其中b1,b2...bk两两互质。但不互质的时候怎么办呢?不互质就解决不了了,因为bi和其他b1*b1*...*bk/bi中有公因数。所以中国剩余定理就解不了这种问题。

  然后来看看扩展的(虽然两者没什么关系)。先从简单入手,假设现在只有两个方程:nΞa1(mod b1),nΞa2(mod b2),它们可以写成:n=a1+k1*b1,n=a2+k2*b2。也就是:a1+k1*b1=a2+k2*b2。然后就得到一个不定方程:bi*k1-b2*k2=a2-a1。于是我们可以用扩欧来求出k1,而k1是刚才关于n的方程里的一个元素,我们当然希望能够通过一个式子得到所有的k1。所以我们还要求出k1的通项。设j1,j2为不定方程b1*j1-b2*j2=gcd(b1,b2)的一组特解,两边为了变成k1,k2的方程的形式,设t=(a2-a1)/gcd(b1,b2),把两边同时乘上t,得到bi*j1*t-b2*j2*t=a2-a1,就得到了k1=j1*t,k2=j2*t。于是我们进一步得到了n的通项:n=a1+b1*j1*t,n=a2+b2*j2*t。于是,每两个解的差值Δ|b1*t,Δ|b2*t。通过观察,两个式子的右边有不同的项:b1,b2,要想能够整除,则Δ必须整除b1,b2;把t=(a2-a1)/gcd(b1,b2)带进去,因为b1*b2/gcd(b1,b2)是b1,b2的最小公倍数,所以能被b1,b2整除,所以能被Δ整除。也就是:Δ|lcm(b1,b2)。然后我们就求出了n的通项:n=b1*k1+x1+k*lcm(b1,b2),就可以把另一个等式扔掉了。再转化为同余方程:nΞbi*k1+x1(mod lcm(b1,b2)),每次求出k1即可。也就是说,我们通过上面的转化,把两个同余方程合并成了一个。

  再推广一下,如果有m个同余方程。。。那一个个合并不就好了吗?时间复杂度就是O(mlog),由于每次的不定方程不同,所以log函数的自变量也就不同,这里不好说。凑合着看吧。


  下面是喜闻乐见的模板题(要么打龟速乘,要么开int128):

  P4777 【模板】扩展中国剩余定理(EXCRT)

  以及代码(丑得一批):

#include <cstdio>
#pragma GCC diagnostic error "-std=c++14"
#pragma GCC target("avx")
#pragma GCC optimize(3)
#pragma GCC optimize("Ofast")
#define maxn 100001
using namespace std;
typedef __int128 ll;

inline ll read(){
    register ll x(0), f=1; register char c(getchar());
    while(c<'0'||'9'<c){ if(c=='-') f=-1; c=getchar(); }
    while('0'<=c&&c<='9')
        x=(x<<1)+(x<<3)+(c^48), c=getchar();
    return x*f;
}

ll a[maxn], b[maxn];

void ex_gcd(ll a, ll b, ll &x, ll &y, ll &d){
    if(!b) x=1,y=0,d=a;
    else{
        ll x1,y1;
        ex_gcd(b,a%b,x1,y1,d);
        x=y1,y=x1-a/b*y1;
    }
}

int main(){
    ll n=read();
    for(register ll i=1;i<=n;i++) b[i]=read(),a[i]=read();
    ll lcm=b[1],prex=a[1];
    for(register ll i=2;i<=n;i++){
        ll lcm_a=((a[i]-prex)%b[i]+b[i])%b[i],x,y,k=lcm,gcd;
        ex_gcd(lcm,b[i],x,y,gcd);
        x=(x*lcm_a/gcd%(b[i]/gcd)+(b[i]/gcd))%(b[i]/gcd);
        lcm=lcm*b[i]/gcd,prex=(prex+k*x)%lcm;
    }
    printf("%lld\n",(long long)((prex%lcm+lcm)%lcm));
    return 0;
}

  学完这个终于可以去搞扩展lucas了(噗)

原文地址:https://www.cnblogs.com/akura/p/10730975.html