zoukankan      html  css  js  c++  java
  • 拓展欧几里德(例题:同余方程)

    例题 同余方程

    题目描述

    求关于 (x) 的同余方程 (a x equiv 1 pmod {b}) 的最小正整数解。
    输入输出格式
    输入格式:

    一行,包含两个正整数 (a,b),用一个空格隔开。

    输出格式:

    一个正整数 (x_0)​,即最小正整数解。输入数据保证一定有解。

    输入输出样例
    输入样例#1:

    3 10
    输出样例#1:

    7
    说明

    【数据范围】
    对于 100%的数据,(2 ≤a, b≤ 2,000,000,000)

    NOIP 2012 提高组 第二天 第一题

    问题转化

    • 题目问的是满足 (ax mod b = 1) 的最小正整数 (x) 。( (a,b) 是正整数),但是不能暴力枚举 (x) ,会超时。

    • 把问题转化一下。观察 (ax mod b = 1) ,它的实质是 (ax + by = 1) :这里 (y) 是我们新引入的某个整数,并且似乎是个负数才对。这样表示是为了用扩展欧几里得算法。我们将要努力求出一组 (x,y) 来满足这个等式。稍微再等一下——

    • 问题还需要转化。扩展欧几里得是用来求 (ax + by = gcd(a,b)) 中的未知数的,怎么牵扯到等于 (1) 呢?

    • 原理是,方程 (ax + by = m) 有解的必要条件是 (m mod gcd(a,b) = 0)

    • 这个简单证一下。

    由最大公因数的定义,可知 (a)(gcd(a,b)) 的倍数,且 (b)(gcd(a,b)) 的倍数,

    (x,y) 都是整数,就确定了 (ax + by)(gcd(a,b)) 的倍数,

    因为 (m = ax + by) ,所以 (m) 必须是 (gcd(a,b)) 的倍数,

    那么 (m mod gcd(a,b) = 0)

    • 可得出在这道题中,方程 (ax + by = 1) 的有解的必要条件是 (1 mod gcd(a,b) = 0) 。可怜的 (gcd(a,b)) 只能等于 (1) 了。这实际上就是 (a,b) 互质。

    扩展欧几里得

    • 前提:知道普通欧几里得算法(辗转相除法)。

    (下面字母挺多,希望你耐心地慢慢地读~)

    • 我们拿到了一组 (a,b) 。设 (G = gcd(a, b)) 。那么目标是求出满足 (ax + by = G) (①) 的整数 (x)(y) 。其中 (x) 应当是满足条件的最小正整数,即答案,而 (y) 是辅助答案。

    • (注意,虽然刚刚已经证明本题的 (gcd(a,b)) 必然等于 (1) ,但是扩展欧几里得算法本身过程求的是 ax + by = gcd(a,b)$ 的解。它正好非常适合本题,不过我们要按照它求解的过程去做)

    • 如果我们先前已经求出了另一组数 (x_2, y_2) ,它们满足这么一个式子: (bx_2 + (a mod b)y_2 = G) (②),则此时结合①②一定有:

    (ax + by = bx_2 + (a mod b)y_2)​ (③)

    • 可见,在这个“如果”实现的时候,我们的目标变成了“求出满足上式的 (x)(y) ”。

    • 其中 (a,b,x_2,y_2)​ 都已知, (x,y) 待求。因为未知数比方程更多,所以没有唯一解。我们先求出一组必然存在的解,最后将在“答案处理”时转为最小解。

    • 怎么求呢?取模运算是 (a mod b = a - b×(a/b)) ,所以方程③实际上是:

    (ax + by = bx_2 + (a-b×(a/b))y_2)

    (Rightarrow ax + by = bx_2 + ay_2 - b × (a/b)y_2)

    (Rightarrow ax + by = ay_2 + b(x_2-(a/b)y_2))

    • 看上面这个方程,一组必然存在的解出现了:

    (x = y_2, y = x_2 - (a/b)y_2) (④)

    • 可见,我们只要求出 (x_2,y_2) ,就能得出正确的 (x,y) 。问题是 (x_2,y_2) 怎么求。

    现在我们手上是 (b,a mod b) 这两个系数,而目标是求出 (x_2)​ 和 (y_2) 满足:

    (bx_2 + (a mod b)y_2 = G)(②)

    把①和②对比一下:

    (ax + by = G) (①)

    原方程中的系数 (a) 变成了②中的系数 (b) ,原方程中的系数 (b) 变成了②中的 (a mod b)而已。

    • 所以,把新的方程也看作 (ax + by = G) 的形式(只是系数 (a)(b) 的具体数值改变了)。然后按照上面的一模一样下来(其实都只是推导过程),我们发现,最好有 (x_3,y_3) 来支撑 (x_2, y_2)

    • 再一模一样下来,我们又需要 (x_4,y_4) 来支撑 (x_3, y_3)

    ……

    这个递归中 (a, b) 不断被替代为 (b, a mod b),这个替换方式与普通欧几里得是一样的,所以最后会出现 (a_n = G, b_n = 0)

    这时要直接返回了,我们需要一组 (x_n,y_n) 满足 (a_nx_n + b_ny_n = G)(⑤)。

    然而该层的 (a_n = G, b_n = 0) 。所以只要⑤左边取 (x_n = 1) ,这个方程就妥妥的成立了。

    (最后一层的 (y_n) 建议取 0 。然而由于 (b = 0) ,就算返回其它数值,方程也一定成立。但这样的程序容易出错,因为 (y_n)​ 在回溯时滚雪球式增长,容易数值越界。

    最后一层结束后,就开始返回,直到最上层。每一层可以轻松地根据下层的 (x{k+1},y_{k+1}) 求出当前层的 (x_k, y_k)

    整个过程就是:以辗转相除的方式向下递进,不断缩小系数,保证会出现有确定解的最后一层。
    答案处理

    一个遗留问题:我们将要求出来的 (x,y) 仅仅保证满足 (ax+by=1) ,而 (x) 不一定是最小正整数解。有可能太大,也有可能是负数。这依然可以解决,那就是, (x) 批量地减去或加上 (b) ,能保证 (ax+by=1),因为:

    (ax+by=1)

    (ax+by+k×ba−k×ba=1)

    (a(x+kb)+(y−ka)b=1)

    1.显然这并不会把方程中任何整数变成非整数。

    2.“加上或减去 (b) ”不会使 (x) 错过任何解。可以这么理解:

    已经求出一组整数 (x,y) 使得 (ax+by=1),也就是 ((1−ax)/b=y)(y) 是整数,可见目前 (1−ax)(b) 的倍数。

    现在想改变 (x) 并使得方程仍然成立。已知 (a,b) 互质,假若 (x) 的变化量( (Δx))不是 (b) 的倍数,则 (1−ax) 的变化量( (−a×Δx) )也不是 (b) 的倍数,这会使得 (1−ax) 不再是 (b) 的倍数,则 (y) 不是整数了。

    仅当 (x) 的变化量是 (b) 的倍数时, (1−ax) 能保持自己是 (b) 的倍数,此时就出现新的解了。

    因此到最后,如果 (x) 太小就不断加 (b) 直到大于等于 (0)$ ,太大则一直减 (b) ,直到最小正整数解。也就是这么写:

    x = (x % b + b) % b;//括号中取模再加,可以处理负数
    

    代码

    推导中的所有 (x,y) 共用全局变量 (long long) 传递也很方便。

    #include <bits/stdc++.h>
    #define int long long //全局变量 long long
    using namespace std;
    int x,y;////目前方程真正的解 
    void exgcd(int a,int b){////当前目的:求解 ax + by = gcd(a, b) 这么一个方程
    	if(b==0){ //在 b = 0 时方程还要成立, 使 x = 1, y = 0 ,必然成立 
    		x=1;y=0;//建议返回0。
    		return;
    	}
    	exgcd(b,a%b);//先求下一个方程的解 
    	//现在我们已经拿到了下一个方程的解x, y
    	int res=x;//暂时存一下x
    	x=y;
    	y=res-(a/b)*y;
    }
    signed main(){
    	int a,b;
    	scanf("%lld%lld",&a,&b);
    	exgcd(a,b);
    	printf("%lld
    ",(x%b+b)%b);//我们求出来的x必然满足方程,但不一定是最小正整数解,所以要进行答案处理
    	return 0;
    }
    

    (当然可能《进阶指南》上或许有更清楚的解说)

  • 相关阅读:
    Ubuntu配置sublime text 3的c编译环境
    ORA-01078错误举例:SID的大写和小写错误
    linux下多进程的文件拷贝与进程相关的一些基础知识
    ASM(四) 利用Method 组件动态注入方法逻辑
    基于Redis的三种分布式爬虫策略
    Go语言并发编程总结
    POJ2406 Power Strings 【KMP】
    nyoj 会场安排问题
    Server Tomcat v7.0 Server at localhost was unable to start within 45 seconds. If the server requires more time, try increasing the timeout in the server editor.
    Java的String、StringBuffer和StringBuilder的区别
  • 原文地址:https://www.cnblogs.com/Lour688/p/13228948.html
Copyright © 2011-2022 走看看