【NOI2016】循环之美【莫比乌斯反演】【整除分块】【杜教筛】【类杜教筛】

传送门

题意：给定n,m,kn,m,kn,m,k,求1≤x≤n,1≤y≤m1\leq x\leq n,1\leq y\leq m1≤x≤n,1≤y≤m时xyx \over yyx中数值不同的纯循环小数或整数的个数。

n≤109,m≤109,k≤2×103n \leq 10^9,m\leq10^9,k\leq2\times10^3n≤109,m≤109,k≤2×103

显然只需要考虑最简分数即gcd⁡(x,y)=1\gcd(x,y)=1gcd(x,y)=1的情况

容（bu）易（yong）证明，xyx\over yyx满足题意当且仅当gcd⁡(y,k)=1\gcd(y,k)=1gcd(y,k)=1

所以

Ans=∑i=1N∑j=1M[gcd⁡(i,j)=1][gcd⁡(j,k)=1]Ans=\sum_{i=1}^N\sum_{j=1}^M[\gcd(i,j)=1][\gcd(j,k)=1]Ans=i=1∑Nj=1∑M[gcd(i,j)=1][gcd(j,k)=1]

换下顺序

=∑j=1M[gcd⁡(j,k)=1]∑i=1N[gcd⁡(i,j)=1]=\sum_{j=1}^M[\gcd(j,k)=1]\sum_{i=1}^N[\gcd(i,j)=1]=j=1∑M[gcd(j,k)=1]i=1∑N[gcd(i,j)=1]

把后面反演掉,前面不管

=∑j=1M[gcd⁡(j,k)=1]∑i=1N∑d∣i,d∣jμ(d)=\sum_{j=1}^M[\gcd(j,k)=1]\sum_{i=1}^N\sum_{d\mid i,d\mid j}\mu(d)=j=1∑M[gcd(j,k)=1]i=1∑Nd∣i,d∣j∑μ(d)

枚举ddd

=∑d=1N⌊Nd⌋μ(d)∑d∣jj≤M[gcd⁡(j,k)=1]=\sum_{d=1}^N\lfloor \frac{N}{d}\rfloor\mu(d)\sum_{d\mid j}^{j\leq M}[\gcd(j,k)=1]=d=1∑N⌊dN⌋μ(d)d∣j∑j≤M[gcd(j,k)=1]

换成枚举jjj是ddd的多少倍

=∑d=1N⌊Nd⌋μ(d)∑j=1⌊Md⌋[gcd⁡(jd,k)=1]=\sum_{d=1}^N\lfloor \frac{N}{d}\rfloor\mu(d)\sum_{j=1}^{\lfloor\frac{M}d{}\rfloor}[\gcd(jd,k)=1]=d=1∑N⌊dN⌋μ(d)j=1∑⌊dM⌋[gcd(jd,k)=1]

gcd⁡\gcdgcd拆开

=∑d=1N⌊Nd⌋μ(d)∑j=1⌊Md⌋[gcd⁡(j,k)=1][gcd⁡(d,k)=1]=\sum_{d=1}^N\lfloor \frac{N}{d}\rfloor\mu(d)\sum_{j=1}^{\lfloor\frac{M}d{}\rfloor}[\gcd(j,k)=1][\gcd(d,k)=1]=d=1∑N⌊dN⌋μ(d)j=1∑⌊dM⌋[gcd(j,k)=1][gcd(d,k)=1]

=∑d=1N⌊Nd⌋μ(d)[gcd⁡(d,k)=1]∑j=1⌊Md⌋[gcd⁡(j,k)=1]=\sum_{d=1}^N\lfloor \frac{N}{d}\rfloor\mu(d)[\gcd(d,k)=1]\sum_{j=1}^{\lfloor\frac{M}d{}\rfloor}[\gcd(j,k)=1]=d=1∑N⌊dN⌋μ(d)[gcd(d,k)=1]j=1∑⌊dM⌋[gcd(j,k)=1]

设

f(n,k)=∑i=1n[gcd⁡(i,k)=1]f(n,k)=\sum_{i=1}^n[\gcd(i,k)=1]f(n,k)=i=1∑n[gcd(i,k)=1]

由于kkk只有200020002000,并且gcd⁡(i,k)=gcd⁡(i%k,k)\gcd(i,k)=\gcd(i\%k,k)gcd(i,k)=gcd(i%k,k)，所以瞎预处理一下就可以算出来

这样

Ans=∑d=1N⌊Nd⌋μ(d)[gcd⁡(d,k)=1]f(⌊Md⌋,k)Ans=\sum_{d=1}^N\lfloor \frac{N}{d}\rfloor\mu(d)[\gcd(d,k)=1]f(\lfloor\frac{M}{d}\rfloor,k)Ans=d=1∑N⌊dN⌋μ(d)[gcd(d,k)=1]f(⌊dM⌋,k)

这是个整除分块的形式，我们只需要想办法算出

g(n,k)=∑i=1nμ(i)[gcd⁡(i,k)=1]g(n,k)=\sum_{i=1}^n\mu(i)[\gcd(i,k)=1]g(n,k)=i=1∑nμ(i)[gcd(i,k)=1]

这个可以反演

g(n,k)=∑i=1nμ(i)∑d∣i,d∣kμ(d)g(n,k)=\sum_{i=1}^n\mu(i)\sum_{d\mid i,d\mid k}\mu(d)g(n,k)=i=1∑nμ(i)d∣i,d∣k∑μ(d)

枚举ddd

g(n,k)=∑d∣kμ(d)∑d∣ii≤nμ(i)g(n,k)=\sum_{d\mid k}\mu(d)\sum_{d\mid i}^{i\leq n}\mu(i)g(n,k)=d∣k∑μ(d)d∣i∑i≤nμ(i)

枚举倍数

g(n,k)=∑d∣kμ(d)∑i=1⌊nd⌋μ(id)g(n,k)=\sum_{d\mid k}\mu(d)\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\mu(id)g(n,k)=d∣k∑μ(d)i=1∑⌊dn⌋μ(id)

前方高能

观察μ(id)\mu(id)μ(id),发现如果gcd⁡(i,d)≠1\gcd(i,d)\neq1gcd(i,d)=1,μ(id)\mu(id)μ(id)一定等于000

所以……我们可以丢一个[gcd⁡(i,d)=1][\gcd(i,d)=1][gcd(i,d)=1]进去

g(n,k)=∑d∣kμ(d)∑i=1⌊nd⌋μ(id)[gcd⁡(i,d)=1]g(n,k)=\sum_{d\mid k}\mu(d)\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\mu(id)[\gcd(i,d)=1]g(n,k)=d∣k∑μ(d)i=1∑⌊dn⌋μ(id)[gcd(i,d)=1]

由于μ\muμ是积性函数，有这个限制就可以拆开了，然后把μ(d)\mu(d)μ(d)丢前面去

g(n,k)=∑d∣kμ2(d)∑i=1⌊nd⌋μ(i)[gcd⁡(i,d)=1]g(n,k)=\sum_{d\mid k}\mu^2(d)\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\mu(i)[\gcd(i,d)=1]g(n,k)=d∣k∑μ2(d)i=1∑⌊dn⌋μ(i)[gcd(i,d)=1]

后面那个有没有很眼熟？

g(n,k)=∑d∣kμ2(d)g(⌊nd⌋,d)g(n,k)=\sum_{d\mid k}\mu^2(d)g(\lfloor\frac{n}{d}\rfloor,d)g(n,k)=d∣k∑μ2(d)g(⌊dn⌋,d)

然后就可以递归了

边界：

当n=0n=0n=0时，g(n,k)=0g(n,k)=0g(n,k)=0

当k=1k=1k=1时，我们发现无法递归

所以写个杜教筛算一下就可以了

复杂度O(能过)O(能过)O(能过)

#include <iostream>
#include <cstdio>
#include <cstring>
#include <cctype>
#include <utility>
#include <map>
using namespace std;
typedef long long ll;
const int N=1e6;
int n,m,k;
int np[N+5],pl[N],cnt;
int mu[N+5],sum[N+5];
void init()
{np[1]=mu[1]=1;for (int i=2;i<=N;i++){if (!np[i]) pl[++cnt]=i,mu[i]=-1;for (int j=1,x;(x=i*pl[j])<=N;j++){np[x]=1;if (i%pl[j]==0){mu[x]=0;break;}else mu[x]=-mu[i];}}for (int i=1;i<=N;i++) sum[i]=sum[i-1]+mu[i];
}
map<int,int> ms;
int getms(int n)
{if (n<=N) return sum[n];map<int,int>::iterator p;if ((p=ms.find(n))!=ms.end()) return p->second;int ans=1;for (int l=2,r;l<=n;l=r+1){r=n/(n/l);ans-=(r-l+1)*getms(n/l);}ms.insert(make_pair(n,ans));return ans;
}
int x[2005];
int gcd(int a,int b){return b? gcd(b,a%b):a;}
inline int f(const int& n){return n/k*x[k]+x[n%k];}
map<pair<int,int>,int> sg;
int g(int n,int k)
{if (n==0) return 0;if (k==1) return getms(n);map<pair<int,int>,int>::iterator p;pair<int,int> pi=make_pair(n,k);if ((p=sg.find(pi))!=sg.end()) return p->second;int ans=0;for (int i=1;i*i<=k;i++)if (k%i==0){if (mu[i]) ans+=g(n/i,i);if (i*i<k&&mu[k/i]) ans+=g(n/(k/i),k/i);}sg.insert(make_pair(pi,ans));return ans;
}
int main()
{init();scanf("%d%d%d",&n,&m,&k);for (int i=1;i<=k;i++) x[i]=x[i-1]+(gcd(i,k)==1);ll ans=0;for (int l=1,r;l<=n&&l<=m;l=r+1){r=min(n/(n/l),m/(m/l));ans+=(ll)(n/l)*f(m/l)*(g(r,k)-g(l-1,k));}printf("%lld\n",ans);return 0;
}

【NOI2016】循环之美【莫比乌斯反演】【整除分块】【杜教筛】【类杜教筛】相关推荐

luogu P3455 [POI2007]ZAP-Queries （莫比乌斯反演 + 整除分块）
整理的算法模板合集: ACM模板题目传送门本题中数据为5e4,我们只需要筛一次5e4就行了. 双倍经验的P4450 双亲数中数据达到了1e6,我们直接筛1e6的莫比乌斯函数有点不可取,因为只有一组 ...
欧拉函数+狄利克雷卷积+莫比乌斯函数+莫比乌斯反演+整除分块+杜教筛
Powered by:NEFU AB-IN 文章目录欧拉函数狄利克雷卷积莫比乌斯函数莫比乌斯反演 P3455 [POI2007]ZAP-Queries 整除分块 P2522 [HAOI2011 ...
BZOJ 1101 Luogu P3455 POI 2007 Zap (莫比乌斯反演+数论分块)
BZOJ 1101 Luogu P3455 POI 2007 Zap (莫比乌斯反演+数论分块) 手动博客搬家: 本文发表于20171216 13:34:20, 原地址https://blog.csd ...
luogu P1587 [NOI2016] 循环之美
https://www.luogu.com.cn/problem/P1587 首先思考我们要求的是什么? {xkly}={xy}\large \begin{Bmatrix}\frac{xk^l}{y} ...
P2522 HAOI2011 Problem b [莫比乌斯反演,数论分块]
P2522 HAOI2011 题意对于给出的n个询问,每次求有多少个数对(x,y)(x,y)(x,y),满足a≤x≤ba≤x≤ba≤x≤b,c≤y≤dc≤y≤dc≤y≤d,且gcd(x,y)=kgc ...
HDU 6833 莫比乌斯反演 + 数论分块
给定下列式子: ∑a1=1n∑a2=1n...∑ax=1n(∏j=1xajk)f(gcd⁡(a1,a2...,ax))⋅gcd⁡(a1,a2...,ax)\sum_{a_1=1}^n\sum_{a_2 ...
[NOI2016]循环之美(杜教筛)
首先要求每个数互不相等,故有$x\perp y$. 可以发现$\frac{x}{y}$在$k$进制下为纯循环小数的充要条件为$x\cdot k^{len}\equiv x(mod\ y)$,即$y\p ...
Noi2016 循环之美
题意相当于求如下式子 $$ \sum _{i=1}^{n} \sum_{j=1}^{m} [gcd(i,j)=1][gcd(j,k)=1] \\ k <= 1000 ; n,m <= 10 ...
数论 —— 莫比乌斯反演
[反演] 假设我们手头有个数列 F,通过某种变换 H,可以得到函数 G.,即: 但现在只有函数 G,需要求 F,那么我们就需要寻找一种变换 ,使得 G 在经过这种变换后能够获得 F,这个过程即为反演, ...
【专题】莫比乌斯反演
问题引入: 给定整数 N 和 M.求满足且为质数的点对的个数. 数据范围: 接下来会见到以下内容: 莫比乌斯函数莫比乌斯函数的线性筛迪利克雷卷积介绍莫比乌斯反演整除分块杜教筛介绍 ...

【NOI2016】循环之美【莫比乌斯反演】【整除分块】【杜教筛】【类杜教筛】

【NOI2016】循环之美【莫比乌斯反演】【整除分块】【杜教筛】【类杜教筛】相关推荐

最新文章

热门文章