从模运算到扩展欧几里得算法,解线性不定方程,同余方程和逆元

先导知识

约数(因数)

如果一个数$A$能被另外一个数$B$整除$A\%B=0$
那么称$B$为$A$的约数,$A$为$B$的倍数

质数(素数)

只能被$1$和其本身整除的数

质约数(质因数,素因数)

能够整出正整数$A$的质数,即约数中为质数的数。

最大公因数(最大公因子,最大公约数)

两个数的公共约数中最大的那个。
特别的,如果两个数的公约数只有$1$,那么这两个数叫做互质整数

最大公约数和最小公倍数的关系

最大公因数的一般求法

最大公因数的一般求法即定义法,从两个数中最小的那个数开始枚举到$1$,如果枚举到的那个数能够同时整除两个数,则终止循环输出那个数。
代码实现如下

#include <bits/stdc++.h>
using namespace std;
int main()
{
    int a,b;
    cin>>a>>b;
    for(int i=min(a,b);i>=1;i--)
    {
        if(a%i==0&&b%i==0){
        cout<<i<<endl;
        break;}
        
    }
    return 0; 
}

但是由于这种算法的时间复杂度很大,所以我们要用到下面要讲到的欧几里得算法。

欧几里得算法(Euclidean algorithm)

欧几里得算法又称辗转相除法,表示为
$$GCD(a,b)=GCD(b,a\%b)$$

欧几里得算法的证明

设$l$为$A$,$B$的最大公因数
$$L=GCD(A,B)$$
则A和B可以表示为:
$$A=Q*L$$
$$B=W*L$$
则有
$$GCD(A,B)=GCD(Q*L,W*L)$$
因为$W$和$Q\%W$互质,所以则有
$$GCD(A,B)=GCD(Q*L,W*L)=GCD(W*L,(Q\%W)*L)=GCD(B,A\%B)$$

当然如果看不懂也没关系,毕竟我们研究的是计算机科学,简单的来说,我们只需要知道不断交换两个数,在$B=0$时,$A$的值就是$A$和$B$的最大公倍数了,通过$C++$的递归可以很简洁的实现。

int gcd(int a,int b)
{
    if(b==0)return 0;
    else gcd(b,a%b)
}

扩展欧几里得算法(Extended Euclidean algorithm)

扩展欧几里得算法可以用来求解二元线性不定方程,在介绍扩展欧几里得算法之前我们先了解几个概念。

裴(pei)蜀(贝祖)定理(Bézout's lemma)

对于任何整数$a$,$b$和他们的最大公约数$d$,关于未知数$x$和$y$的线性不定方程为称为裴蜀等式:
$$ax+by=m$$
使裴蜀等式有解的定理被称为裴蜀定理:
$ax+by=m$有整数解时,当且仅当$m$是$d$的倍数,且有解时必然无穷个整数解,每组解$x$,$y$都被称为裴蜀数。
简单的来说,要使得$$ax+by=m$$有整数解,$m$必须是$a$和$b$的最大公约数的倍数。并且一旦有解那么就有无数个解。
例如:
$12x+42y=6$有解,因为$GCD(42,12)=6$,而$6$是$m$的倍数,我们可以找到其中的一个解,例如:
$$x=4,y=-1$$
只要我们稍加改动,就可以得到一个新的解,例如:
$$x=11,y=-3$$
那么我们该如何改动呢?答案很简单,我们只需要先使用最大公倍数$d$进行约分:
$$2x+7y=1$$
每当$x$增加$y$的系数的值,那么$y$就应该降低$x$的系数值,例如$x$上升一个$y$的系数,由$4$变成了$7$,那么很明显,另外一项不可能再增加,应该减小才能维持平衡,那么该降低多少呢?前面$x$上升了一个$y$的系数$7$,那么这里$y$应该下降一个$x$的系数,也就是$2$,那么最后的结果是$x=11,y=-3$,于是我们就得到了一个新解。
另外,特别的对于方程$ax+by=1$有整数解当且仅当整数$a$和$b$互素,因为上面讲了,只有两个数互素他们的最大公倍数才是$1$.

扩展欧几里得算法的代码实现

扩展欧几里得算法就是再求$a$,$b$的最大公约数$m=gcd(a,b)$的同时求出裴蜀等式的一个解
那么如何实现扩展欧几里得算法呢?首先我们可以来看先来看看最终状态,那么欧几里得算法的最终状态是什么呢?
很显然,和之前的欧几里得算法一样:
$$a=gcd,b=0时$$。
那么此时的裴蜀等式是什么样的呢:
$$ gcd*x+0*y=gcd$$
很显然,此时$x=1$,$y$可以为任意值。但是此时的$x=1$,仅为这种情况下的解,此时的$a$和$b$已经早就不是之前的了。但是我们可以利用递归的性质,层层倒推回去。所以我们此时需要找到递推式,也就是递推关系。
我们可以来看看相邻两种状态下的裴蜀等式:
(1)$$a*x+b*y=gcd$$
(2)$$b*x1+(a\%b)*y1=gcd$$
同时通过余数的性质我们可以得到
(3)$$a\%b=a-(a/b)*b$$

需要注意的是这里的除法是C++中的除法,即两个整数相除只保存整数部分,小数部分舍去

那么将式$3$代入式$2$可以得到
$$gcd=b *x1+(a-(a/b)*b)*y1$$
$$gcd=b*x1+a*y1-(a/b)*b*y1$$
(4)
$$gcd=a*y1+b*(x1-a/b*y1)$$
等式$4$等价于等式$3$,所以将等式$4$与等式$1$相比较可以得到以下递推式
$$x=y1 $$
$$y=x1-a/b*y1$$
其中$x1$和$y1$分别为下一个状态时求出的不定方程的解,得到递推式后利用递归的性质可以简单实现。
需要注意的是由于我们要求的是$a*x+b*y=m$的解,而最后我们计算出的却是$a*x+b*y=d$所以结果应该乘以一个$\frac{m}{d}$才是最终结果。

int x, y;//使用全局变量储存上一轮返回的x和y
int ext_gcd(int a, int b)
{
    if (b == 0)
    {
        x = 1;
        y= 0;//解除b=0时的一种解
        return a;
    }
    int re=ext_gcd(b, a % b);//计算下一轮
    int x1 = x;//暂时保存x的值用于计算y
    x = y;
    y = x1 - a / b * y;
    return re;//返回

}
int bdfc(long a, long b, long m)//计算不定方程解的函数
{
    int d = ext_gcd(a, b);//计算GCD
    if (m % d != 0)cout << "无解" << endl;//如果GCD不是倍数则无解
    long n = m / d;//解乘以M/D
    x *= n;
    y*= n;
    return d;
}

使用扩展欧几里得算法求解同余方程

同余方程长下面这个样子
$$ a x \equiv 1 \pmod {b}$$
脑壳痛,弱鸡我高数重修了,先来说一下这个三线等号是啥意思吧,这个三线等号(≡)叫做同余符号。这式子的意义为对于两个整数$ax$,$1$,他们除以整数$b$所得的余数相等,那么则称$ax$,$1$对于模$b$同余。
由于在这里
$$1 mod b=1$$
那么求解该同余方程也就是求
$$ax mod b=1$$
要解决这个问题我们需要了解一个定理
定理:

方程$ax+by=c$ 与方程$a x \equiv 1 \pmod {b} $ 是等价的,有整数解的充要条件为$gcd(a,b)|c$ (也就是$c$是$a$和$b$最小公倍数的倍数)

有了这个定理,那么问题就转化为了求解以下方程:
$$ax+by=1$$
转化为了我们所熟悉的线性不定方程的解,利用扩展欧几里得定理可以轻松解决,但是需要注意的是,题目要求的是最小整数解,但是我们解出来得x并不一定就是正数,也不一定为正,所以我们需要用到上面讲到的裴蜀定理由解出一个解推出最小整数解。那么如何由一个解推出另外一个解呢?答案是加上或者减去一个$b/gcd(a,b)$.
下面我们来看看洛谷的一道题:

P1082 同余方程
题目描述
求关于$x$的同余方程 $a x \equiv 1 \pmod {b}$ 的最小正整数解。
输入格式
一行,包含两个正整数 $a,ba,b$,用一个空格隔开。
输出格式
一个正整数 即最小正整数解。输入数据保证一定有解。
输入输出样例
输入
$3 10$
输出
$7$

代码实现:

#include <iostream>
using namespace std;
int x,y;
int ext_gcd(int a,int b)//扩展欧几里得定理
{
    if(b==0)
    {
        x=1;
        y=0;
        return a;
    }
    int res=ext_gcd(b,a%b);
    int x1=x;
    x=y;
    y=x1-a/b*y;
    return res;

}
void bdfc(int a,int b,int m)//求二元线性不等式
{
    int d=ext_gcd(a,b);//计算最大公倍数
    if(m%d!=0)cout<<-1<<endl;//根据裴蜀定义,有解当且仅当m是gcd(a,b)的倍数
    int e=m/d;//每个解都应该乘以m/d
    x*=e;
    y*=e;
}
int main()
{
    int a,b;
    cin>>a>>b;
    bdfc(a,b,1);
    while(x>0)x-=b;//获得最小整数解
    while(x<0)x+=b;//这里gcd(a,b)=1所以可以直接减去b,否则则需要减去b/gcd(a,b)
    cout<<x<<endl;
    return 0;
}

特殊的同余方程-逆元

逆元是一类特殊的同余方程,满足$ax \pmod 1$形式的同余方程被称之为一整数a对同余n之模逆元。
例如:
求整数3对同余11的逆元x写作
$${\displaystyle 3x\equiv 1{\pmod {11}}}$$
于普通同余方程相比,与逆元等价的线性不定方程为
$$ax+nb=1$$
可以看见,逆元的$c$固定为了$1$,那么结合裴蜀定理可以得到$a$和$n$,必然互质。
利用逆元的一些特殊的性质我们可以$解决部分扩欧无法解决的问题,以题说明:

HDU1576 A/B
要求$(A/B)\%9973$,但由于$A$很大,我们只给出$n$($n=A\%9973$)(我们给定的$A$必能被$B$整除,且$gcd(B,9973) = 1$)。
$Input$:
数据的第一行是一个$T$,表示有T组数据。
每组数据有两个数$n$(0 <= n < 9973$)和$$B$$($$1 <= B <= 10^9$)。
$Output$
对应每组数据输出$(A/B)\%9973$。
$Sample Input$
$2$
$1000 53$
$87 123456789$
$Sample Output$
$7922$
$6060$

题目告诉我们要求的是$(A/B)\%9973$的结果,但是题目并没告诉我们$A$的值,给我们的是$n$($n=A\%9973$)和$B$和$9973$的最大公约数,那么我们该如何求呢?
首先我们这里要讲以下模的运算
$$(a + b) \% m == ((a \% m) + (b \% m)) \% m$$
$$(a - b) \% m == ((a \% m) - (b \% m)) \% m$$
$$(a * b) \% m == ((a \% m) * (b \% m)) \% m$$
模的运算存在以上性质,但是需要注意的是以上性质并不包含除法
$$\color{red}{(a / b) \% m == ((a \% m) / (b \% m)) \% m} $$
这是错误的,那么这如何把$(A/B)\%9973$拆开,把$n=A\%9973$的性质利用上去呢,答案是使用逆元。
既然我们没有除法分配律,那么我们可以把除法改写成乘法。利用逆元可以做到。
逆元具有初中学过的倒数的性质,他们互为逆运算,实际上,逆元就是模的逆运算.
在这里,我们边可以将其改写为
$$(A*B^{-1})\%9973$$
$B^{-1}$为B相对于模9973的逆元,此后再用模的乘法分配律将其展开为
$$((A \% 9973) * ( B^{-1} )) \% 9973$$
由于$A \% 9973$已知,所以我们只需要求出B的逆元,最后模上$9973$即可解决咯嗷。
对于逆元,我们可以像解普通同余方程那样利用扩展欧几里得算法去解。

#include <iostream>
using namespace std;
int x, y;
int ext_gcd(int a, int b)
{
    if (b == 0)
    {
        x = 1;
        y = 0;
        return a;
    }
    int res = ext_gcd(b, a % b);
    int _x = x;
    x = y;
    y = _x - a / b * y;
    return res;
}
int bdfc(int a, int b, int c)
{
    int d = ext_gcd(a, b);
    if (c % d != 0)return -1;
    int e = c / d;
    x *= e;
    y *= e;
}
int main()
{
    int T;
    cin >> T;
    while (T--)
    {
        int n, b;
        cin >> n >> b;
        bdfc(b, 9973, 1);


        while (x < 0)
        {
            x += 9973;
            y -= b;
        }
        while (x > 0&&x-9973>0)
        {
            x -=9973;
            y += b;
        }
        //cout << x << " " << y<<" ";
         cout<<(x*n)%9973<<endl;
        //cout << b * x + 9973 * y;

    }
    return 0;
}

刷题时间

/***                                                                          
 *          .,:,,,                                        .::,,,::.          
 *        .::::,,;;,                                  .,;;:,,....:i:         
 *        :i,.::::,;i:.      ....,,:::::::::,....   .;i:,.  ......;i.        
 *        :;..:::;::::i;,,:::;:,,,,,,,,,,..,.,,:::iri:. .,:irsr:,.;i.        
 *        ;;..,::::;;;;ri,,,.                    ..,,:;s1s1ssrr;,.;r,        
 *        :;. ,::;ii;:,     . ...................     .;iirri;;;,,;i,        
 *        ,i. .;ri:.   ... ............................  .,,:;:,,,;i:        
 *        :s,.;r:... ....................................... .::;::s;        
 *        ,1r::. .............,,,.,,:,,........................,;iir;        
 *        ,s;...........     ..::.,;:,,.          ...............,;1s        
 *       :i,..,.              .,:,,::,.          .......... .......;1,       
 *      ir,....:rrssr;:,       ,,.,::.     .r5S9989398G95hr;. ....,.:s,      
 *     ;r,..,s9855513XHAG3i   .,,,,,,,.  ,S931,.,,.;s;s&BHHA8s.,..,..:r:     
 *    :r;..rGGh,  :SAG;;G@BS:.,,,,,,,,,.r83:      hHH1sXMBHHHM3..,,,,.ir.    
 *   ,si,.1GS,   sBMAAX&MBMB5,,,,,,:,,.:&8       3@HXHBMBHBBH#X,.,,,,,,rr    
 *   ;1:,,SH:   .A@&&B#&8H#BS,,,,,,,,,.,5XS,     3@MHABM&59M#As..,,,,:,is,   
 *  .rr,,,;9&1   hBHHBB&8AMGr,,,,,,,,,,,:h&&9s;   r9&BMHBHMB9:  . .,,,,;ri.  
 *  :1:....:5&XSi;r8BMBHHA9r:,......,,,,:ii19GG88899XHHH&GSr.      ...,:rs.  
 *  ;s.     .:sS8G8GG889hi.        ....,,:;:,.:irssrriii:,.        ...,,i1,  
 *  ;1,         ..,....,,isssi;,        .,,.                      ....,.i1,  
 *  ;h:               i9HHBMBBHAX9:         .                     ...,,,rs,  
 *  ,1i..            :A#MBBBBMHB##s                             ....,,,;si.  
 *  .r1,..        ,..;3BMBBBHBB#Bh.     ..                    ....,,,,,i1;   
 *   :h;..       .,..;,1XBMMMMBXs,.,, .. :: ,.               ....,,,,,,ss.   
 *    ih: ..    .;;;, ;;:s58A3i,..    ,. ,.:,,.             ...,,,,,:,s1,    
 *    .s1,....   .,;sh,  ,iSAXs;.    ,.  ,,.i85            ...,,,,,,:i1;     
 *     .rh: ...     rXG9XBBM#M#MHAX3hss13&&HHXr         .....,,,,,,,ih;      
 *      .s5: .....    i598X&&A&AAAAAA&XG851r:       ........,,,,:,,sh;       
 *      . ihr, ...  .         ..                    ........,,,,,;11:.       
 *         ,s1i. ...  ..,,,..,,,.,,.,,.,..       ........,,.,,.;s5i.         
 *          .:s1r,......................       ..............;shs,           
 *          . .:shr:.  ....                 ..............,ishs.             
 *              .,issr;,... ...........................,is1s;.               
 *                 .,is1si;:,....................,:;ir1sr;,                  
 *                    ..:isssssrrii;::::::;;iirsssssr;:..                    
 *                         .,::iiirsssssssssrri;;:.                      
 */

P1516 青蛙的约会

题目描述
两只青蛙在网上相识了,它们聊得很开心,于是觉得很有必要见一面。它们很高兴地发现它们住在同一条纬度线上,于是它们约定各自朝西跳,直到碰面为止。可是它们出发之前忘记了一件很重要的事情,既没有问清楚对方的特征,也没有约定见面的具体位置。不过青蛙们都是很乐观的,它们觉得只要一直朝着某个方向跳下去,总能碰到对方的。但是除非这两只青蛙在同一时间跳到同一点上,不然是永远都不可能碰面的。为了帮助这两只乐观的青蛙,你被要求写一个程序来判断这两只青蛙是否能够碰面,会在什么时候碰面。
我们把这两只青蛙分别叫做青蛙A和青蛙B,并且规定纬度线上东经0度处为原点,由东往西为正方向,单位长度1米,这样我们就得到了一条首尾相接的数轴。设青蛙A的出发点坐标是x,青蛙B的出发点坐标是y。青蛙A一次能跳m米,青蛙B一次能跳n米,两只青蛙跳一次所花费的时间相同。纬度线总长L米。现在要你求出它们跳了几次以后才会碰面。
输入格式:
输入只包括一行5个整数x,y,m,n,L其中0<x≠y < =2000000000,0 < m、n < =2000000000,0 < L < =2100000000。
输出格式:
输出碰面所需要的天数,如果永远不可能碰面则输出一行"Impossible"。

解题思路:


首先我们把图画出来,两只青蛙初始位置为a和b,他们一起向西跳(正方向),a每次能跳m,b每次能跳n,设跳数为x,那么我们可以写出这两只青蛙关于离远点距离的方程。
$$y1=(a+mx)\%L$$
$$y2=(b+nx)\%L$$
这里之所以需要让其和L取模是因为这是一个环,到达0点就重新计算了。看图可以很轻松的理解。
那么什么时候相遇呢?很显然y1=y2时,那么可以得到:
$$(a+mx)\%L=(b+nx)\%L$$
是不是有点眼熟,是的,其可以写做如下形式:
$$a +mx \equiv b+nx \pmod {L} $$
我们只需要解出这个方程,x即为所求。
那么该如何解呢?很明显,这和同余方程的一般形式不同,所以不能直接用定理去转化为二元不定方程来解,但是我们可以稍微的代换一下其中的变量
(1)$$(a+mx)=Lt1+余数$$
(3)$$(b+nx)=Lt2+余数$$
使1式和2式相减可以得到
$$(m-n)x+(a-b)=L(t1-t2)$$
我们可以将t1-t2看作一个未知数,设y=t1-t2,那么原式为
$$(m-n)x+Ly=(a-b)$$
m,n,a,b都是已知的,所以就转换为了一个线性不定方程,使用扩展欧几里得定理可以解决。

#include <iostream>
using  namespace std;
long _x,_y;
long ext_gcd(long a,long b)
{
    if(b==0)
    {
        _x=1;
        _y=0;
        return a;
    }
    long ans=ext_gcd(b,a%b);
    long tmp=_x;
    _x=_y;
    _y=tmp-a/b*_y;
    return ans;
}
void bdfc(long a,long b,long c)
{
    long gcd=ext_gcd(a,b);
    if(c%gcd!=0) {
        cout<<"Impossible"<<endl;
        return;
    }
    long e=c/gcd;
    _x*=e;
    _y*=e;
    while(_x<0)_x+=b;
   while(_x>0&&_x-b>0)_x-=b;
    //cout<<_x<<" "<<_y<<endl;
    cout<<_x<<" "<<endl;
}
int main()
{
    long x,y,m,n,L;
    cin>>x>>y>>m>>n>>L;
    long a=m-n;
    long b=L;
    long c=y-x;
    //cout<<a<<" "<<b<<" "<<c;
    bdfc(a,b,c);

    return 0;
}

emmm,很迷WA了两个点

不知道错在了哪里,如果有大佬可以帮忙看看。看洛谷评论说的是题目错了,脑壳痛,先放放吧,等几天再来解决

最后笔者想吐槽一下这两只青蛙。。。这两只青蛙是智障吧。。。两个蛙都向原点走不好吗,都要朝着同一方向

参考资料

1.【日】秋叶拓哉,【日】岩田阳一 【日】北川宜稔 《挑战程序设计竞赛 第二版》人民邮电出版社 ISBN 978-7-115-32010-0
2.IOI国家集训队论文-2009-《欧几里得算法的应用》-金斌
3.模逆元-维基百科,自由的百科全书
4.[模运算中的倒数 - 乘法逆元](http://conw.net/archives/6/

版权信息

知识共享许可协议
本作品采用知识共享署名-相同方式共享 4.0 国际许可协议进行许可。

Last modification:September 28th, 2019 at 09:15 pm
如果觉得我的文章对你有用,请随意赞赏

Leave a Comment