zoukankan      html  css  js  c++  java
  • 浅析中国剩余定理(从CRT到EXCRT))

    前置知识

    1.

    a%b=d,c%b=e,

    则(a+c)%b=(d+e)%b(正确性在此不加证明)

    2.

    a%b=1,则(d( imes)a)%b=d%b(正确性在此不加证明)

    下面先看一道题(改编自曹冲养猪):

    烤绿鸟的故事

    题目描述:

    mian包是一个贪吃的孩子,这天,他买了一堆绿鸟吃。当然他的妈妈并不想让他吃太多食物(因为那样会发胖),为了避免老妈的唠叨,他决定不告诉他的妈妈绿鸟数量,而是将绿鸟的数量x用以下式子来描述

    [egin{cases}x≡b_1 (mod a_1)\x≡b_2 (mod a_2)\...\x≡b_n (mod a_n)end{cases} ]

    (a_1),(a_2)......(a_n)两两互质)
    由于他的妈妈数学不好,于是就来向你求助了,请你求出mian包最少买了多少烤绿鸟

    输入输出格式

    输入格式:

    第一行包含一个整数n (n <= 10) – 告诉妈妈的式子数,接下来n行,每行两个整数ai, bi( bi <= ai <= 1000)

    首先,我们把这道题简化成下面的图

    样例为:

    3
    3 1
    5 1
    7 2
    

    很显眼吧?怎么做呢?

    1.妈妈,我会暴力

    这就很简单了

    我们从第二个开始枚举

    初始值为a1+b1;

    每次自加当前的已经满z足的(a_1)**(a_2)......(a_{i-1})的乘积(为什么我接下来说)

    如果当前值已经满足%(a_i)=(b_i)

    则当前已经成立

    ps:上式为什么成立呢?

    以上面的例子为例,看下面我的图

    当已经满足前两个等式时,我们增加(a_1)*(a_2)

    那么我们会增加15只绿鸟

    看图:

    很显然,15%3=0,15%5=0

    所以增加15的话一定满足当前的等式

    自然,暴力就解决了

    暴力代码:

    #include<iostream>
    #include<cstdio>
    #include<algorithm>
    using namespace std;
    long long a,b,c,d,e,f,g,kkk,n;
    struct zzzz{
        int x,y,z;
    }ltt[15];
    bool cmp(zzzz s,zzzz t)
    {
        return s.x>t.x;
    }
    int main()
    {
        cin>>n;
        kkk=1;
        for(a=1;a<=n;a++)
        {
            cin>>ltt[a].x>>ltt[a].y;
        }
        sort(ltt+1,ltt+n+1,cmp);
        kkk=ltt[1].x;
        g=ltt[1].y;
        f+=2;
        while(f<=n)//枚举每个等式
        {
            while(1)
            {
                if(g%ltt[f].x==ltt[f].y)
                {
                    kkk=kkk*ltt[f].x;
                    f++;
                    break;
                }
                g=g+kkk;//自加,直到满足当前等式
            }
        }
        cout<<g;
    }
    

    诚然,这种方法好想也好写,但是,我们要注意一个问题

    看下面这个式子:

    [egin{cases}x≡1 (mod 3)\x≡1 (mod 5)\x≡19260817 (mod 1000000007)end{cases} ]

    好吧,你自加去吧,保证TLE

    这时,我们就需要一个高效的算法

    前置知识3:

    每一个crt方程组的解在模(a_1 imes a_2 imes …… imes a_n)意义下有且只有1个

    为什么呢?因为每个解之间的差是所有a的乘积,加数不管对序列中的哪个数取模都是0呗

    2.关于暴力的扩展

    先说点别的东西

    由上面的暴力我们可以看出,该样例crt解法的本质是从5和7的公倍数中找一个除以3余1的数(n_1),从3和7的公倍数中找一个除以5余1的数(n_2),从3和5的公倍数中找一个除以7余2的数(n_3),再将三个数相加得到解。

    那么我们可不可以说,更多项的答案也可以这样的得来吗?

    答案是肯定的

    我们设{(a_1 imes a_2 imes ...... imes a_n)}为mul

    然后按照暴力的方法,我们可以得到这样一个式子:

    (k_1 imesfrac{(mul)}{a1})+(k_2 imesfrac{(mul)}{a2})+……+(k_n imesfrac{(mul)}{an})

    并满足((k_i imesfrac{(mul)}{a_i})%(a_i)=(b_i))

    正确性?看下面:

    首先,因为(a_i)两两互质,所以我们可以知道,(frac{(mul)}{a_i})一定不是(a_i)的倍数

    (frac{(mul)}{a_j})(j( e)i)一定是({a_i})的倍数

    所以,我们知道上式中除了第 (i) 项其他所有项 ≡ 0(mod({a_i})

    所以,原式 ≡ (k_i imesfrac{(mul)}{a_i})(mod (a_i))

    然后,我们又根据一开始那个式子的约束条件,得出(k_i imesfrac{(mul)}{a_i})%(a_i)=(b_i)

    所以,整个式子满足对(a_i)的约束条件

    我们将这个推广开来,即可证得原式是crt方程组的一个解

    然后根据前置知识三,即求得原始的最小解

    你说了这么多,但该怎么实现呢?

    我们设(frac{mul}{a_i})(m_i)

    那么,原式就变成了:

    (k_1 imes m_1 + k_2 imes m_2+......+k_n imes m_n)

    (k_i imes m_i)(b_i) (mod (a_i))

    ** 这是什么?**

    很多人表示一脸懵逼,不知道该怎么求???

    好,让我们仔细看一下

    我们由前置定理2可得:

    如果(k imes m_i) ≡ 1 (mod (a_i))

    那么当(k_i=b_i imes k)

    (k_i imes m_i)(b_i) (mod (a_i))

    那么怎么求 (k) 呢?

    我们知道有一个东西叫逆元。。

    什么是逆元?

    请出门左转以前的一篇日报

    浅谈模质数意义下的乘法逆元

    这里,我们就可以用逆元来解k咯

    代码里我用的是exgcd哦

    代码是重点:

    #include<iostream>
    #include<cstdio>
    #include<cstring>
    #include<algorithm>
    #define rii register int i
    #define rij register int j
    #define int long long
    using namespace std;
    int a[15],b[15],n;
    void exgcd(int a,int b,int &ls,int &x,int &y)//求逆元
    {
        if(b==0)
    	{
            ls=a;
            x=1;
    		y=0;
        }
        else
    	{
            exgcd(b,a%b,ls,y,x);
            y-=(a/b)*x;
        }
    }
    signed main()
    {
    	int mod=1,ls,y,x=0;
        scanf("%lld",&n);
        for(rii=1;i<=n;i++)
        {
    		scanf("%lld%lld",&a[i],&b[i]);  	
        }
        for(rii=1;i<=n;i++)
    	{
    		mod*=a[i];
    	}
        for(rii=1;i<=n;i++)
    	{
            int kkk=mod/a[i];
            exgcd(a[i],kkk,ls,ls,y);
            x=(x+y*kkk*b[i])%mod;
        }
    	int ans=(x+mod)%mod;
        printf("%lld",ans);
    } 
    

    毒瘤 Imagine(模板excrt出题人):我的(a_i)不互质!

    好吧,我们进入下一步,

    excrt

    难的不是一点半点,好像和上一个已经没什么关系了

    因为大家注意一点,上面理论的基础是(a_i)互质

    那么这里这里不互质了我们怎么办呢?

    (说暴力的站墙角)

    但我们又没有好的解法......

    看来,我们不得不学着暴力两两合并了

    当然,我们不能像暴力一样通过自加两两合并

    我们要用一个好玩的东西——(exgcd)

    对于原方程的两个式子,我们可以化成这样(请跟着我的公式):

    [egin{cases}x≡b_1 (mod a_1)\x≡b_2 (mod a_2)end{cases} ]

    (a_1 imes y+b_1=a_2 imes z+b_2)

    这就是一眼看穿的问题了

    所以,这里就能用上面exgcd了

    顺次两两合并

    每次得出一个当前的解

    合并很简单,

    我们首先利用exgcd求出一个解

    然后利用这个解,新构建一个新的方程带入计算

    怎么构建新的解呢?

    我们目前有两个方程,一个是由前面的方程合并来的,我们称它为方程1

    另一个是我们这一步要合并的方程,称为方程2

    [egin{cases}x≡b (mod a)---1\x≡b_i (mod a_i)---2end{cases} ]

    我们利用(exgcd)可以求出这里的x

    现在的关键就在于如何合并(a)(a_i)

    我们知道,(x)+(lcm)(a),(a_i))仍然满足上面的这个方程组

    (因为(lcm)(a),(a_i))%(a)=0,(lcm)(a),(a_i))%(a_i)=0)

    我们就可以得出一个新的方程了

    代码:

    #include<iostream>
    #include<cstdio>
    #include<cstring>
    #include<algorithm>
    #define rii register int i
    #define rij register int j
    #define int long long
    using namespace std;
    int a[100005],b[100005],n,cj,ans,mod;
    int exgcd(int a,int b,int &x,int &y)
    {	
        if(b==0)
        {
            x=1;
            y=0;
            return a;
        }
        int gcd=exgcd(b,a%b,x,y);
        int zcq=x;
        x=y;
        y=zcq-a/b*y;
        return gcd;
    }
    int pw(int l,int r,int p)
    {
        int an=0;
        while(r>0)
        {
            if(r&1)
            {
                an=(an+l)%p;
            }
            l*=2;
            l%=p;
            r/=2;
        }
        return an;
    }
    signed main()
    {
        scanf("%lld",&n);
        for(rii=1;i<=n;i++)
        {
            scanf("%lld%lld",&a[i],&b[i]);
        }
        cj=a[1];
        ans=b[1];
        int x,y;
        for(rii=2;i<=n;i++)
        {
            int zcq=cj;
            int st=a[i];
            int ltt=(b[i]-ans%st+st)%st;
            int gys=exgcd(zcq,st,x,y);//求解
            int kkk=st/gys;
            x=pw(x,ltt/gys,kkk);
            ans+=cj*x;
            cj*=kkk;
            ans+=cj;
            ans%=cj;//防溢出(同时确保答案最小)
        }
        cout<<ans;
    }
    

    会了后别忘了写NOI2018屠龙勇士哦~

    特别感谢luogu曹冲养猪和excrt题解作者,他们给了我很多思路

    如果出锅,请私信联系(评论我不一定看的到),万分感谢

  • 相关阅读:
    ES7 cat API的小结
    zabbix5.0 使用elasticsearch7.6按日期索引存储历史数据
    Archlinux爬坑指南
    ArchLinux安装常用软件QQ、TIM、微信等常用软件(三)
    ArcnLinux安装KDE桌面环境(二)
    ArchLinux安装步骤(一)
    DDD领域驱动及落地方案
    Text Classification with Keras
    Mattermost Server安装及配置AD/LADP
    使用队列问题
  • 原文地址:https://www.cnblogs.com/ztz11/p/9681123.html
Copyright © 2011-2022 走看看