poj-2480(数论+积性函数性质)

男娘i 2022-06-10 14:23 167阅读 0赞

Description

Longge is good at mathematics and he likes to think about hard mathematical problems which will be solved by some graceful algorithms. Now a problem comes: Given an integer N(1 < N < 2^31),you are to calculate ∑gcd(i, N) 1<=i <=N.

"Oh, I know, I know!" Longge shouts! But do you know? Please solve it.

**Input**

Input contain several test case.   
A number N per line.

**Output**

For each N, output ,∑gcd(i, N) 1<=i <=N, a line

**Sample Input**

2
    6

**Sample Output**

3
    15

这个大牛写得很好！

[点击打开链接][Link 1]

题目题意:题目给我们一个数字n，让我们计算1到n的每个数i 与n的gcd(i,n)的和

题目分析:这道题目得用到一种数论函数（积性函数）

积性函数指对于所有互质的整数a和b有性质f(ab)=f(a)f(b)的 [数论函数][Link 2] 。

百度一下我们可以知道很多常用的积性函数，就这道题目来说，我们只用知道gcd(i,k)是积性函数

即k=n\*m   gcd(i,k)=gcd(i,n)\*gcd(i,m)

∑gcd(i, N)这个函数也是积性函数。

我们对N质因子分解

N=p1^a1\*p2^a2\*p3^a3\*\*\*\*\*pn^an

这根据积性函数的性质:∑gcd(i, N)=∑gcd(i, p1^a1)\*∑gcd(i, p2^a2)\*∑gcd(i, p3^a3)\*\*\*\*\*\*\*\*∑gcd(i, pn^an)

现在我们就是求每一个∑gcd(i, pk^ak)。

很明显gcd(i,pk^ak)只可能是pk^at次方(0<=at<=ak)

假设gcd(i,pk^ak)=pk^at 那么i的个数为phi\[pk^ak/pk^at\]

我们拿gcd(i,a)=b证明

即i=k\*b  a=p\*b 使得gcd(k\*b,p\*b)=b, k与p互质的对数就是phi\[a/b\]

则我们的到来一个式子的答案：

∑gcd(i, pk^ak)=pk^0\*phi\[pk^(ak-0)\]+pk^1\*phi\[pk^(ak-1)+.....+pk^ak\*phi\[pk^(ak-ak)\]

我然后就写了代码，然后TLE了！

其实还可以变形:

对于phi\[pk^ak\],因为pk是质数的原因，那么能与不互质的一定是pk的倍数

则phi\[pk^ak\]=pk^ak-pk^(ak-1)

然后代入：

我们得到了

pk^ak\*(1+ak\*(1-1/pk)

再进一步化简：

pk^ak+ak\*pk^ak-ak\*p^(ak-1)

然后就可以做了.

代码如下:

#include<iostream>
    #include<cstdio>
    #include<cmath>
    #include<cstring>
    #define ll long long
    using namespace std;
    
    const int maxn=1e6+10;
    int prime[maxn],cnt;
    bool vis[maxn];
    void get_prime()//打表
    {
        memset (vis,true,sizeof (vis));
        vis[1]=false;
        for (int i=2;i<maxn;i++) {
            if (vis[i]) prime[cnt++]=i;
            for (int j=0;j<cnt&&i*prime[j]<maxn;j++) {
                vis[i*prime[j]]=false;
                if (i%prime[j]==0) break;
            }
        }
    }
    int factor[1005][2],fcnt;
    void get_factor(int n)//分解质因子
    {
        memset (factor,0,sizeof (factor));
        fcnt=0;
        for (int i=0;i<cnt&&prime[i]*prime[i]<=n;i++) {
            if (n%prime[i]==0) {
                factor[fcnt][0]=prime[i];
                while (n%prime[i]==0) {
                    factor[fcnt][1]++;
                    n=n/prime[i];
                }
                fcnt++;
            }
        }
        if (n!=1) {
            factor[fcnt][0]=n;
            factor[fcnt++][1]=1;
        }
    }
    
    int main()
    {
        get_prime();
        int n;
        while (scanf("%d",&n)!=EOF) {
            get_factor(n);
            ll ans=1;
            for (int i=0;i<fcnt;i++) {
                ll res=pow(factor[i][0],factor[i][1]);//套入公式里面
                ans*=(res+factor[i][1]*res-factor[i][1]*res/factor[i][0]);
            }
            printf("%lld\n",ans);
        }
        return 0;
    }

[Link 1]: http://lydws.blog.163.com/blog/static/22621105120152265175340/
[Link 2]: https://baike.baidu.com/item/%E6%95%B0%E8%AE%BA%E5%87%BD%E6%95%B0