zoukankan      html  css  js  c++  java
  • Problem b 莫比乌斯反演+枚举除法的取值

    莫比乌斯反演+枚举除法的取值
    第二种形式:
    f(n)表示gcd(x,y)=n的数量。
    F(n)表示gcd(x,y)是n的倍数的数量。
    /**
    题目:Problem b
    链接:https://vjudge.net/contest/178455#problem/G
    题意:对于给出的 n 个询问,每次求有多少个数对 (x,y) ,
    满足 a ≤ x ≤ b , c ≤ y ≤ d ,且 gcd(x,y) = k , gcd(x,y) 函数为 x 和 y 的最大公约数。
    1≤n≤50000,1≤a≤b≤50000,1≤c≤d≤50000,1≤k≤50000
    思路:
    首先容斥:ans = solve(b,d,k)-solve(b,c-1,k)-solve(a-1,d,k)+solve(a-1,c-1,k);
    
    solve(n,m,k)表示x在[1,n],y在[1,m] gcd(x,y)==k的对数。
    
    定义:
    f(n)表示gcd(x,y)=n的数量。
    F(n)表示gcd(x,y)是n的倍数的数量。
    
    如何求F(n)?
    
    F(n) = (x/n) * (y/n); 要加括号,因为这是取整之后的乘积
    
    根据定义用第二种形式:f(n) = sigma(mu[d/n]*F(d)) (n|d)
    
    这样只要枚举k的倍数一直到min(n,m)就可以了。可是如果k=1,那么枚举一次就是O(N);总复杂度为O(N*N);
    
    
    实际上可以继续优化;
    
    solve(n,m,k)等价于solve(n/k,m/k)表示x在[1,n/k],y在[1,m/k],gcd(x,y)==1的对数。
    
    由于x/i,x/(i+1),x/(i+2)...x/(i+t)存在连续相同的结果,也就是这段区间[l,r]内(n/i)*(m/i)的结果是相同的;
    
    这样i在[l,r] 范围内的(n/i)*(m/i)*mu[i];就等价于 (n/i)*(m/i)*(sum[r]-sum[l-1]); sum表示mu的前缀和。
    
    所以这里可以快速处理。复杂度为sqrt(N); 总时间复杂度为N*sqrt(N);
    
    参考:https://wenku.baidu.com/view/fbec9c63ba1aa8114431d9ac.html
    
    */
    #include <cstdio>
    #include <cstring>
    #include <algorithm>
    #include <set>
    #include <iostream>
    #include <vector>
    #include <map>
    using namespace std;
    typedef long long LL;
    #define ms(x,y) memset(x,y,sizeof x)
    typedef pair<int, int> P;
    const LL INF = 1e10;
    const int mod = 1e9 + 7;
    const int maxn = 5e4 + 100;
    int prime[maxn], tot, not_prime[maxn];
    int mu[maxn], sum[maxn];
    void init()
    {
        mu[1] = 1;
        tot = 0;
        for(int i = 2; i < maxn; i++){
            if(!not_prime[i]){
                prime[++tot] = i;
                mu[i] = -1;
            }
            for(int j = 1; prime[j]*i<maxn; j++){
                not_prime[prime[j]*i] = 1;
                if(i%prime[j]==0){
                    mu[prime[j]*i] = 0;
                    break;
                }
                mu[prime[j]*i] = -mu[i];
            }
        }
        for(int i = 1; i < maxn; i++) sum[i] = sum[i-1]+mu[i];
    }
    LL solve(int n,int m)
    {
        LL ans = 0;
        if(n>m) swap(n,m);
        int last;
        for(int i = 1; i <= n; i=last+1){
            last = min(n/(n/i),m/(m/i));
            ans += (LL)(sum[last]-sum[i-1])*(n/i)*(m/i);
        }
        return ans;
    }
    int main()
    {
        //freopen("in.txt","r",stdin);
        int T;
        int a, b, c, d, k;
        init();
        cin>>T;
        while(T--)
        {
            scanf("%d%d%d%d%d",&a,&b,&c,&d,&k);
            printf("%lld
    ",solve(b/k,d/k)-solve(b/k,(c-1)/k)-solve((a-1)/k,d/k)+solve((a-1)/k,(c-1)/k));
        }
        return 0;
    }
  • 相关阅读:
    10 个深恶痛绝的 Java 异常。。
    为什么公司宁愿 25K 重新招人,也不给你加到 20K?原因太现实……
    推荐一款代码神器,代码量至少省一半!
    Spring Cloud Greenwich 正式发布,Hystrix 即将寿终正寝。。
    hdu 3853 LOOPS(概率 dp 期望)
    hdu 5245 Joyful(期望的计算,好题)
    hdu 4336 Card Collector(期望 dp 状态压缩)
    hdu 4405 Aeroplane chess(概率+dp)
    hdu 5036 Explosion(概率期望+bitset)
    hdu 5033 Building (单调栈 或 暴力枚举 )
  • 原文地址:https://www.cnblogs.com/xiaochaoqun/p/7360155.html
Copyright © 2011-2022 走看看