Petrozavodsk Programming Camp Winter 2018 Day 6: Grand Prix of Gomel K. Kids Aren‘t Alright（容斥好题）

Code92007

于 2024-06-27 01:13:40 发布

阅读量363

点赞数 5

分类专栏：组合数学(容斥原理) 文章标签：容斥

本文链接：https://blog.csdn.net/Code92007/article/details/139897554

版权

组合数学(容斥原理) 专栏收录该内容

51 篇文章 1 订阅

订阅专栏

题目

19313번: Kids Aren't Alright

求非空集合的数量，满足gcd=1，lcm=m(m<=1e18)，答案取模998244353

思路来源

稲葉廻群、b站绿导师原谅你了容斥视频等

Grand Prix of Gomel (ITMO Day 5) - EOJ Wiki

题解1

先pollard_pho找因子，

然后f_x 表示 gcd=1, lcm=x 的组个数

我们仅计算 m 的因子的 f，那么考虑知道 x 所有真因子的时候怎么计算 f_m

首先直接要求每个元素都是 m 的因子，有 2^d(m) 种方案

然后考虑容斥掉，就是这里面会出现 gcd=x, lcm=y 的集合，且 (x,y)!=(1,m)

那么这种集合的数量恰好等于 f_(y/x)，也就是把 f_m 减去所有 d(m/i)*f(i) 应该就行了

(x,y)的数量是满足y/x=i，y|m 的个数，即 x|(m/i)，也就是d(m/i)

题解2

这个题解写的略抽象，想了很久，没太懂对称是啥意思，后来才想明白

对于每个质因子p，我们要用全集，减掉gcd是pi的倍数的集合，也要减掉lcm不是pi^ai的集合

每个质因子有两个集合，1e18最多有15个不同质因子，最多有30个集合，暴力枚举是O(2^30)的

实际是因为，只考虑gcd是pi的倍数时，方案是 $2^{d(m/p)}-1$ ，

而只考虑lcm不是pi^ai时，说明是pi^ai删了至少一个质因子pi，方案数也是 $2^{d(m/p)}-1$

这俩是一样的集合，所以可以把这俩合并成一个集合，

用选了其中0个/选了其中1个/选了其中2个去做区分，

选了两个就除以p两次，就优化到O(3^15)了

代码1（O(d(n)^2，但常数略小)

#include<bits/stdc++.h>
#include<iostream>
#include<cstdio>
#include<vector>
#include<map>
#include<queue>
#include<set>
using namespace std;
#define rep(i,a,b) for(int i=(a);i<=(b);++i)
#define per(i,a,b) for(int i=(a);i>=(b);--i)
typedef long long ll;
typedef double db;
typedef pair<int,int> P;
#define fi first
#define se second
#define pb push_back    
#define dbg(x) cerr<<(#x)<<":"<<x<<" ";
#define dbg2(x) cerr<<(#x)<<":"<<x<<endl;
#define SZ(a) (int)(a.size())
#define sci(a) scanf("%d",&(a))
#define pt(a) printf("%d",a);
#define pte(a) printf("%d\n",a)
#define ptlle(a) printf("%lld\n",a)
#define debug(...) fprintf(stderr, __VA_ARGS__)
using namespace std;
typedef long long ll;
typedef __int128 LL;
const int N=114514,mod=998244353;
int two[N],mp[N];
ll ksm(ll a,ll b,ll mod) {
    ll ans=1;
    for(;b;a=(LL)a*a%mod,b/=2) {
        if(b&1) ans=(LL)ans*a%mod;
    }
    return ans;
}
inline ll mul(ll x,ll y,ll Mod){
    ll r=x*y-Mod*(ll)(1.L/Mod*x*y);
    return r-Mod*(r>=Mod)+Mod*(r<0);
}
ll gcd(ll x,ll y){return !y?x:gcd(y,x%y);}
namespace MR{
    const ll pr[7]={2,325,9375,28178,450775,9780504,1795265022};
    inline bool isprime(ll x){
        if(x<3||x%2==0)return x==2;
        int r=0;
        ll k=x-1;
        for(;!(k&1);k>>=1,++r);
        for(int T=7,j;T--;){
            ll p=pr[T];
            if(p%x==0)continue;
            ll v=ksm(p,k,x);
            if(v==1)continue;
            for(j=1;j<=r;++j,v=mul(v,v,x))
                if(v==x-1)break;
            if(j>r)return 0;
        }
        return 1;
    }
}
namespace PR{
    vector<ll>ans;
    inline ll calc(ll n){
        if(!(n&1))return 2;
        ll x=0,y=0,z=0,p=1,q,g;
        for(int i=0;(i&255)||(g=gcd(p,n))==1;\
            ++i,x=mul(x,x,n)+1,y=mul(y,y,n)+1,y=mul(y,y,n)+1){
            if(x==y)x=z++,y=mul(x,x,n)+1;
            q=mul(p,x-y+n,n);
            if(q)p=q;
        }
        return g;
    }
    inline void split(ll n){
        if(n==1)return;
        if(MR::isprime(n)){
            ans.push_back(n);
            return;
        }
        ll d=calc(n);
        while(d==1||d==n)d=calc(n);
        split(d),split(n/d);
    }
}
using namespace MR;
using namespace PR;
namespace Factor{
    typedef long long ll;
    static const int M=20;// https://oeis.org/A066150
    int num[M],d[M],c;
    ll p[M],v[N];
    vector<int>F[N];
    map<ll,int>fac;
    void init(){
        c=0;
        fac.clear();
    }
    // 质因数分解，视情况可以换为pollard_pho
    void cal_f(ll x){
        for(auto &d:ans){
            while(x%d==0)x/=d,fac[d]++;
        }
        if(x>1)fac[x]++;
    }
    // x=prod{v in x}，获取x的每个质因子值p，出现个数num
    void get_p(ll x){
        cal_f(x);
        for(auto &w:fac){
            p[c]=w.first;
            num[c++]=w.second;
        }
    }
    // x=prod{v in x}, 获取x的每个因子值v，因子个数前缀积d
    void merge(vector<int>&x,int p,vector<int>&z){
        int i=0,j=0,sz=x.size();
        while(i<sz && j<sz){
            int pos=x[j]+p;
            if(v[x[i]]<v[pos])z.pb(x[i++]);
            else if(v[x[i]]>v[pos])z.pb(pos),j++;
            else z.pb(x[i++]),j++;
        }
        while(i<sz){
            if(z.back()!=x[i])z.pb(x[i]);
            i++;
        }
        while(j<sz){
            int pos=x[j]+p;
            if(z.back()!=pos)z.pb(pos);
            j++;
        }
    }
    void get_d(ll x){
        //printf("x:%d\n",x);
        get_p(x);
        v[0]=d[0]=1;
        F[0].pb(0);
        for(int i=0;i<c;++i){
            d[i+1]=d[i]*(num[i]+1);
            for(int j=d[i];j<d[i+1];++j){
                v[j]=v[j-d[i]]*p[i];
                merge(F[j-d[i]],d[i],F[j]);
                //printf("j:%d v:%lld sz:%d\n",j,v[j],SZ(F[j]));
            }
        }
    }
};
using namespace Factor;
ll m;
int sol(int p){
    if(~mp[p])return mp[p];
    int sz=SZ(F[p]);
    int res=(two[sz]+mod-1)%mod;
    int p2=sz-1;
    //printf("p:%d v:%lld fac:\n",p,v[p]);
    //for(auto &v:F[p])printf("%d ",v);puts("");
    for(int i=0;i<sz;++i){
        int p1=F[p][i],sz1=SZ(F[p1]);
        while(p2>=0 && v[F[p][p2]]!=v[p]/v[p1])p2--;
        if(p1==0)continue;
        //printf("p:%d p1:%d p2:%d\n",p,p1,p2);
        res=(res+mod-1ll*sz1*sol(F[p][p2])%mod)%mod;
    }
    //printf("p:%d res:%d\n",p,res);
    return mp[p]=res;
}
int main(){
    two[0]=1;
    rep(i,1,N-5)two[i]=2ll*two[i-1]%mod;
    srand(time(0));
    scanf("%lld",&m);
    split(m);
    sort(ans.begin(),ans.end());
    ans.erase(unique(ans.begin(),ans.end()),ans.end());
    get_d(m);
    memset(mp,-1,sizeof mp);
    printf("%d\n",sol(d[c]-1));
    return 0;
}

代码2（3^n容斥）

#include <bits/stdc++.h>
#define ff first
#define ss second
using namespace std;

typedef long long ll;
typedef pair<int, int> pii;
typedef pair<ll, ll> pll;

typedef unsigned long long ull;
namespace pollard_rho {
	ull modmul(ull a, ull b, ull M) {
		ll ret = a * b - M * ull(1.L / M * a * b);
		return ret + M * (ret < 0) - M * (ret >= (ll)M);
	}
	ull modpow(ull b, ull e, ull mod) {
		ull ans = 1;
		for (; e; b = modmul(b, b, mod), e /= 2)
			if (e & 1) ans = modmul(ans, b, mod);
		return ans;
	}
	bool isPrime(ull n) {
		if (n < 2 || n % 6 % 4 != 1) return (n | 1) == 3;
		ull A[] = {2, 325, 9375, 28178, 450775, 9780504, 1795265022},
			s = __builtin_ctzll(n-1), d = n >> s;
		for (ull a : A) {   // ^ count trailing zeroes
			ull p = modpow(a%n, d, n), i = s;
			while (p != 1 && p != n - 1 && a % n && i--)
				p = modmul(p, p, n);
			if (p != n-1 && i != s) return 0;
		}
		return 1;
	}
	ull pollard(ull n) {
		auto f = [n](ull x) { return modmul(x, x, n) + 1; };
		ull x = 0, y = 0, t = 30, prd = 2, i = 1, q;
		while (t++ % 40 || __gcd(prd, n) == 1) {
			if (x == y) x = ++i, y = f(x);
			if ((q = modmul(prd, max(x,y) - min(x,y), n))) prd = q;
			x = f(x), y = f(f(y));
		}
		return __gcd(prd, n);
	}
	vector<ull> factor(ull n) {
		if (n == 1) return {};
		if (isPrime(n)) return {n};
		ull x = pollard(n);
		auto l = factor(x), r = factor(n / x);
		l.insert(l.end(), r.begin(), r.end());
		return l;
	}
}

const int MOD = 998244353;

vector<pair<ull, int>> P;

ll mpow(ll a, ll x) {
	ll r = 1;
	while (x) {
		if (x & 1) r = r * a % MOD;
		a = a * a % MOD;
		x /= 2;
	}
	return r;
}

int ans = 0;
void dfs(int i, int m, int c) {
	if (i == P.size()) {
		ans = (ans + (mpow(2, m) - 1) * c) % MOD;
		return;
	}
	int k = min(P[i].ss, 2);
	for (int s = 0; s <= k; s++) {
		dfs(i + 1, 1ll * m * (P[i].ss - s + 1) % (MOD - 1), 1ll * c * ((s == 1) ? MOD - 2 : 1) % MOD);
	}
}
int main() {
	cin.tie(0); ios_base::sync_with_stdio(0);
	ll n;
	cin >> n;

	vector<ull> V = pollard_rho::factor(n);
	map<ull, int> mp;
	for (ull x : V) mp[x]++;

	for (auto [p, c] : mp) P.push_back({p, c});

	dfs(0, 1, 1);

	cout << ans << '\n';
}

Code92007

关注

5
点赞
踩
5

收藏

觉得还不错? 一键收藏
打赏
0
评论
Petrozavodsk Programming Camp Winter 2018 Day 6: Grand Prix of Gomel K. Kids Aren‘t Alright（容斥好题）

那么这种集合的数量恰好等于 f_(y/x)，也就是把 f_m 减去所有 d(m/i)*f(i) 应该就行了。每个质因子有两个集合，1e18最多有15个不同质因子，最多有30个集合，暴力枚举是O(2^30)的。然后考虑容斥掉，就是这里面会出现 gcd=x, lcm=y 的集合，且 (x,y)!对于每个质因子p，我们要用全集，减掉gcd是pi的倍数的集合，也要减掉lcm不是pi^ai的集合。求非空集合的数量，满足gcd=1，lcm=m(m<=1e18)，答案取模998244353。
复制链接

扫一扫