zoukankan      html  css  js  c++  java
  • 【题解】HAOI2011Problem b

      第一次接触莫比乌斯反演,总之脑子都快要炸掉了……好难啊!本蒟蒻表示根本无法理解呜呜呜呜呜……不过在机房DL的帮助下总算是磕磕绊绊的A掉了这一题:

      这道题目要我们的求的是:(1)ΣiΣj [gcd(i,j)==k], 区间范围内的限定我们都可以利用容斥来解决,所以默认为所求函数值的(i, j)取值范围为(1~n,1~m)。注意到gcd(i, j) == k比较的不好求解,我们把k除到外面去,得到:(2)ΣiΣj [gcd(i,j)==1](i∈n/k, j∈m/k); 这里将[gcd(i,j)==1]替换为单位元ε(ε(i) = i == 1 ? 1 : 0)。根据狄利克雷卷积有μ * 1 = ε(*为卷积符号),所以我们得到式子(3)ΣiΣjΣd|gcd(i,j)μ(d)。这里是考虑对于每一对i,j有哪几个d与之对应,我们将d提取出来,反之枚举符合d|gcd(i,j)的数对。由此有(4)Σd μ(d) * (n / (d*k)) * (m / (d * k))(d∈min(n, m) / (d * k))。

      但是利用最后这个式子求解的复杂度依然太高了,不过因为n/i的取值最多只有2*sqrt(n)种(若i<sqrt(n),i的取值只有sqrt(n)种,若i>=sqrt(n),n/i<=sqrt(n),所以一共最多只有2*sqrt(n)种),所以在很大的范围内实际上(n / (d*k)) * (m / (d * k))都是一个常数,我们就可以利用μ(d)的前缀和来实现快速求解。

      设当前枚举到的d,pos = min(n / (n / d), m / (m / d)); 这个pos就是从当前位置一直到达pos,(n / (d*k)) * (m / (d * k))都是一个常数,只需要计算一次之后利用前缀和求解。

      贴代码,再一次被自己蠢哭,智商堪忧怎么办啊QAQ(向数学组DL低头,伏地%%%,跪地%%%!!!)

    #include <bits/stdc++.h>
    using namespace std;
    #define maxn 65000
    int tot, T, maxx = 60000, sum[maxn], pri[maxn], Mu[maxn];
    bool is_prime[maxn];
    
    int read()
    {
        int x = 0, k = 1;
        char c;
        c = getchar();
        while(c < '0' || c > '9') { if(c == '-') k = -1; c = getchar(); }
        while(c >= '0' && c <= '9') x = x * 10 + c - '0', c = getchar();
        return x * k;
    }
    
    int Get_Mu()
    {
        Mu[1] = 1;
        for(int i = 2; i <= maxx; i ++)
        {
            if(!is_prime[i]) pri[++ tot] = i, Mu[i] = -1;
            for(int j = 1; j <= tot; j ++)
            {
                if(i * pri[j] > maxx) break;
                is_prime[i * pri[j]] = true;
                if(!(i % pri[j]))
                {
                    Mu[i * pri[j]] = 0;
                    break;
                }
                else Mu[i * pri[j]] = - Mu[i];
            }
        }
        for(int i = 1; i <= maxx; i ++)
            sum[i] = sum[i - 1] + Mu[i];
    }
    
    int cal(int n, int m)
    {
        if(n > m) swap(n, m);
        int ans = 0, pos;
        for(int i = 1; i <= n; i = pos + 1)
        {
            pos = min(n / (n / i), m / (m / i));
            ans += (sum[pos] - sum[i - 1]) * (n / i) * (m / i);
        }
        return ans;
    }
    
    int main()
    {
        Get_Mu();
        T = read();
        while(T --)
        {
            int a = read(), b = read(), c = read(), d = read(), k = read();
            a --, c --;
            a /= k, b /= k, c /= k, d /= k;
            int ans = cal(a, c) + cal(b, d) - cal(a, d) - cal(b, c);
            printf("%d
    ", ans);
        }
        return 0;
    } 
  • 相关阅读:
    redis可编译
    不要用Serverzoo 提供的CloudLinux 的五大原因 Linode 強大VPS 資源為你解密
    linux加载指定目录的so文件
    超级rtmp服务器和屌丝wowza
    标准IO: 文件的打开与关闭函数 fopen & fclose
    《gdb调试之基础篇》
    linux信号Linux下Signal信号太详细了,终于找到了
    【干货】Chrome插件(扩展)开发全攻略
    斯坦福开源无Bug的随机计算图Certigrad
    心跳包:告诉别人,我还活着
  • 原文地址:https://www.cnblogs.com/twilight-sx/p/8470967.html
Copyright © 2011-2022 走看看