Summer Training day4 Code hdu5212 莫比乌斯反演

来源:互联网 发布:知乎日报 离线下载 编辑:程序博客网 时间:2024/06/15 06:16

WLD likes playing with codes.One day he is writing a function.Howerver,his computer breaks down because the function is too powerful.He is very sad.Can you help him? 

The function: 


int calc 

   
  int res=0; 
   
  for(int i=1;i<=n;i++) 
     
    for(int j=1;j<=n;j++) 
     
    { 
       
      res+=gcd(a[i],a[j])*(gcd(a[i],a[j])-1); 
       
      res%=10007; 
     
    } 
   
  return res; 

}
Input
There are Multiple Cases.(At MOST 1010

For each case: 

The first line contains an integer N(1N10000)N(1≤N≤10000)

The next line contains NN integers a1,a2,...,aN(1ai10000)a1,a2,...,aN(1≤ai≤10000)
Output
For each case: 

Print an integer,denoting what the function returns.
Sample Input
51 3 4 2 4
Sample Output
64          
Hint
gcd(x,y) means the greatest common divisor of x and y.


方法1:莫比乌斯 反演

统计1-10000中每个数出现的次数并存在num数组里面

从先枚举i和j变成先枚举gcd

其中sum数组表示的含义是

sum(i) = num(i)+num(2i)+...+num(ni)


代码:

#include <iostream>#include <cstdio>#include <cstring>using namespace std;typedef long long LL;const int MAX = 10007;const int MOD = 10007;int nums[MAX];int sums[MAX];const int maxn = 222232;  int mu[maxn],a[maxn],prime[maxn],cnt[maxn],num[maxn];  bool vis[maxn];  int n,pnum;  int aa,b,c,d,k;void mobeius(int N)  {      pnum=0;      vis[1]=mu[1]=1;      for(int i=2;i<=N;i++)      {          if(!vis[i])          {              mu[i]=-1;              prime[pnum++]=i;          }          for(int j=0;j<pnum;j++)          {              if(i*prime[j]>N)break;              vis[i*prime[j]]=1;//ɸµôºÏÊý              if(i%prime[j]==0)              {                  mu[i*prime[j]]=0;                  break;//±£Ö¤ºÏÊýʹÓÃ×îСµÄËØÊýɸµôµÄ              }              mu[i*prime[j]]=-mu[i];          }      }  }  int main(){mobeius(MAX);int N;while(scanf("%d",&N) != EOF){memset(nums,0,sizeof(nums));memset(sums,0,sizeof(sums));for(int i = 0 ;i < N;i++){int tmp;scanf("%d",&tmp);nums[tmp] ++;}for(int i = 1;i <= 10000;i++){for(int j = i;j <= 10000;j += i){sums[i] += nums[j];}sums[i] %= MOD;}LL ans = 0;for(int n = 1;n <= 10000;n++){LL g = 0;for(int j = n;j <= 10000;j += n){g = (g + mu[j/n] * sums[j] % MOD * sums[j]) % MOD;}ans = (ans + g * n % MOD * (n-1)) % MOD;}printf("%lld\n",ans);}return 0;}

方法2:


依然同方法一一样处理,只不过在到达这一步以后我们可以方便的求出后面一半的值

即这个式子的值

我们定义上面的式子为F[k]

那么我们如果求F的时候将k从大到小递减,那么F[i] = sum(i)*sum(i) - F[2i] -F[3i] -.....

这里sum(i)的定义与方法一中的定义一致


代码:

#include <iostream>#include <cstdio>#include <algorithm>#include <cstring>using namespace std;#define int long longconst int MOD = 10007;const int MAX = 10006;int N;int num[MAX];int F[MAX];main(){while(scanf("%lld",&N) != EOF){int ma = 0;memset(num,0,sizeof(num));memset(F,0,sizeof(F));for(int i = 0;i < N;i++){int tmp;scanf("%lld",&tmp);num[tmp] ++;ma = max(ma,tmp);}int ans = 0;for(int i = ma;i >= 1;i--){int tmp = 0;int sum = 0;for(int j = i;j <= ma;j += i){sum = (sum + num[j]) % MOD;tmp = (F[j] + tmp) % MOD;}F[i] = (sum * sum % MOD - tmp + MOD )% MOD;ans += F[i] * i % MOD * (i-1) % MOD;}printf("%lld\n",ans % MOD);}return 0;}



原创粉丝点击