zoukankan      html  css  js  c++  java
  • 拓展中国剩余定理(exCRT)摘要

    清除一个误区

    虽然中国剩余定理和拓展中国剩余定理只差两个字,但他俩的解法相差十万八千里,所以会不会CRT无所谓

    用途

    求类似$$egin{cases}x equiv b_{1}pmod{a_{1}} \x equiv b_{2}pmod{a_{2}} \...\x equiv b_{n}pmod{a_{n}} \ end{cases}$$的线性同余方程组的解

    具体过程

    假设现在我们只有两个同余方程$$x equiv b_{1}pmod{a_{1}}$$ $$x equiv b_{2}pmod{a_{2}}$$

    改写它们得$$x=a_{1}k_{1}+b_{1}                  ①$$ $$x=a_{2}k_{2}+b_{2}$$

    联立$$a_{1}k_{1}+b_{1}=a_{2}k_{2}+b_{2}$$

    所以$$a_{1}k_{1}=b_{2}-b_{1}+a_{2}k_{2}$$

    记$$c=gcd(a_{1},a_{2})$$ $$d_{1}=frac{a_{1}}{c}$$ $$d_{2}=frac{a_{2}}{c}$$

    两边同除以$c$得$$d_{1}k_{1}=frac{b_{2}-b_{1}}{c}+d_{2}k_{2}$$

    即$$d_{1}k_{1} equiv frac{b_{2}-b_{1}}{c}\%d_{2}pmod{d_{2}}$$

    显然,方程组有解的条件是$c|(b_{2}-b_{1})$

    为什么要在这个时候取模$\%d_{2}$?因为之后取模的时候模数和被模数已经不互质了

    记$d_{1}$在模$d_{2}$意义下的逆元为$p$,两边同乘$p$得$$k_{1} equiv frac{p(b_{2}-b_{1})}{c}\%d_{2} pmod{d_{2}}$$

    即$$k_{1} = frac{p(b_{2}-b_{1})}{c}\%d_{2}+d_{2}*y$$

    带入①式得$$x=frac{p(b_{2}-b_{1})}{c}\%d_{2}*a_{1}+d_{2}a_{1}y+b_{1}$$

    即$$x equiv frac{p(b_{2}-b_{1})}{c}\%d_{2}*a_{1}+b_{1}pmod{a_{1}d_{2}}$$

    注意,上式中的$a_{1}$万万不可乘到$\%d_{2}$前面

    所以我们就弄出了一个新的同余方程。这个方程也可以写成$$x equiv b pmod a$$其中$$b=(inv(frac{a_{1}}{(a_{1},a_{2})},frac{a_{2}}{(a_{1},a_{2})})*frac{b_{2}-b_{1}}{(a_{1},a_{2})}\%frac{a_{2}}{(a_{1},a_{2})}*a_{1}+b_{1})$$ $$a=frac{a_{1}a_{2}}{(a_{1},a_{2})}$$

    其中$inv(m,n)$表示$m$在模$n$意义下的逆元,$(m,n)$表示$m$和$n$的最大公约数

    然后把新求出的方程和后面的方程联立,迭代求解,直到只剩一个方程

    精度问题

    毫无疑问的是,在求解方程组的时候很容易爆$long long$,那该怎么办呢?

    __int128 快(龟)速(速)乘(乘)是个不错的主意。

    何为快速乘?

    就是一种和快速幂类似的东西,利用二进制在log时间内求出两个数相乘对另一个数取模的结果并且不会爆long long

    但要注意的是,千万不要让负数出现在快速乘里面,尤其是被拆分的那个乘数,会死循环的

    1 int mul(int x,int k,int mod) {
    2     ll ans=0;
    3     while(k) {
    4         if(k&1)(ans+=x)%=mod;
    5         k>>=1;
    6         (x<<=1)%=mod;
    7     }
    8     return ans;
    9 }
    快速乘

    题目

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

     1 #include<bits/stdc++.h>
     2 using namespace std;
     3 #define INF 0x7fffffff
     4 #define ME 0x7f
     5 #define FO(s) freopen(s".in","r",stdin);freopen(s".out","w",stdout)
     6 #define fui(i,a,b,c) for(int i=(a);i<=(b);i+=(c))
     7 #define fdi(i,a,b,c) for(int i=(a);i>=(b);i-=(c))
     8 #define fel(i,a) for(register int i=hd[a];i;i=dg[i].nxt)
     9 #define ll long long
    10 #define int long long
    11 #define MEM(a,b) memset(a,b,sizeof(a))
    12 #define maxn 100010
    13 ll n;
    14 ll a[maxn],b[maxn],a0,b0;
    15 template<class T>
    16 inline T read(T &n){
    17     n=0;int t=1;double x=10;char ch;
    18     for(ch=getchar();!isdigit(ch)&&ch!='-';ch=getchar());(ch=='-')?t=-1:n=ch-'0';
    19     for(ch=getchar();isdigit(ch);ch=getchar()) n=n*10+ch-'0';
    20     if(ch=='.') for(ch=getchar();isdigit(ch);ch=getchar()) n+=(ch-'0')/x,x*=10;
    21     return (n*=t);
    22 }
    23 template<class T>
    24 T write(T n){
    25     if(n<0) putchar('-'),n=-n;
    26     if(n>=10) write(n/10);putchar(n%10+'0');return n;
    27 }
    28 template<class T>
    29 T writeln(T n){write(n);putchar('
    ');return n;}
    30 int exgcd(int a,int b,int &x,int &y){
    31     if(!b)return a+((x=1)*(y=0));
    32     int xx,yx,z;z=exgcd(b,a%b,xx,yx);
    33     x=yx,y=xx-a/b*yx;return z;
    34 }
    35 int gcd(int a,int b){if(!b)return a;return gcd(b,a%b);}
    36 int mul(int x,int k,int mod){
    37     ll ans=0;
    38     while(k){
    39         if(k&1)(ans+=x)%=mod;
    40         k>>=1;
    41         (x<<=1)%=mod;
    42     }
    43     return ans;
    44 }
    45 signed main(){
    46     read(n);
    47     fui(i,1,n,1)read(a[i]),read(b[i]);
    48     a0=a[1],b0=b[1];
    49     fui(i,2,n,1){
    50         int c=gcd(a0,a[i]),d1=a0/c,d2=a[i]/c,p,x,y;
    51         int a1=a0,a2=a[i],b1=b0,b2=b[i];
    52         if(((b2-b1)/c*c)^(b2-b1))return 0*writeln(-1);
    53         exgcd(d1,d2,p,y);((p%=d2)+=d2)%=d2;
    54         x=mul(p,(((b2-b1)/c)%d2+d2)%d2,d2)*a1+b1;
    55         y=a1*d2;b0=x,a0=y;((b0%=a0)+=a0)%=a0;
    56     }
    57     int x,y;exgcd(1,a0,x,y);
    58     ((x%=a0)+=a0)%=a0;
    59     x=mul(x,b0,a0);
    60     writeln(x);
    61     return 0;
    62 }
    AC代码

     好像还有NOI2018屠龙勇士?23333先等我去A掉它

    作者:A星际穿越
    我的博客写得这么烂,应该不会有人想转载的吧
    如果要转载的话,请在文章显眼处标明作者和出处谢谢
  • 相关阅读:
    关于Django 报错 ImportError: cannot import name RegexUrlResolver解决
    Git版本控制
    java Labmda表达式
    java注解&反射
    maven
    数据库
    为什么java中子类重写父类的方法时声明抛出异常不能比父类范围大
    【Linux】Linux基本命令
    阿里云服务器linux环境搭建SSM项目(一)--Linux环境配置jdk和Tomcat.md
    阿里云服务器linux环境搭建SSM项目(二)--linux环境配置mysql5.7.md
  • 原文地址:https://www.cnblogs.com/Axjcy/p/9801119.html
Copyright © 2011-2022 走看看