Min_25筛(EES筛)法
Extended Eratosthenes Sieve 参考链接
∀2≤i≤n\forall \ 2\le i\le n∀ 2≤i≤n,我们可以将iii分为两类
初始化S(n)=f(1)S(n)=f(1)S(n)=f(1)
另外,若kkk的最大质因子次幂>1>1>1,S(n)+=f(k)S(n)+=f(k)S(n)+=f(k),此时kkk就是一个第二类数。
- 如果我们dfsdfsdfs质因子来得到每个kkk,通常可以由积性函数的定义与性质简单地计算出f(k)f(k)f(k)。另外需要注意的是,如果L>nkL>\frac n kL>kn,则继续递归对答案的贡献为0,此时需要及时break,否则影响时间复杂度
- 如果我们可以O(1)O(1)O(1)地求出∑L<p≤nkf(p)\sum_{L<p\le \frac n k}f(p)∑L<p≤knf(p),那么上面过程的时间复杂度是O(满足kL<n的k的个数)O(满足kL< n的k的个数)O(满足kL<n的k的个数),当n≤1013n\le 10^{13}n≤1013时,时间复杂度为O(n34logn)O(\frac {n^{\frac 3 4}}{ log n})O(lognn43)
- 设g(i)=∑1≤p≤if(p)g(i)=\sum_{1\le p\le i}f(p)g(i)=∑1≤p≤if(p) ,现在问题只剩下了求∑L<p≤nkf(p)=g(⌊nk⌋)−g(L)\sum_{L<p\le \frac n k}f(p)=g(\lfloor \frac n k\rfloor)-g(L)∑L<p≤knf(p)=g(⌊kn⌋)−g(L) 。 由于⌊nk⌋\lfloor \frac n k\rfloor⌊kn⌋只有O(n)O(\sqrt{n})O(n)种,L≤nL\le \sqrt{n}L≤n也只有O(n)O(\sqrt{n})O(n)种,因此我们只需要计算ggg的O(n)O(\sqrt{n})O(n)项。
- 在题设里提到了f(p)f(p)f(p)是一个关于ppp的多项式,即f(p)=∑aipkipisprimef(p)=\sum a_ip^{k_i} \quad p \ is \ primef(p)=∑aipkip is prime ,我们对于每个iii,假设f(p)=pkif(p)=p^{k_i}f(p)=pki,最后乘上系数累加就可以得到ansansans
- 现在的问题是求g(i)=∑1≤p≤if(p)g(i)=\sum_{1\le p\le i}f(p)g(i)=∑1≤p≤if(p),注意这里ppp是素数,因此fff的非质数项的结果是不影响答案的,我们强行规定f(n)=nkif(n)=n^{k_i}f(n)=nki使得fff成为一个完全积性函数。
对于每个我们可能用到的g(i)g(i)g(i),我们只会在遍历不超过n\sqrt{n}n的质数时访问到,因此每个iii贡献的时间复杂度为O(ilogi)=O(ilogi)O(\frac {\sqrt{i} }{log\sqrt{i}})=O(\frac {\sqrt{i} }{logi})O(logii)=O(logii)
例题1 DIVCNTK
定义σ(n)=n的因子数\sigma(n)=n的因子数σ(n)=n的因子数 ,求∑i=1nσ(ik)mod264n,k1e10\sum_{i=1}^n\sigma(i^k) \ mod \ 2^{64} \quad \ n,k \ 1e10∑i=1nσ(ik) mod 264 n,k 1e10
//SPOJ DIVCNTK - Counting Divisors (general)
//Author : Feynman1999 9.27.2018
//f(1)=1
//f(p)=k+1
//f(p^e)=ek+1
#include<bits/stdc++.h>
using namespace std;
typedef unsigned long long u64;
u64 n,M,k;
//pre预处理后是2~i的p^0的和 p是素数
//hou是2~n/i的p^0的和
//同理,一个题目可能出现p^1 p^2等需要维护
vector<u64> pre,hou,primes;// 这里res是n/枚举的数
u64 dfs(u64 res, int last, u64 f){//最大质因子是prime[last-1] 但将1放在外面值显然一样u64 t=(res > M ? hou[n/res] : pre[res])-pre[primes[last]-1];u64 ret= t*f*(k+1);//这里需修改for(int i=last;i<(int) primes.size();++i){int p = primes[i];if((u64)p*p > res) break;for(u64 q=p,nres=res,nf=f*(k+1);q*p<=res;q*=p){//nf需修改ret += dfs (nres/=p,i+1,nf);//枚举更大的数nf += f*k;//继续枚举当前素数,指数大于1时,指数每加1,nf+=f*k ,k是系数ret += nf;//指数大于1时,记上贡献}}return ret;
}
u64 solve(u64 n){M=sqrt(n);pre.clear();pre.resize(M+1);hou.clear();hou.resize(M+1);primes.clear();primes.reserve(M+1);for(int i=1;i<=M;++i){pre[i]=i-1;hou[i]=n/i-1;}for(int p=2;p<=M;++p){if(pre[p]==pre[p-1]) continue;primes.push_back(p);const u64 q=(u64)p*p,m=n/p,pnt=pre[p-1];const int mid=M/p;const int End=min((u64)M,n/q);for(int i=1;i<=mid;++i) hou[i]-=hou[i*p]-pnt;for(int i=mid+1;i<=End;++i) hou[i]-=pre[m/i]-pnt;for(int i=M;i>=q;--i) pre[i]-=pre[i/p]-pnt;}primes.push_back(M+1);return n>1 ? 1+dfs(n,0,1) : 1;
}
int main()
{//freopen("in.txt","r",stdin);ios::sync_with_stdio(false);int t;cin>>t;while(t--){cin>>n>>k;cout<<solve(n)<<endl;}return 0;
}
例题2 神犇和蒟蒻
1<=N<=1E9
答案mod 1e9+7
显然A=1,我们只要计算B
//phi(1^2)=1
//phi(p^2)=p^2-p
//phi((p^e)^2)=phi((p^{e-1})^2)*p*p;
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int mod=1e9+7;
const int ni6=166666668;
const int ni2=500000004;
ll n,M;
vector<int> pre[3],hou[3],primes;inline int add(const int x, const int v) {return x + v >= mod ? x + v - mod : x + v;
}
inline int dec(const int x, const int v) {return x - v < 0 ? x - v + mod : x - v;
}//这里res是n/枚举的数
int dfs(ll res, int last, ll f){//最大质因子是prime[last-1] 但将1放在外面值显然一样int t=dec((res > M ? hou[2][n/res] : pre[2][res]),pre[2][primes[last]-1]);int ret= (ll)t*f%mod;//这里需修改for(int i=last;i<(int) primes.size();++i){int p = primes[i];if((ll)p*p > res) break;const int p2=(ll)p*p%mod;for(ll q=p,nres=res,nf=f*p%mod*(p-1)%mod;q*p<=res;q*=p){//nf需修改ret = add(ret,dfs(nres/=p,i+1,nf));//枚举更大的数nf = nf*p2%mod;//继续枚举当前素数,指数大于1时,指数每加1,nf=nf*p*p;ret =add(ret,nf);//指数大于1时,记上贡献}}return ret;
}inline int ff(ll x){x%=mod;return x*(x+1)%mod*ni2%mod;
}inline int fff(ll x){x%=mod;return x*(x+1)%mod*(2*x+1)%mod*ni6%mod;
}int solve(ll n){M=sqrt(n);for(int i=0;i<3;++i){pre[i].clear();pre[i].resize(M+1);hou[i].clear();hou[i].resize(M+1);}primes.clear();primes.reserve(M+1);for(int i=1;i<=M;++i){pre[0][i]=i-1;hou[0][i]=n/i-1;pre[1][i]=dec(ff(i),1);;hou[1][i]=dec(ff(n/i),1);pre[2][i]=dec(fff(i),1);hou[2][i]=dec(fff(n/i),1);}for(int p=2;p<=M;++p){if(pre[0][p]==pre[0][p-1]) continue;primes.push_back(p);const ll q=(ll)p*p;const int m=n/p,pnt0=pre[0][p-1],pnt1=pre[1][p-1],pnt2=pre[2][p-1];const int mid=M/p;const int End=min((ll)M,n/q);for(int i=1;i<=mid;++i){hou[0][i]=dec(hou[0][i],dec(hou[0][i*p],pnt0));hou[1][i]=dec(hou[1][i],dec(hou[1][i*p],pnt1)*(ll)p%mod);hou[2][i]=dec(hou[2][i],dec(hou[2][i*p],pnt2)*q%mod);}for(int i=mid+1;i<=End;++i){hou[0][i]=dec(hou[0][i],dec(pre[0][m/i],pnt0));hou[1][i]=dec(hou[1][i],dec(pre[1][m/i],pnt1)*(ll)p%mod);hou[2][i]=dec(hou[2][i],dec(pre[2][m/i],pnt2)*q%mod);}for(int i=M;i>=q;--i){pre[0][i]=dec(pre[0][i],dec(pre[0][i/p],pnt0));pre[1][i]=dec(pre[1][i],dec(pre[1][i/p],pnt1)*(ll)p%mod);pre[2][i]=dec(pre[2][i],dec(pre[2][i/p],pnt2)*q%mod);}}//cout<<clock()<<endl;primes.push_back(M+1);for (int i = 1; i <= M; i++) {pre[2][i] = dec(pre[2][i], pre[1][i]);//p^p-phou[2][i] = dec(hou[2][i], hou[1][i]);}return n>1 ? add(dfs(n,0,1),1) : 1;
}int main()
{//freopen("in.txt","r",stdin);ios::sync_with_stdio(false);cin>>n;cout<<1<<endl;cout<<solve(n)<<endl;return 0;
}
例题3 APS2 一道不是积性函数的题,说明了非积性函数的可行性
定义f(n)=n的最小质因子f(n)=n的最小质因子f(n)=n的最小质因子,求∑i=1nf(i)mod264\sum_{i=1}^nf(i) \ mod \ 2^{64}∑i=1nf(i) mod 264 ,1≤N≤12345678910111≤N≤12345678910111≤N≤1234567891011
//SPOJ DIVCNTK - Counting Divisors (general)
//f(1)=0
//f(p)=p
//f(p^e)=p
#include<bits/stdc++.h>
using namespace std;
typedef unsigned long long u64;
u64 n,M;
vector<u64> pre[2],hou[2],primes;
u64 ff(u64 A, u64 B){return (A+B)%2?(B-A+1)/2*(A+B):(A+B)/2*(B-A+1);
}
u64 dfs(u64 res, int last, u64 f){u64 ret,t;if(f>0){t = (res > M ? hou[0][n/res] : pre[0][res]) - pre[0][primes[last]-1];//有多少个素数//直接-id也行 因为pre[0][primes[last]-1]=idret = t * f * 1;//每个的权值这里就是f}else ret = hou[1][1];//第一次f=0 计算所有的素数的贡献for(int i=last;i<(int) primes.size();++i){int p = primes[i];if((u64)p*p > res) break;for(u64 q=p,nres=res,nf=(f==0?p:f);q*p<=res;q*=p){//f==0表明初始化ret += dfs (nres/=p,i+1,nf);ret += nf;//f==0 nf就是p 否则就是f}}return ret;
}
u64 solve(u64 n){M=sqrt(n);for(int i=0;i<2;++i){pre[i].clear();pre[i].resize(M+1);hou[i].clear();hou[i].resize(M+1);}primes.clear();primes.reserve(M+1);for(int i=1;i<=M;++i){pre[0][i]=i-1;pre[1][i]=ff(2,i);hou[0][i]=n/i-1;hou[1][i]=ff(2,n/i);}for(int p=2;p<=M;++p){if(pre[0][p]==pre[0][p-1]) continue;primes.push_back(p);const u64 q=(u64)p*p,m=n/p,pnt0=pre[0][p-1],pnt1=pre[1][p-1];const int mid=M/p;const int End=min((u64)M,n/q);for(int i=1;i<=mid;++i){hou[0][i]-=hou[0][i*p]-pnt0;hou[1][i]-=(hou[1][i*p]-pnt1)*p;}for(int i=mid+1;i<=End;++i){hou[0][i]-=pre[0][m/i]-pnt0;hou[1][i]-=(pre[1][m/i]-pnt1)*p;}for(int i=M;i>=q;--i){pre[0][i]-=pre[0][i/p]-pnt0;pre[1][i]-=(pre[1][i/p]-pnt1)*p;}}primes.push_back(M+1);return n>1 ? dfs(n,0,0) : 0;
}
int main()
{//freopen("in.txt","r",stdin);ios::sync_with_stdio(false);int t;cin>>t;while(t--){cin>>n;cout<<solve(n)<<endl;}return 0;
}
Min_25筛(EES筛)法相关推荐
- 经典中的经典之——筛选法求素数(埃氏筛 | 线性筛)
题目描述 统计小于非负整数n的质数数量 浑水摸鱼之蛮力验证法 直接上代码 bool is_zen(int x) {int i = 2;while (i * i <= x) {if (x % i ...
- [质数筛] 质数筛算法详解
今天给大家讲解质数筛这个算法. 在信息竞赛中,我们总是会遇到很多判断质数的题目,那么在这里就由我来给大家讲解一下质数筛算法(这里所有讲的算法都是基于筛出从 1 1 1 到 n n n 之间的素数的算法 ...
- 信息学奥赛一本通 2040:【例5.7】筛选法找质数 (普通筛 线性筛)
[题目链接] ybt 2040:[例5.7]筛选法找质数 [题目考点] 1. 普通筛法找质数(埃拉托色尼筛法) 基本思想:把从2到N的一组正整数从小到大按顺序排列.从中依次删除2的倍数.3的倍数.5的 ...
- 素数筛(筛选法求素数)
求素数 Problem Description 求小于n的所有素数的数量.(素数筛概念) Input 多组输入,输入整数n(n<1000000),以0结束. Output 输出n以内所有素数的个 ...
- 埃氏筛 线性筛(欧拉筛) 算法解析
埃氏晒 埃拉托斯特尼筛法,简称埃氏晒,是一种用来求自然数n以内的全部素数. 他的基本原理是,如果我们要获得小于n的所有素数,那就把不大于根号n的所有素数的倍数剔除. 埃氏晒的原理很容易理解,一个合数, ...
- 素数的线性筛法java,埃氏筛 线性筛(欧拉筛) 算法解析
埃氏晒 埃拉托斯特尼筛法,简称埃氏晒,是一种用来求自然数n以内的全部素数. 他的基本原理是,如果我们要获得小于n的所有素数,那就把不大于根号n的所有素数的倍数剔除. 埃氏晒的原理很容易理解,一个合数, ...
- 洲阁筛/Min25筛
1:完整的Min25筛学习方案 2:虽然不是但我也挂 3:网上捞了好久才捞到的纯递推写法Min25筛 终于看懂了..... 就是充分利用合数一定有小于等于n\sqrt nn的质因子和积性函数的性质来 ...
- 【BZOJ1607】轻拍牛头,筛一筛
1607: [Usaco2008 Dec]Patting Heads 轻拍牛头 Time Limit: 3 Sec Memory Limit: 64 MB Submit: 1404 Solved: 7 ...
- 素数筛线性筛详细详解(个人总结思路超长版)
一.埃氏筛 由于传统的sqrt写法比较简单,直接进行相除看是否能整除即可,这里不想过多的介绍此种方法.大多数比赛中这种写法也只能用于判断少量数据或无需大表即可通过的题目.这里从此种埃氏筛开始介绍. 此 ...
最新文章
- C++builder XE 安装控件 及输出路径
- 职业中专的计算机综合应用,职业中专计算机教学的思考
- 数据结构-深度优先遍历和广度优先遍历(漫画)
- [问题已处理]-[jenkins]-Jenkins 反向代理有误
- 在idea里如何实现Git项目回滚
- C眼看J - 初窥JAVA
- 漫步数学分析番外五(下)
- linux:C++的socket编程
- 第一个C#程序—C#基础回顾
- 如何提高Eclipse的运行速度 之总结
- iis出现HTTP 错误 403.14 - Forbidden Web问题
- 一步一步教你搭建外卖cps小程序
- requests-BeautifulSoup爬取美女贴吧图片
- 工业交换机的管理方式有哪些?
- 微信小程序《难忘便签》开发记录
- oracle自带的缓存,ORACLE缓存表与ORACLE缓存
- <_main__.类名 object at 0x0000000002A7CEB8>
- logging日志管理
- hiho一下 第139周 买零食 动态规划
- 天大青医堂第十期报告会之一