• 模乘逆元与孙子定理


    孙子定理也称为中国剩余定理。

    《孙子算经》卷下第二十六题(“物不知数”问题):有物不知其数,三三数之剩二,五五数之剩三,七七数之剩二。问物几何?

    孙子定理讲的是求解一元线性同余方程组的方法。


    有人指出,孙子定理有以下5种解法:

    1.枚举法
    2.解不定方程法
    3.逐级满足法
    4.化为相同除数的同余式法
    5.典经同余式方程组解法

    据说最为简洁快速的是第4种方法。

    这里介绍的程序为第5种解法。

    要得到解x,首先需要计算模乘逆元。模乘逆元定义为:

    满足 ab≡1(mod m),称b为a模乘逆元。

    程序中给出了两种求模乘逆元的方法,一种是穷举法,此法对于模值较大时是一种低效率的方法;另外一种是利用扩展欧几里得算法求逆元。

    求得模乘逆元之后,就可以算出同余方程组的解。

    “物不知数”问题的一个解为23(最小值解),其解系为x=23+105k(k>=0)(105=3*5*7)。

    #include <stdio.h>
    
    // 递推法实现扩展欧几里德算法
    long exgcd(long a, long b, long *x, long *y)
    {
        long x0=1, y0=0, x1=0, y1=1;
        long r, q;
        *x=0;
        *y=1;
    
        r = a % b;
        q = (a - r) / b;
        while(r)
        {
            *x = x0 - q * x1;
            *y = y0 - q * y1;
            x0 = x1;
            y0 = y1;
            x1 = *x;
            y1 = *y;
    
            a = b;
            b = r;
            r = a % b;
            q = (a - r) / b;
        }
        return b;
    }
    
    // 扩展欧几里德算法求逆元
    long minv(long a, long p)
    {
        long x, y;
    
        exgcd(a, p, &x, &y);
    
        return x<0 ? x+p : x;
    }
    
    // 试探法求逆元
    long minv2(long a, long p)
    {
        long y=1, t;
        int i;
    
        if(a < 0) {
            a = a % p;
            a += p;
        }
    
        for(i=1; i<p; i++) {
            t = a * i;
            if(t % p == 1) {
                y = i;
                return y;
            }
        }
    
        return y;
    }
    
    int main(void)
    {
        printf("a=%d m=%d x=%ld x=%ld
    ", 65, 83, minv(65, 83), minv2(65, 83));
        printf("a=%d m=%d x=%ld x=%ld
    ", 11663, 103, minv(11663, 103), minv2(11663, 103));
        printf("a=%d m=%d x=%ld x=%ld
    ", 11227, 107, minv(11227, 107), minv2(11227, 107));
        printf("a=%d m=%d x=%ld x=%ld
    ", 11021, 103, minv(11021, 109), minv2(11021, 109));
    
        // 孙子算经
        long a[] = {2, 3, 2};
        long m[] = {3, 5, 7};
    
        int i, size= sizeof(a)/sizeof(long);
        long bm=1, subm[size], x=0;
    
        for(i=0; i<size; i++)
            bm *= m[i];
        for(i=0; i<size; i++)
            subm[i] = bm / m[i];
        for(i=0; i<size; i++) {
            x += subm[i] * minv(subm[i], m[i]) * a[i];
            x %= bm;
        }
    
        printf("x=%ld
    ", x);
    
        return 0;
    }

    程序运行结果:

    a=65 m=83 x=23 x=23
    a=11663 m=103 x=73 x=73
    a=11227 m=107 x=40 x=40
    a=11021 m=103 x=100 x=100
    x=23


  • 相关阅读:
    Winfrom 动画实现
    Android-SD卡相关操作
    Android-动态权限获取
    Java 常用知识点
    无锁队列的实现
    稳定的快排
    设计模式
    map的线程安全
    win 消息
    memecpy源码
  • 原文地址:https://www.cnblogs.com/tigerisland/p/7564911.html
Copyright © 2020-2023  润新知