详解--乘法逆元

怼烎@ 2022-06-11 07:48 364阅读 0赞

今天我们来探讨逆元在**ACM-ICPC**竞赛中的应用，逆元是一个很重要的概念，必须学会使用它。

对于正整数![20140613102654328][]和![20140613102712781][]，如果有![20140613102734984][]，那么把这个同余方程中![20140613102856531][]的最小正整数解叫做![20140613102654328][]模![20140613102712781][]的逆元。

逆元一般用扩展欧几里得算法来求得，如果![20140613102712781][]为素数，那么还可以根据费马小定理得到逆元为![20140613103413828][]。

推导过程如下

![20140613104248984][]

求现在来看一个逆元最常见问题，求如下表达式的值**（已知![20140613104752312][]）**

![20140613104619203][]

当然这个经典的问题有很多方法，最常见的就是扩展欧几里得，如果![20140613102712781][]是素数，还可以用费马小定理。

但是你会发现费马小定理和扩展欧几里得算法求逆元是有局限性的，它们都会要求![20140613102654328][]与![20140613102712781][]互素。实际上我们还有一

种通用的求逆元方法，适合所有情况。公式如下

![20140613105646406][]

现在我们来证明它，已知![20140613104752312][]，证明步骤如下

![20140613112023328][]

接下来来实战一下，看几个关于逆元的题目。

**题目：[http://poj.org/problem?id=1845][http_poj.org_problem_id_1845]**

**题意：**给定两个正整数![20140613112543656][]和![20140613112553093][]，求![20140613112622265][]的所有因子和对**9901**取余后的值。

**分析：**很容易知道，先把![20140613112543656][]分解得到![20140613113059203][]，那么得到![20140613113255234][]，那么![20140613112622265][]

的所有因子和的表达式如下

![20140613132246375][]

所以我们有两种做法。第一种做法是二分求等比数列之和。

**代码：**

#include <iostream>
     #include <string.h>
     #include <stdio.h>
    using namespace std;
    typedef long long LL;
    const int N = 10005;
    const int MOD = 9901;
    
    bool prime[N];
    int p[N];
    int cnt;
    
    void isprime()
    {
        cnt = 0;
        memset(prime,true,sizeof(prime));
        for(int i=2; i<N; i++)
        {
            if(prime[i])
            {
                p[cnt++] = i;
                for(int j=i+i; j<N; j+=i)
                    prime[j] = false;
            }
        }
    }
    
    LL power(LL a,LL b)
    {
        LL ans = 1;
        a %= MOD;
        while(b)
        {
            if(b & 1)
            {
                ans = ans * a % MOD;
                b--;
            }
            b >>= 1;
            a = a * a % MOD;
        }
        return ans;
    }
    
    LL sum(LL a,LL n)
    {
        if(n == 0) return 1;
        LL t = sum(a,(n-1)/2);
        if(n & 1)
        {
            LL cur = power(a,(n+1)/2);
            t = (t + t % MOD * cur % MOD) % MOD;
        }
        else
        {
            LL cur = power(a,(n+1)/2);
            t = (t + t % MOD * cur % MOD) % MOD;
            t = (t + power(a,n)) % MOD;
        }
        return t;
    }
    
    void Solve(LL A,LL B)
    {
        LL ans = 1;
        for(int i=0; p[i]*p[i] <= A; i++)
        {
            if(A % p[i] == 0)
            {
                int num = 0;
                while(A % p[i] == 0)
                {
                    num++;
                    A /= p[i];
                }
                ans *= sum(p[i],num*B) % MOD;
                ans %= MOD;
            }
        }
        if(A > 1)
        {
            ans *= sum(A,B) % MOD;
            ans %= MOD;
        }
        cout<<ans<<endl;
    }
    
    int main()
    {
        LL A,B;
        isprime();
        while(cin>>A>>B)
            Solve(A,B);
        return 0;
    }

第二种方法就是用等比数列求和公式，但是要用逆元。用如下公式即可

![20140613105646406][]

因为![20140613114857000][]可能会很大，超过**int**范围，所以在快速幂时要二分乘法。

**代码：**

#include <iostream>
     #include <string.h>
     #include <stdio.h>
    using namespace std;
    typedef long long LL;
    const int N = 10005;
    const int MOD = 9901;
    
    bool prime[N];
    int p[N];
    int cnt;
    
    void isprime()
    {
        cnt = 0;
        memset(prime,true,sizeof(prime));
        for(int i=2; i<N; i++)
        {
            if(prime[i])
            {
                p[cnt++] = i;
                for(int j=i+i; j<N; j+=i)
                    prime[j] = false;
            }
        }
    }
    
    LL multi(LL a,LL b,LL m)
    {
        LL ans = 0;
        a %= m;
        while(b)
        {
            if(b & 1)
            {
                ans = (ans + a) % m;
                b--;
            }
            b >>= 1;
            a = (a + a) % m;
        }
        return ans;
    }
    
    LL quick_mod(LL a,LL b,LL m)
    {
        LL ans = 1;
        a %= m;
        while(b)
        {
            if(b & 1)
            {
                ans = multi(ans,a,m);
                b--;
            }
            b >>= 1;
            a = multi(a,a,m);
        }
        return ans;
    }
    
    void Solve(LL A,LL B)
    {
        LL ans = 1;
        for(int i=0; p[i]*p[i] <= A; i++)
        {
            if(A % p[i] == 0)
            {
                int num = 0;
                while(A % p[i] == 0)
                {
                    num++;
                    A /= p[i];
                }
                LL M = (p[i] - 1) * MOD;
                ans *= (quick_mod(p[i],num*B+1,M) + M - 1) / (p[i] - 1);
                ans %= MOD;
            }
        }
        if(A > 1)
        {
            LL M = MOD * (A - 1);
            ans *= (quick_mod(A,B+1,M) + M - 1) / (A - 1);
            ans %= MOD;
        }
        cout<<ans<<endl;
    }
    
    int main()
    {
        LL A,B;
        isprime();
        while(cin>>A>>B)
            Solve(A,B);
        return 0;
    }

其实有些题需要用到![20140613122746828][]模![20140613122821000][]的所有逆元，这里![20140613122821000][]为奇质数。那么如果用快速幂求时间复杂度为![20140613123226296][]，

如果对于一个**1000000**级别的素数![20140613122821000][]，这样做的时间复杂度是很高了。实际上有![20140613123603656][]的算法，有一个递推式如下

![20140613123955765][]

它的推导过程如下，设![20140613124932968][]，那么

![20140613125025984][]

对上式两边同时除![20140613125121234][]，进一步得到

![20140613125208000][]

再把![20140613125307921][]和![20140613125324765][]替换掉，最终得到

![20140613123955765][]

初始化![20140613125520578][]，这样就可以通过递推法求出![20140613122746828][]模奇素数![20140613122821000][]的所有逆元了。

另外![20140613122746828][]模![20140613122821000][]的所有逆元值对应![20140613122746828][]中所有的数，比如![20140613125913078][]，那么![20140613130144843][]对应的逆元是![20140613130219921][]。

**题目：[http://www.lydsy.com/JudgeOnline/problem.php?id=2186][http_www.lydsy.com_JudgeOnline_problem.php_id_2186]**

**题意：**求![20140613132413593][]中![20140613132447468][]互质的数的个数，其中![20140613132544453][]。

**分析：**因为![20140613132544453][]，所以![20140613132643546][]，我们很容易知道如下结论

**   对于两个正整数![20140613132805796][]和![20140613132820031][]，如果![20140613132820031][]是![20140613132805796][]的倍数，那么![20140613132916781][]中与![20140613132805796][]互素的数的个数为![20140613132958781][]**

**    ** 本结论是很好证明的，因为![20140613133123640][]中与![20140613132805796][]互素的个数为![20140613133215171][]，又知道![20140613133304375][]，所以

结论成立。那么对于本题，答案就是

![20140613133456343][]

其中![20140613133637687][]为小于等于![20140613133717765][]的所有素数，先筛选出来即可。由于最终答案对一个质数取模，所以要用逆元，这里

求逆元就有技巧了，用刚刚介绍的递推法预处理，否则会**TLE**的。

**代码：**

#include <iostream>
     #include <string.h>
     #include <stdio.h>
     #include <bitset>
    using namespace std;
    typedef long long LL;
    const int N = 10000005;
    
    bitset<N> prime;
    
    void isprime()
    {
        prime.set();
        for(int i=2; i<N; i++)
        {
            if(prime[i])
            {
                for(int j=i+i; j<N; j+=i)
                    prime[j] = false;
            }
        }
    }
    
    LL ans1[N],ans2[N];
    LL inv[N];
    
    int main()
    {
        isprime();
        int MOD,m,n,T;
        scanf("%d%d",&T,&MOD);
        ans1[0] = 1;
        for(int i=1; i<N; i++)
            ans1[i] = ans1[i-1] * i % MOD;
        inv[1] = 1;
        for(int i=2;i<N;i++)
        {
            if(i >= MOD) break;
            inv[i] = (MOD - MOD / i) * inv[MOD % i] % MOD;
        }
        ans2[1] = 1;
        for(int i=2; i<N; i++)
        {
            if(prime[i])
            {
                ans2[i] = ans2[i-1] * (i - 1) % MOD;
                ans2[i] = ans2[i] * inv[i % MOD] % MOD;
            }
            else
            {
                ans2[i] = ans2[i-1];
            }
        }
        while(T--)
        {
            scanf("%d%d",&n,&m);
            LL ans = ans1[n] * ans2[m] % MOD;
            printf("%lld\n",ans);
        }
        return 0;
    }

接下来还有一个关于逆元的有意思的题目，描述如下

![20140613140842765][]

**证明：**由

![20140613141419609][]

其中![20140613141505140][]

所以只需要证明![20140613141749171][]，而我们知道![20140613141943718][]模![20140613142014640][]的逆元对应全部

![20140613141943718][]中的所有数，既是单射也是满射。

所以进一步得到

![20140613142541046][]

证明完毕！

[20140613102654328]: https://img-blog.csdn.net/20140613102654328
[20140613102712781]: /images/20220611/5e196e5fcc9c41a4a98656edebf73269.png
[20140613102734984]: https://img-blog.csdn.net/20140613102734984
[20140613102856531]: https://img-blog.csdn.net/20140613102856531
[20140613103413828]: /images/20220611/bbb718e13c414e60a110c834bcb67cc4.png
[20140613104248984]: /images/20220611/4452ff03098749838ff1b7027f9b4e42.png
[20140613104752312]: /images/20220611/fa362e62926d4d48a2d5793789ecc1e0.png
[20140613104619203]: /images/20220611/df2765d1c0524510bf22d768e5d9ed48.png
[20140613105646406]: /images/20220611/59101309fbea4cafb0d82e63f12e5655.png
[20140613112023328]: /images/20220611/41457505d0a54a02b2e0b1fb5305395d.png
[http_poj.org_problem_id_1845]: http://poj.org/problem?id=1845
[20140613112543656]: https://img-blog.csdn.net/20140613112543656
[20140613112553093]: https://img-blog.csdn.net/20140613112553093
[20140613112622265]: /images/20220611/aa19ce7d748f4e528e62695ecc8bbe15.png
[20140613113059203]: https://img-blog.csdn.net/20140613113059203
[20140613113255234]: https://img-blog.csdn.net/20140613113255234
[20140613132246375]: /images/20220611/d8a45e4a8c954a6dbdc3822d64be9c67.png
[20140613114857000]: /images/20220611/60728ff40baa4aa58842906adfeecfda.png
[20140613122746828]: https://img-blog.csdn.net/20140613122746828
[20140613122821000]: /images/20220611/b66d3c8838aa452a949584ae3b1a252c.png
[20140613123226296]: /images/20220611/3f4e85dcd51d4a98af4cde96471b219e.png
[20140613123603656]: /images/20220611/6ff4f81dcb5a4f2e92ad3595d7c1b304.png
[20140613123955765]: /images/20220611/0d68db06c55e4c2d968c59c7dc7a9e1f.png
[20140613124932968]: /images/20220611/487ed2cb29e74f5caa1c4272c0a27d7a.png
[20140613125025984]: /images/20220611/251699e15aa24aeebeba38f9135864a9.png
[20140613125121234]: /images/20220611/9a71b4e371624e54b32fab8e4393de36.png
[20140613125208000]: /images/20220611/826d162175d4464aa0785ef1eff16ff8.png
[20140613125307921]: https://img-blog.csdn.net/20140613125307921
[20140613125324765]: /images/20220611/ae2b1a663e9942e88e6f5a0be469e194.png
[20140613125520578]: https://img-blog.csdn.net/20140613125520578
[20140613125913078]: https://img-blog.csdn.net/20140613125913078
[20140613130144843]: https://img-blog.csdn.net/20140613130144843
[20140613130219921]: /images/20220611/cdeb1312cce2455daf5d4db39d4cb5d8.png
[http_www.lydsy.com_JudgeOnline_problem.php_id_2186]: http://www.lydsy.com/JudgeOnline/problem.php?id=2186
[20140613132413593]: https://img-blog.csdn.net/20140613132413593
[20140613132447468]: https://img-blog.csdn.net/20140613132447468
[20140613132544453]: /images/20220611/521552404a424fa083d6f6cbbb4fd811.png
[20140613132643546]: /images/20220611/13bba8f3c5ac4ca09f1d394a659f196b.png
[20140613132805796]: https://img-blog.csdn.net/20140613132805796
[20140613132820031]: https://img-blog.csdn.net/20140613132820031
[20140613132916781]: https://img-blog.csdn.net/20140613132916781
[20140613132958781]: /images/20220611/0691deb00d454c56a795974c1b52459a.png
[20140613133123640]: https://img-blog.csdn.net/20140613133123640
[20140613133215171]: https://img-blog.csdn.net/20140613133215171
[20140613133304375]: /images/20220611/7b1995aeb6be463d8382f7e14ef808da.png
[20140613133456343]: /images/20220611/14ec9b864dff4457bf06981720278c0a.png
[20140613133637687]: https://img-blog.csdn.net/20140613133637687
[20140613133717765]: /images/20220611/e20079af04e7419db4cd8748ff29f064.png
[20140613140842765]: /images/20220611/7231944e50d34fa6b1e29c8389ffd718.png
[20140613141419609]: /images/20220611/ee5e1574021b4bccb07ca38e605a9b4f.png
[20140613141505140]: /images/20220611/9c8b223907b94d86a1159e5de97d5b1c.png
[20140613141749171]: https://img-blog.csdn.net/20140613141749171
[20140613141943718]: /images/20220611/8692566efb9445829f360e93f0a9e965.png
[20140613142014640]: /images/20220611/2fdb5e836b4b4f6b8c78b7ef3533d291.png
[20140613142541046]: /images/20220611/75b4aa2d08614825aca871582ac26b21.png