zoukankan      html  css  js  c++  java
  • poj2480

    题意:给出n(2^31),求∑gcd(i, n) 1<=i <=n。

    分析:首先所有与n互质的数字x都满足gcd(x,n)=1。我们先计算这种等于1的情况,恰好是n的欧拉函数phi(n)。我们可以将上述情况视为跟n最大公约数为1的情况,现在我们将其推广至最大公约数为p的情况。对于对于所有满足gcd(x,n)=p(p为常数)的x,他们与n拥有相同的gcd,那么他们同时除以p之后,就会变得和n/p互质,这种数字x有phi(n/p)个,这些x的和为p*phi(n/p)个。所以我们要计算∑gcd(i, n) 1<=i <=n,只需要根据gcd的值不同,分类进行计算即可,总结成公式:∑p*phi(n/p)(p是n的约数)

    对于这个公式我们有一个快速计算的方法,我们先想,若要计算∑p,我们可以用公式(x1^0+x1^1+...+x1^c1)*(x2^0+x2^1+...+x2^c2)*...*(xm^0+xm^1+...+xm^cm)。xi是n的质因数,ci是n中含有xi的个数。同理,根据欧拉函数的公式,每个欧拉函数值也和各个质因子有关,求∑phi(n/p)=∑phi(p)=(1 + x1^0 * (x1 - 1) + x1^1 * (x1- 1) + ...+x1^(c1-1) * (x1- 1))*...*(...)。然后由于二者都是质因子幂相乘的形式,所以可以把两者综合起来,对于每个括号里的每一项,都是约数因式和欧拉函数因式的乘积的形式,即p的因式要乘以phi(n/p)的因式,两因式的质因子指数是互补的,最终公式变为(x1^(c1-1)*(c1-1) + x1^c1)*...。

    View Code
    #include <iostream>
    #include <cstdlib>
    #include <cstdio>
    #include <cstring>
    #include <cmath>
    using namespace std;
    
    int n;
    
    bool is[(1 << 16) + 5];
    int prm[(1 << 16) + 5];
    
    int getprm(int n)
    {
        int i, j, k = 0;
        int s, e = (int) (sqrt(0.0 + n) + 1);
        memset(is, 1, sizeof(is));
        prm[k++] = 2;
        is[0] = is[1] = 0;
        for (i = 4; i < n; i += 2)
            is[i] = 0;
        for (i = 3; i < e; i += 2)
            if (is[i])
            {
                prm[k++] = i;
                for (s = i * 2, j = i * i; j < n; j += s)
                    is[j] = 0;
            }
        for (; i < n; i += 2)
            if (is[i])
                prm[k++] = i;
        return k;
    }
    
    long long power(int p, int cnt)
    {
        if (cnt < 0)
            return 1;
        long long ret = 1;
        long long temp = p;
        while (cnt)
        {
            if (cnt & 1)
                ret *= temp;
            cnt >>= 1;
            temp = temp * temp;
        }
        return ret;
    }
    
    long long cal(int p, int cnt)
    {
        return power(p, cnt - 1) * cnt * (p - 1) + power(p, cnt);
    }
    
    int main()
    {
        //freopen("t.txt", "r", stdin);
        int m = getprm(1 << 16);
        long long ans;
        while (~scanf("%d", &n))
        {
            ans = 1;
            for (int i = 0; i < m && prm[i] * prm[i] <= n; i++)
            {
                int cnt = 0;
                while (n % prm[i] == 0)
                {
                    n /= prm[i];
                    cnt++;
                }
                ans *= cal(prm[i], cnt);
            }
            if (n != 1)
                ans *= cal(n, 1);
            printf("%lld\n", ans);
        }
        return 0;
    }

     

  • 相关阅读:
    Windows Azure 网站开发Stacks支持
    AzureDev 社区活动获奖者公布
    Android 改变窗口标题栏的布局
    cocos2d-x游戏开发系列教程-超级玛丽01-前言
    cocos2dx进阶学习之CCObject
    基于visual Studio2013解决算法导论之055拓扑排序
    查看某文件夹内文件大小
    vmstat命令
    uname 命令
    iostat命令
  • 原文地址:https://www.cnblogs.com/rainydays/p/2745672.html
Copyright © 2011-2022 走看看