题目描述
我们知道,从区间[L,H](L和H为整数)中选取N个整数,总共有(H-L+1)^N种方案。小z很好奇这样选出的数的最大公约数的规律,他决定对每种方案选出的N个整数都求一次最大公约数,以便进一步研究。然而他很快发现工作量太大了,于是向你寻求帮助。你的任务很简单,小z会告诉你一个整数K,你需要回答他最大公约数刚好为K的选取方案有多少个。由于方案数较大,你只需要输出其除以1000000007的余数即可。
输入输出格式
输入格式:
输入一行,包含4个空格分开的正整数,依次为N,K,L和H。
输出格式:
输出一个整数,为所求方案数。
输入输出样例
说明
样例解释
所有可能的选择方案:(2, 2), (2, 3), (2, 4), (3, 2), (3, 3), (3, 4), (4, 2), (4, 3), (4, 4)
其中最大公约数等于2的只有3组:(2, 2), (2, 4), (4, 2)
对于100%的数据,1<=N,K<=10^9,1<=L<=H<=10^9,H-L<=10^5
看到这题之后,我的第一感觉竟然是排列组合。。。。。
我想的是把k除了之后,只要是任意两个数的gcd为1的话,那整个的n个数的gcd
就为1,可想了好久发现这并不好处理,而且还不知道如何去重,还是看题解吧qwq
看到题解,才发现自己也太蒟了,但是做完这题之后感觉队反演有了更深的理解
反演的F(n) ,两个数的gcd的 d的倍数的个数是: n/d * m/d,
然后就可以推广到 n维, m/d ^ n
然后这题就是求L - R 之间的gcd为d的倍数的个数,那就是 R/d - L/d;
然后就是当L/d 为0 ,也就是d >L 的时候特判一下(一开始把这个给漏掉了。。。)
最后就是整出分块的时候需要mu的前缀和,1e9的数据量上杜教筛就行了
这题太容易爆long long 了,多mod几次没坏处
1 #include"bits/stdc++.h"
2 using namespace std;
3 typedef long long ll;
4 #define int ll
5
6 const int nn =2000100;
7 int maxn=2000000; const int mod = 1000000007;
8 int n,k,L,H;
9
10 int ptot;int prime[nn]; int vis[nn];
11
12 int mu[nn];int sum[nn];
13
14 unordered_map<int,int >mp;
15
16
17 int Mu()
18 { mu[1]=1;
19 for (int i=2;i<=maxn;i++)
20 {
21 if(!vis[i])
22 {
23 prime[++ptot]=i; mu[i]=-1;
24 }
25 for (int j=1;j<=ptot&&1ll*i*prime[j]<=1l*maxn;j++)
26 {
27 vis[i*prime[j]]=1;
28 if(i%prime[j]==0)
29 {
30 break;
31 }
32 else
33 {
34 mu[i*prime[j]]=-mu[i];
35 }
36 }
37
38 }
39 // for (int i=1;i<=20;i++)cout<<mu[i]<<" ";
40 for(int i=1;i<=maxn;i++)sum[i]=sum[i-1]+mu[i];
41
42
43 }
44
45 ll djs(int x)
46 {
47 if(x<=maxn)return sum[x];
48 if(mp[x])return mp[x];
49 int ans=1;
50 for(int l=2,r;l<=x;l=r+1)
51 {
52 r=x/(x/l);
53 ans-=1ll*(r-l+1)*djs(x/l)%mod;
54 ans%=mod;
55
56 }
57 ans=(ans%mod+mod)%mod;
58 return mp[x]=ans;
59
60 }
61
62 ll poww(ll a,ll b)
63 {
64 ll res=1; a%=mod;
65 for (;b;b>>=1,a=a*a%mod)if(b&1)res=res*a%mod;
66 return res;
67 }
68 void work()
69 {
70 ll ans=0;
71 H/=k; --L/=k;
72 // cout<<"MU"<<sum[1]<<" "<<sum[2]<<endl;
73 for (int l=1,r;l<=H;l=r+1)
74 {
75
76
77 if(L/l==0)
78 r=H/(H/l);
79 else r=min(L/(L/l),H/(H/l));
80 // cout<<l<<" "<<r<<endl;
81 // cout<<H/l<<" "<<L/l<<endl;
82 // cout<<djs(r)-djs(l-1)<<" "<<poww((H/l-L/l),n)<<endl;
83 ans+=( (djs(r)-djs(l-1))*1ll * poww((H/l-L/l),n))%mod;
84 ans%=mod;
85 }
86 ans=(ans%mod+mod)%mod; //这里要mod
87
88 cout<<ans<<endl;
89 return ;
90
91 }
92
93 signed main()
94 {
95 Mu();
96 //cout<<djs(1e9);
97 cin>>n>>k>>L>>H;
98 work();
99 }