博客
关于我
强烈建议你试试无所不能的chatGPT,快点击我
SPOJ PGCD (mobius反演 + 分块)
阅读量:7128 次
发布时间:2019-06-28

本文共 3666 字,大约阅读时间需要 12 分钟。

转载请注明出处,谢谢    by---cxlove

题意 :求满足gcd(i , j)是素数(1 <= i <= n && 1 <= j <= m)二元组(i , j)个数。

很值得总结的题。。。

首先得会一点前提东西 。。。先简单说下Mobius反演,就是偏序集上的容斥原理。

定义

F(n) = sigma (G(d))   d | n

那么G(n) = sigma (F(d) * miu (n / d))  d | n

还有另外一个表达形式

F(n) = sigma (G(d))  n | d

G(n) = sigma (F(d) * miu (d / n))  n | d

miu(i)是莫比乌斯函数

                1                (i = 1)

miu(i) =    (-1) ^ k       (i是square free number,k是i的素因子个数)

                0                (其它情况)

然后 就是要会求一个东西 GCD (i , j) = k 的数量。 (可以拿SPOJ   VLATTICE  先练练)

令F(k)表示GCD(i , j) >= k的数量。    这个东西 比较好求 就是(n / k) * (m / k)

然后令G(k)表示GCD(i , j) = k的数量。这个东西不 好求,就可以反演出来了

显然F(k) = sigma (G[d])  k | d

则G(k) = sigma (G[d] * miu (d / k))  k | d

这样就可以枚举d做出来了。。。。

对于这题,分成三个阶段。。。

stage 1 :打个素数表,然后预处理出miu函数。那么直接枚举素数p,求gcd(i , j) = p的数量,然后累加就行了。

预处理可以用线性筛选搞到O(n),每次查询是O(prime_count(n) * n)。显然会TLE。。。。

stage 2 :整理1中的式子,求的是sigma(F(d) * miu(d / p))  p为素数。

那么其中的F(d)只和d有关。。。可以枚举d,那么就要求出P[d] = sigma(miu(d / p))了p是d的某一个素因子。。。

如果d = p1 ^ k1 * p2 ^ k2 * …… pm ^ km

1、如果ki全为1,那么d / pi是square free number,而且素因子个数为m - 1。答案就是((m - 1) & 1 ? -1 : 1) * m

2、如果ki中只有1个为2,那么如果pi的指数为1,那么d / pi肯定不是square free number,那么miu值为0,不用管。如果pi的指数为2的话,那么d / pi肯定是square free number,素因子数个还为m,答案就是(m & 1 ? -1 : 1)

3、其余情况d / pi肯定不是square free number,miu值全为0,答案就是0.

处理完之后,对于查询只需要枚举d,然后 sigma(F(d) * P(d))

预处理同样可以用线性筛选完成O(n),每次查询O(n)。但是还是TLE。。。好SXBK。

PS:线性筛选太强大了。。。记录了每一个合数的最小素因子,保证每个合数最多被筛一次,保证了线性复杂度。

而且有了最小素因子minumfac[],就可以记录素因子个数totnum[]以及不同的素因子个数diffnum[]。就可以处理出P[]。

stage 3:又是经典的分块,考虑F(d) = (n / d) * (m / d),这个东西只有log()级别种,相邻的部分值是一样的。

对于这种分块不熟的可以做CQOI2007 余数之和sum

这样的话预处理O(n),查询O(sqrt(n))。。。终于过了

PS:SPOJ把MLE判成RE,然后 我就RE了10多次。。。。。线性筛选的时候少了一句话,又TLE好多次。。。

 

#include 
#include
#include
#include
#include
#include
#define pii pair
using namespace std;typedef long long LL;const int N = 10000005;int cnt = 0 , prime[N] , minumfac[N];int totnum[N] , diffnum[N] , p[N];bool isnotprime[N];// int miu[N];void init () { isnotprime[0] = isnotprime[1] = true; cnt = 0 ; // miu[1] = 1; for (int i = 2 ; i < N ; i ++) { if (!isnotprime[i]) { prime[cnt ++] = i; // miu[i] = -1; } for (int j = 0 ; j < cnt && i * prime[j] < N ; j ++) { isnotprime[i * prime[j]] = true; minumfac[i * prime[j]] = prime[j]; // if (i % prime[j] == 0) miu[i * prime[j]] = 0; // else miu[i * prime[j]] = - miu[i]; if (i % prime[j] == 0) break; } } p[1] = 0; for (int i = 2 ; i < N ; i ++) { if (!isnotprime[i]) totnum[i] = diffnum[i] = 1; else { totnum[i] = totnum[i / minumfac[i]] + 1; diffnum[i] = diffnum[i / minumfac[i]] + ((i / minumfac[i]) % minumfac[i] != 0); } if (totnum[i] == diffnum[i]) p[i] = (((totnum[i] - 1) & 1) ? -1 : 1) * totnum[i]; else if (totnum[i] - diffnum[i] == 1) p[i] = (diffnum[i] & 1) ? -1 : 1; else p[i] = 0; p[i] += p[i - 1]; }}int main () { #ifndef ONLINE_JUDGE freopen ("input.txt" , "r" , stdin); // freopen ("output.txt" , "w" , stdout); #endif int t; scanf ("%d" , &t); init (); while (t --) { int n , m ; LL ans = 0; scanf ("%d %d" , &n , &m); if (n > m) swap (n , m); /* // stage 1: // make talbe : O(n) // for every query : O(n * prime_count(n)) for (int i = 0 ; i < cnt ; i ++) { int p = prime[i]; for (int j = p ; j <= n ; j += p) { #define f(a , b , p) ((a / p) * (b / p)) ans += miu[j / p] * f(n , m , j); } } */ /* // stage 2: // make table : O(n) // for every query : O(n) for (int i = 1 ; i <= n ; i ++) { #define f(d) ((n / d) * (m / d)) ans += p[i] * f(i); } */ // stage 3: // make table : O(n) // for every query : O(sqrt (n)) for (int i = 2 , next; i <= n ; i = next) { int d1 = n / i , d2 = m / i; int next_n = n / d1 + 1 , next_m = m / d2 + 1; next = min (next_m , next_n); ans += 1LL * (p[next - 1] - p[i - 1]) * d1 * d2; } printf ("%lld\n" , ans); } return 0;}

 

 

你可能感兴趣的文章
Kali Linux 2019.1 发布,Metasploit 更新到 5.0 版本
查看>>
【mysql的设计与优化专题(1)】ER图,数据建模与数据字典
查看>>
Jibo’s Name: How did we pick it?
查看>>
device's media capture mechanism,利用input:file调用设备的照相机/相册、摄像机、录音机...
查看>>
BroadLink:三款新品力求无障碍人机交互,三大平台分三期对外开放 ...
查看>>
掌门1对1获3.5亿美元E-1轮融资,华人文化产业基金、中金甲子基金等投资 ...
查看>>
Unity中的通用对象池
查看>>
ORA-00600: internal error code, arguments: [16703], [1403], [28], [...
查看>>
忆芯科技发布新一代国产主控芯片STAR1000P!4月完成量产版本 ...
查看>>
如何用条码标签打印软件实现商品价签制定会员价 ...
查看>>
如何轻松实现个性化推荐系统
查看>>
Mysql高级查询 内连接和外连接详解
查看>>
基于AWS的电子商务网站架构——Web前端
查看>>
基于险企传统资源优势的“一核三环”规划——互联网平台建设
查看>>
社交网络:有意义的不仅是邓巴数
查看>>
MySQL优化案例
查看>>
02 贝叶斯算法 - 案例一 - 鸢尾花数据分类
查看>>
场景数据互为表里!畅想2027,保险行业发展愿景
查看>>
hibernate4整合spring3出现java.lang.NoClassDefFoundError: [Lorg/hibernate/engine/FilterDefinition;...
查看>>
港科大教授权龙:三维视觉重新定义人工智能安防
查看>>