hdu1695 容斥原理 莫比乌斯反演

给定两个数b，d，问[1,b]和[1,d]区间上有多少对互质的数。（x,y）和（y,x）算一个。
对于[1,b]部分，用欧拉函数直接求。对于大于b的部分，求n在[1,b]上有多少个互质的数，用容斥原理。
主要学习容斥原理的写法，本题使用DFS。容斥原理复杂度比较高，是指数复杂度。
输出长整型不能用lld，必须用I64d。

#include<stdio.h>
#include<iostream>
#include<stdlib.h>
#include<string.h>
#include<algorithm>
using namespace std;
typedef long long ll;
const int maxn = 1e5 + 7;
int a, b, c, d, k;
bool isPrime[maxn];
ll euler[maxn];
ll sumeuler[maxn];
int p[maxn][20], ps[maxn];
void initPrimes(){
    memset(isPrime, 1, sizeof(isPrime));
    for (int i = 1; i < maxn; i++){
        euler[i] = i;
    }
    for (int i = 2; i < maxn; i++){
        if (isPrime[i]){
            for (int j = i; j < maxn; j += i){
                isPrime[j] = false;
                p[j][ps[j]++] = i;
                euler[j] = euler[j] * (i - 1) / i;
            }
        }
    }
    sumeuler[1] = 1;
    for (int i = 2; i < maxn; i++){
        sumeuler[i] = sumeuler[i - 1] + euler[i];
    }
}
ll dfs(int ind, int n, int x){//容斥原理核心代码
    ll s = 0;
    for (int i = ind; i < ps[x]; i++){
        s += n / p[x][i] - dfs(i + 1, n / p[x][i], x);
    }
    return s;
}
ll huzhi(int n, int x){//0~n之间，与x互质的数字个数
    return n - dfs(0, n, x);
}
int main(){
    freopen("in.txt", "r", stdin);
    initPrimes();
    int caseCount;
    scanf("%d", &caseCount);
    for (int caseid = 1; caseid <= caseCount; caseid++){
        scanf("%d%d%d%d%d", &a, &b, &c, &d, &k);
        ll ans;
        if (k == 0)
            ans = 0;
        else{
            if (b > d)swap(b, d);
            b /= k, d /= k;
            ans = sumeuler[b];
            for (int i = b + 1; i <= d; i++){
                ans += huzhi(b, i);
            }
        }
        printf("Case %d: %I64d\n", caseid, ans);
    }
    return 0;
}

容斥原理的另一种写法：

int calc(int n,int m)//n < m,求1-n内和m互质的数的个数
{
    getFactors(m);
    int ans = 0;
    for(int i = 1;i < (1<<fatCnt);i++)
    {
        int cnt = 0;
        int tmp = 1;
        for(int j = 0;j < fatCnt;j++)
            if(i&(1<<j))
            {
                cnt++;
                tmp *= factor[j][0];
            }
        if(cnt&1)ans += n/tmp;
        else ans -= n/tmp;
    }
    return n - ans;
}

容斥原理项的个数为2的幂次，肯定不会太大，所以一定可以用一个int来表示所有情况。

本题还可以用莫比乌斯反演来解决。

#include <stdio.h>
#include <string.h>
#include <iostream>
#include <algorithm>
#include <vector>
#include <queue>
#include <set>
#include <map>
#include <string>
#include <math.h>
#include <stdlib.h>
#include <time.h>
using namespace std;
const int MAXN = 1e5 + 7;
typedef long long ll;
bool isprime[MAXN];
int prime[MAXN], psize;
int mu[MAXN];
ll euler[MAXN], sumeuler[MAXN];
void init(){
    memset(isprime, false, sizeof(isprime));
    mu[1] = euler[1] = sumeuler[1] = 1;
    psize = 0;
    for (int i = 2; i <= MAXN; i++) {
        if (!isprime[i]) {
            prime[psize++] = i;
            mu[i] = -1;
            euler[i] = i - 1;
        }
        for (int j = 0; j < psize&&i*prime[j] < MAXN; j++) {
            isprime[i * prime[j]] = true;
            if (i % prime[j] == 0) {
                mu[i * prime[j]] = 0;
                euler[i*prime[j]] = euler[i] * prime[j];
                break;
            }
            else {
                mu[i * prime[j]] = -mu[i];
                euler[i*prime[j]] = euler[i] * (prime[j] - 1);
            }
        }
        sumeuler[i] = sumeuler[i - 1] + euler[i];
    }
}
ll eu(int b){
    long long ans2 = 0;
    for (int i = 1; i <= b; i++)
        ans2 += (long long)mu[i] * (b / i)*(b / i);
    return ans2;
}
void debug(){
    for (int i = 1; i < MAXN; i++){ 
        if (eu(i) / 2 != sumeuler[i] - 1)
            cout <<i<<" "<< eu(i) / 2 << " " << sumeuler[i] - 1 << endl;
    }
    exit(0);
}
int main(){
    freopen("in.txt", "r", stdin);
    int T;
    int a, b, c, d, k;
    init();
    //debug();
    scanf("%d", &T);
    int iCase = 0;
    while (T--)
    {
        iCase++;
        scanf("%d%d%d%d%d", &a, &b, &c, &d, &k);
        if (k == 0) {
            printf("Case %d: 0\n", iCase);
            continue;
        }
        b /= k; d /= k;
        if (b > d)swap(b, d);
        long long ans1 = 0;
        for (int i = 1; i <= b; i++)
            ans1 += (long long)mu[i] * (b / i)*(d / i);
        //这里ans2表示重复的部分，这部分可以用欧拉函数直接求，完全不需要for循环，但是提交却是wa，经过本机检测，完全没有问题
        long long ans2 = 0;
        for (int i = 1; i <= b; i++)
            ans2 += (long long)mu[i] * (b / i)*(b / i); 
        ans1 -= ans2 / 2;
        //ans1 -= (sumeuler[b] - 1);
        printf("Case %d: %I64d\n", iCase, ans1);
    }
    return 0;
}

转载于:https://www.cnblogs.com/weiyinfu/p/6369018.html