辗转相除法

辗转相除是求两个数的最大公约数的。

要证这个定理成立，只需要证明

gcd(a, b) = gcd(b, a % b) 就行了

证明：令a % b = r, 所以a = k * b + r, 所以r = a - k * b,假设d为a，b的一个公约数，那么 d|a， d|b,（d|a的意思就是d整除a，也就是a能被d整除），所以a - k * b 也一定能被d整除，即 d|r， 也就是 d|(a % b)， 因此d也是b 和 (a % b)的公约数，因此a,b 的公约数和b, (a%b)的公约数也是一样的，其最大公约数也一定相同，所以gcd(a, b) = gcd(b, a % b);

所以有了这个等式之后，基本上就算完了，还有一步就是怎么到最后求个具体的数，当b等于0时候就可以了，因为最后递归好多还是和原来的那个公约数是相同的，最后有0了，他俩的最大公约数就是他本身了，也就是a了，

代码:

1
2
3
4

long long gcd(long long a,long long b)  //辗转相除求最大公约数
{
    return b?gcd(b,a%b):a;
}

扩展欧几里得算法

算法思想

现在我们知道了 a 和 b 的最大公约数是 gcd ，那么，我们一定能够找到这样的 x 和 y ，使得: a*x + b*y = gcd 这是一个不定方程，有多解是一定的，但是只要我们找到一组特殊的解 x0 和 y0 那么，我们就可以用 x0 和 y0 表示出整个不定方程的通解：

1
2

x = x0 + (b/gcd)*t
y = y0 – (a/gcd)*t

可以通过将通解代入原方程来验证正确性.

为什么不是：

x = x0 + b*t
y = y0 – a*t

那是因为：

b/gcd 是 b 的因子， a/gcd 是 a 的因子是吧？那么，由于 t的取值范围是整数，你说 (b/gcd)*t 取到的值多还是 b*t 取到的值多？同理，(a/gcd)*t 取到的值多还是 a*gcd 取到的值多？那肯定又要问了，那为什么不是更小的数，非得是 b/gcd 和a/gcd ？

注意到：我们令 B = b/gcd ， A = a/gcd ， 那么，A 和 B 一定是互素,若要让两个通解代入原方程后等式成立,A和B必须缩小相同的倍数,即有公约数.所以找不到比A和B更小的数了

现在，我们知道了一定存在 x 和 y 使得 ： a*x + b*y = gcd ， 那么，怎么求出这个特解 x 和 y 呢？只需要在欧几里德算法的基础上加点改动就行了。

我们观察到：欧几里德算法停止的状态是： a= gcd ， b = 0 ，那么，这是否能给我们求解 x y 提供一种思路呢？因为，这时候，只要 a = gcd 的系数是 1 ，那么只要 b 的系数是 0 或者其他值（无所谓是多少，反正任何数乘以 0 都等于 0 但是a 的系数一定要是 1），这时，我们就会有： a*1 + b*0 = gcd

当然这是最终状态，但是我们是否可以从最终状态反推到最初的状态呢？

假设当前我们要处理的是求出 a 和 b的最大公约数，并求出 x 和 y 使得 a*x + b*y= gcd ，而我们已经求出了下一个状态：b 和 a%b 的最大公约数，并且求出了一组x1 和y1 使得： b*x1 + (a%b)*y1 = gcd ， 那么这两个相邻的状态之间是否存在一种关系呢？

我们知道： a%b = a - (a/b)*b（这里的 “/” 指的是整除，例如 5/2=2 , 1/3=0），那么，我们可以进一步得到：

gcd = b*x1 + (a-(a/b)*b)*y1
    = b*x1 + a*y1 – (a/b)*b*y1
    = a*y1 + b*(x1 – a/b*y1)

对比之前我们的状态：求一组 x 和 y 使得：a*x + b*y = gcd ，是否发现了什么？

这里： x = y1 y = x1 – a/b*y1

以上就是扩展欧几里德算法的全部过程

代码:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20

void exgcd( long long int a, long long int b, long long int &x, long long int &y )
{
    if( b== 0 )
    {
        x= 1;
        y= 0;
        return;
    }
    exgcd( b, a% b, x, y );
    long long int t= x;
    x= y;
    y= t- a/ b* y;
    return;
}
//1、先计算Gcd(a, b)，若n不能被Gcd(a, b)整除，则方程无整数解；否则，在方程两边同时除以Gcd(a, b)，得到新的不定方程a' * x + b' * y = n'，此时Gcd(a', b')=1;
//2、求出方程a' * x + b' * y = 1的一组整数解x0, y0，则n' * x0,n' * y0是方程a' * x + b' * y = n'的一组整数解；
//3、根据数论中的相关定理，可得方程a' * x + b' * y = n'的所有整数解为：
//x = n' * x0 + b' * k
//y = n' * y0 - a' * k
//(t为整数)

举例

题目描述：

两只青蛙在网上相识了，它们聊得很开心，于是觉得很有必要见一面。它们很高兴地发现它们住在同一条纬度线上，于是它们约定各自朝西跳，直到碰面为止。可是它们出发之前忘记了一件很重要的事情，既没有问清楚对方的特征，也没有约定见面的具体位置。不过青蛙们都是很乐观的，它们觉得只要一直朝着某个方向跳下去，总能碰到对方的。但是除非这两只青蛙在同一时间跳到同一点上，不然是永远都不可能碰面的。为了帮助这两只乐观的青蛙，你被要求写一个程序来判断这两只青蛙是否能够碰面，会在什么时候碰面。

我们把这两只青蛙分别叫做青蛙A和青蛙B，并且规定纬度线上东经0度处为原点，由东往西为正方向，单位长度1米，这样我们就得到了一条首尾相接的数轴。设青蛙A的出发点坐标是x，青蛙B的出发点坐标是y。青蛙A一次能跳m米，青蛙B一次能跳n米，两只青蛙跳一次所花费的时间相同。纬度线总长L米。现在要你求出它们跳了几次以后才会碰面。

输入只包括一行5个整数x，y，m，n，L，其中x≠y < 2000000000，0 < m、n < 2000000000，0 < L < 2100000000。

输出碰面所需要的跳跃次数，如果永远不可能碰面则输出一行"Impossible"

Sample Input:

1 2 3 4 5

Sample Output:

4

题目分析：

构造方程 (x + m * s) - (y + n * s) = k * l(k = 0, 1, 2,…) 变形为 (n-m) * s + k * l = x - y。即转化为模板题，a * x + b * y = c，是否存在整数解。

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51

#include <stdio.h>
long long gcd(long long a,long long b)  //辗转相除求最大公约数
{
    return b?gcd(b,a%b):a;
}
void exgcd( long long int a, long long int b, long long int &x, long long int &y )
{
    if( b== 0 )
    {
        x= 1;
        y= 0;
        return;
    }
    exgcd( b, a% b, x, y );
    long long int t= x;
    x= y;
    y= t- a/ b* y;
    return;
}
//1、先计算Gcd(a, b)，若n不能被Gcd(a, b)整除，则方程无整数解；否则，在方程两边同时除以Gcd(a, b)，得到新的不定方程a' * x + b' * y = n'，此时Gcd(a', b')=1;
//2、求出方程a' * x + b' * y = 1的一组整数解x0, y0，则n' * x0,n' * y0是方程a' * x + b' * y = n'的一组整数解；
//3、根据数论中的相关定理，可得方程a' * x + b' * y = n'的所有整数解为：
//x = n' * x0 + b' * k
//y = n' * y0 - a' * k
//(t为整数)
int main(  )
{
    long long int x, y, m, n, l;
    while( scanf( "%lld %lld %lld %lld %lld", &x, &y, &m, &n, &l )!= EOF )
    {
        long long int a= n- m, b= l, c= x- y, p, q;
        long long int d= gcd( a, b );
        if( c% d )//a*x + b*y = c 有解的充要条件：c%gcd(a,b)==0
        {
            puts( "Impossible" );
            continue;
        }
        a/= d, b/= d, c/= d;
        exgcd( a, b, p, q );  //
//此时方程的所有解为：x = c*p - b*k,  x与cp关于b同余那么根据最小整数原理，一定存在一个最小的正整数，
//它是 a 关于m 的逆元，而最小的肯定是在（0 , m）之间的，而且只有一个。 当 m 是负数的时候，我们取 m 的绝对值就行了，
//当 x0 是负数的时候，他模上 m 的结果仍然是负数，我们仍然让对b取模，然后结果再加上b就行了
        p*= c;//将整个方程扩大c倍，让右边为c
        long long int ans= p% b;//求最小解
        while( ans< 0 )
        {
            ans+= b;
        }
        printf( "%lld\n", ans );
    }
}

乘法逆元

算法思想

扩展欧几里德可以求解形如 a*x +b*y = c 的通解，但是一般没有谁会无聊到让你写出一串通解出来，都是让你在通解中选出一些特殊的解，比如一个数对于另一个数的乘法逆元

什么叫乘法逆元？

这里，我们称 x 是 a 关于 m 的乘法逆元

这怎么求？可以等价于这样的表达式： a*x + m*y = 1

看出什么来了吗？没错，当gcd(a , m) != 1 的时候是没有解的这也是 a*x + b*y = c 有解的充要条件： c % gcd(a , b) == 0

接着乘法逆元讲，一般，我们能够找到无数组解满足条件，但是一般是让你求解出最小的那组解，怎么做？我们求解出来了一个特殊的解 x0 那么，我们用 x0 % m其实就得到了最小的解了。为什么？

可以这样思考：

x 的通解不是 x0 + m*t 吗？

那么，也就是说， a 关于 m 的逆元是一个关于 m 同余的，那么根据最小整数原理，一定存在一个最小的正整数，它是 a 关于m 的逆元，而最小的肯定是在（0 , m）之间的，而且只有一个，这就好解释了。

可能有人注意到了，这里，我写通解的时候并不是 x0 + (m/gcd)*t ，但是想想一下就明白了，gcd = 1，所以写了跟没写是一样的，但是，由于问题的特殊性，有时候我们得到的特解 x0 是一个负数，还有的时候我们的 m 也是一个负数这怎么办？

当 m 是负数的时候，我们取 m 的绝对值就行了，当 x0 是负数的时候，他模上 m 的结果仍然是负数（在计算机计算的结果上是这样的，虽然定义的时候不是这样的），这时候，我们仍然让 x0 对abs(m) 取模，然后结果再加上abs(m) 就行了，于是，我们不难写出下面的代码求解一个数 a 对于另一个数 m 的乘法逆元：

 1
 2
 3
 4
 5
 6
 7
 8
 9
10

void extend_gcd (ll a , ll b , ll& d, ll &x , ll &y) {
	if(!b){d = a; x = 1; y = 0;}
	else {extend_gcd(b, a%b, d, y, x); y-=x*(a/b);}
}
ll inv(ll a, ll m) { //计算%m下 a的逆。如果不存在逆return -1;
	ll d, x, y;
	m=abs(m);
	extend_gcd(a, m, d, x, y);
	return d == 1 ? (x+m)%m : -1;
}

中国剩余定理

算法思想

问题：今有物不知其数，三三数之剩二，五五数之剩三，七七数之剩二。问物几何？

解法：三三数之剩二，置一百四十；五五数之剩三，置六十三；七七数之剩二，置二十三，并之，得二百三十三。以二百一十减之，即得。凡三三数之剩一，则置七十；五五数之剩一，则置二十一；七七数之剩一，则置十五。一百六以上，一百五减之，即得。

说明白一点就是说，存在一个数x，除以3余2，除以5余三，除以7余二，然后求这个数。上面给出了解法。再明白这个解法的原理之前，需要先知道一下两个定理。

定理1：几个数相加，如果只有一个加数，不能被数a整除，而其他加数均能被数a整除，那么它们的和，就不能被整数a整除。

定理2：二数不能整除，若被除数扩大（或缩小）了几倍，而除数不变，则其余数也同时扩大（或缩小）相同的倍数（余数必小于除数）。

这两个定理浅显易懂就不再赘述了。

先给出求解该问题的具体步骤：

1. 求出最小公倍数

    lcm=3\*5\*7=105
   

2. 求各个数所对应的基础数

   1. 105÷3=35

       35÷3=11......2 //基础数35(满足除3余2)
      

   2. 105÷5=21

       21÷5=4......1
      

      定理2把1扩大3倍得到3，那么被除数也扩大3倍，得到

       21\*3=63//基础数63(满足除5余3)
      

   3. 105÷7=15

       15÷7=2......1
      

      定理2把1扩大2倍得到2，那么被除数也扩大2倍，得到

       15*2=30//基础数30(满足除7余2)
      

3. 把得到的基础数加和（注意：基础数不一定就是正数） 35+63+30=128

4. 减去最小公倍数lcm（在比最小公倍数大的情况下） x=128-105=23

那么满足题意得最小的数就是23了。一共有四个步骤。下面详细解释每一步的原因。

（1）最小公倍数就不用解释了，跳过（记住，这里讨论的都是两两互质的情况）

（2）观察求每个数对应的基础数时候的步骤，比如第一个。105÷3=35。显然这个35是除了当前这个数不能整除以外都能够被其他数整除，就是其他数的最小公倍数。相当于找到了最小的开始值，用它去除以3发现正好余2。那么这个基础数就是35。记住35的特征，可以整除其他数但是不能被3整除，并且余数是2。

再看下一个5对应基础数。21是其他数的最小公倍数，但是不能被5整除，用21除以5得到的余数是1，而要求的数除以5应该是余1的。所以余数被扩大，就得到了相应的基础数63。记住这个数的特征，可以被其他数整除但是被5除应该余三。同理，我们得到了第三个基础数23，那么他的特征就是：可以被其他数整除，但是不能被7整除，并且余数为2。

（3）第三步基础数加和，为什么要这样做呢？利用就是上面提到的定理1。

35+63+30=128。对于3来说，可以把63+30的和看作一个整体，应该他们都可以被3整除。看着上面写出的三个数的特征，运用定理1来说，就是在35的基础上加上一个可以被3整除的倍数，那么得到的结果依然还是满足原先的性质的，就是128除以同样还是余2的。同理，对于5还说，这个数被除之后会剩余3；对于7来说，被除之后剩余2。所以说，我们当前得到的这个数是满足题目要求的一个数。但是这个数是不是最小的，那就不一定了。

（4）应该不能确定是不是最小的数，这个时候就要用到他们的最小公倍数了。最小公倍数顾名思义，一定是一个同时被几个数整除的最小的一个数，所以减去它剩余下来的余数还是符合题意要求的。当然也同样可以运用定理1来解释，只不过是加法变成了减法，道理还是一样的。当然具体要不要减还是要看和lcm的大小关系的。

编码

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16

int CRT(int a[],int m[],int n)
{
    int M = 1;
    int ans = 0;
    for(int i=1; i<=n; i++)//M=m1*m2*m3…………*mn
        M *= m[i];
    for(int i=1; i<=n; i++)
    {
        int x, y;
        int Mi = M / m[i];
        extend_Euclid(Mi, m[i], x, y);//扩展欧几里得求逆元x，因为mi互素，此处直接普通扩展欧几里得就可求得逆元
        ans = (ans + Mi * x * a[i]) % M;
    }
    if(ans < 0) ans += M;
    return ans;
}

我们能够针对每一个mi都能够求出来一个Mi，但是注意基础数可以是负数呀。同时我们看到函数的语句ans = (ans + Mi * x * a[i]) % M;负数取余难道不会受影响吗？答案是不会的。对于负数取余，不同的语言有不同的标准，举个简单例子：-7%3。在C++语言里面，答案是-1，Java语言里面也是-1，但是在Python里面就是2了。但不管是-1还是2都是对了.

(n+ret%n)%n;这是return语句。加了n之后再取余，得到的就是最小的正余数了，真是机智。同样你要是觉得不够保险，可以把上面的语句改成ans = (ans + Mi * x * a[i]+M) % M;这样足够保险了。同样，在其他题目里面凡是遇到取余的操作都需要去这样做，就可以消除因为编译器不同带来的烦恼了。

举例

问题描述

人自出生起就有体力，情感和智力三个生理周期，分别为23，28和33天。一个周期内有一天为峰值，在这一天，人在对应的方面（体力，情感或智力）表现最好。通常这三个周期的峰值不会是同一天。现在给出三个日期，分别对应于体力，情感，智力出现峰值的日期。然后再给出一个起始日期，要求从这一天开始，算出最少再过多少天后三个峰值同时出现。

问题分析

首先我们要知道，任意两个峰值之间一定相距整数倍的周期。假设一年的第N天达到峰值，则下次达到峰值的时间为N+Tk(T是周期，k是任意正整数)。所以，三个峰值同时出现的那一天(S)应满足

  S = N1 + T1*k1 = N2 + T2*k2 = N3 + T3*k3

N1,N2,N3分别为为体力，情感，智力出现峰值的日期， T1，T2,T3分别为体力，情感，智力周期。 我们需要求出k1,k2,k3三个非负整数使上面的等式成立。

思路:

想直接求出k1,k2,k3貌似很难，但是我们的目的是求出S， 可以考虑从结果逆推。根据上面的等式，S满足三个要求：除以T1余数为N1，除以T2余数为N2，除以T3余数为N3。这样我们就把问题转化为求一个最小数，该数除以T1余N1，除以T2余N2,除以T3余N3。这就是著名的中国剩余定理

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58

#include <iostream>
#include <string.h>
#include <stdio.h>

using namespace std;

int a[4], m[4];

void extend_Euclid(int a, int b, int &x, int &y)//扩展欧几里得
{
    if(b == 0)
    {
        x = 1;
        y = 0;
        return;
    }
    extend_Euclid(b, a % b, x, y);
    int tmp = x;
    x = y;
    y = tmp - (a / b) * y;
}

int CRT(int a[],int m[],int n)
{
    int M = 1;
    int ans = 0;
    for(int i=1; i<=n; i++)//M=m1*m2*m3…………*mn
        M *= m[i];
    for(int i=1; i<=n; i++)
    {
        int x, y;
        int Mi = M / m[i];
        extend_Euclid(Mi, m[i], x, y);//扩展欧几里得求逆元x，因为mi互素，此处直接普通扩展欧几里得就可求得逆元
        ans = (ans + Mi * x * a[i]) % M;
    }
    if(ans < 0) ans += M;
    return ans;
}
int main()
{
    int p, e, i, d, t = 1;
    while(cin>>p>>e>>i>>d)
    {
        if(p == -1 && e == -1 && i == -1 && d == -1)
            break;
        a[1] = p;
        a[2] = e;
        a[3] = i;
        m[1] = 23;
        m[2] = 28;
        m[3] = 33;
        int ans = CRT(a, m, 3);
        if(ans <= d)
            ans += 21252;
        cout<<"Case "<<t++<<": the next triple peak occurs in "<<ans - d<<" days."<<endl;
    }
    return 0;
}