HDU 6053(莫比乌斯反演+快速幂)

来源:互联网 发布:ugg和jumbougg知乎 编辑:程序博客网 时间:2024/06/04 19:06

       (ps:声明一下,这是一篇误人子弟的文章2333//手动滑稽)

   首先讲一下题目思路,他要找出满足两个条件的B数组的个数,显然要枚举GCD,然后求出每个(Ai/GCD)的乘积,然后累加,

对于某个数及其倍数能产生的B数组方案个数为:


因为过程中会出现重复,所以要用容斥原理,然后就要用到强大的莫比乌斯反演了,(莫比乌斯反演的推导)。


本题的难点还有一个就是如何快速求得(Ai/gcd)的乘积,我们可以维护一个前缀和数组b[],b[i]表示有多少Ai<=i,

然后我们类似筛法一样的枚举,假设当前枚举的gcd为d,那么j * d - 1 到 (j + 1) * d - 1 除以d的商都是j,结合b[]数组

我们就可以快速求出除以d商为j的Ai有几个了,这样复杂度和筛法差不多,nlogn数量级。

#include  #define ll long long  #define pb push_back  #define fi first  #define se second  #define pi acos(-1)  #define inf 0x3f3f3f3f  #define lson l,mid,rt<<1  #define rson mid+1,r,rt<<1|1  #define rep(i,x,n) for(int i=x;i=x;i--)  using namespace std;  const int mod = 1e9 + 7;  typedef pairP;  const int MAXN=100010;  int gcd(int a,int b){return b?gcd(b,a%b):a;}  int f[MAXN];//f为moebius函数值   ll F[MAXN];   void moebius()  {      f[1] = 1;      for(int i = 2; i < MAXN; i++)      {          if(f[i] == 0)          {              for(int j = i; j < MAXN; j += i)              f[j]++;          }      }      for(int i = 2; i < MAXN; i++)      if(f[i] & 1)      f[i] = -1;      else      f[i] = 1;      for(int i = 2; i * i < MAXN; i++)      {          for(int j = i * i; j < MAXN; j += i * i)          f[j] = 0;      }  }  int a[MAXN];  int b[MAXN * 2];//前缀和数组   ll qmul(int a, int b)  {      ll ans = 1, t = a;      while(b)      {          if(b & 1) ans *= t;          t *= t;          b >>= 1;          ans %= mod;          t %= mod;      }      return ans;  }  int main()  {      int T;      moebius();      cin >> T;      for(int kase = 1; kase <= T; kase++)      {          int n;          scanf("%d", &n);          memset(b, 0, sizeof(b));          for(int i = 0; i < n; i++) scanf("%d", a + i), b[a[i]]++;          for(int i = 0; i < 2 * MAXN; i++) b[i] += b[i - 1];          int up = *min_element(a, a + n);          ll ans = 0;          for(int i = 2; i <= up; i++)//枚举gcd           {              ll num = 1;              for(int j = 1; j * i <= MAXN; j++)              {                  num *= qmul(j, b[(j + 1) * i - 1] - b[j * i - 1]);                  num %= mod;              }              F[i] = num;           }          for(int i = 2; i <= up; i++)          {              ll num = 1;              for(int j = 1; j * i <= MAXN; j++)              ans += f[j] * F[i * j], ans %= mod;           }          printf("Case #%d: %lld\n", kase, ans);      }       return 0;  }