汉口网站建设制作,网站设计工,国内 wordpress,微信 公众号 微网站开发扯淡
min_25筛是由min_25提出的求积性函数前缀和的亚线性算法#xff0c;和一个叫“扩展埃氏筛”的东西有着微妙的关系。
至于是什么关系#xff0c;我也不太清楚#xff0c;反正有人说很像有人说就是一个东西#xff08;雾#xff09;
这段话并不是废话
约定
为了方…扯淡
min_25筛是由min_25提出的求积性函数前缀和的亚线性算法和一个叫“扩展埃氏筛”的东西有着微妙的关系。
至于是什么关系我也不太清楚反正有人说很像有人说就是一个东西雾
这段话并不是废话
约定
为了方便后面描述这里写一些用到的约定和符号表示以免产生恐惧
111被开除正整数籍 也就是说“前缀和”之类的都是从222开始对答案所求的前缀和同样最后手动加111
π(n)\pi(n)π(n)表示1∼n1 \sim n1∼n中质数个数的规模其实问题不大后面就懂了
pip_ipi表示第iii个质数单独的ppp均表示质数
primeprimeprime表示质数集合 minp(i)minp(i)minp(i)表示iii的最小质因子
pcp^cpc表示一个只含一个质因子的数
流程
首先所求的函数f(n)f(n)f(n)需要满足
是个积性函数在质数处的取值f(p)f(p)f(p)是一个关于ppp的多项式。我们把这个多项式拆成若干单项式分别计算在相加这样就变成了f(p)f(p)f(p)的值是一个关于ppp的单项式我们记为f(p)pkf(p)p^kf(p)pkf(pc)f(p^c)f(pc)的值可以快速计算
首先考虑计算这个东西
∑i2n[i∈prime]ik\sum_{i2}^n[i\in prime]i^ki2∑n[i∈prime]ik
注意是iki^kik不是f(i)f(i)f(i)虽然这里还没有区别但它们只是质数位置相等
首先既然叫扩展埃氏筛先回忆一下埃氏筛在干什么开始写下所有正整数然后从小到大如果一个没被筛说明它是质数用它筛掉后面的倍数
我们可以用这个思路计算上面的式子先算出所有的和再用质数筛掉所有合数项
设g(n,j)g(n,j)g(n,j)表示用前jjj个质数筛了之后剩余项的和
g(n,j)∑i2n[i∈primeorminp(i)pj]ikg(n,j)\sum_{i2}^n[i \in prime \quad or\quad minp(i)p_j]i^kg(n,j)i2∑n[i∈primeorminp(i)pj]ik
和就是g(n,π(n))g(n,\pi(\sqrt{n}))g(n,π(n))
注意合数的最小质因子不会超过n\sqrt nn,所以如果pj2np_j^2npj2n有g(n,j)g(n,j−1)g(n,j)g(n,j-1)g(n,j)g(n,j−1)
对剩下情况考虑转移即从g(n,j−1)g(n,j-1)g(n,j−1)去掉被pjp_jpj删掉的项
g(n,j−1)∑i2n[i∈primeorminp(i)≥pj]g(n,j-1)\sum_{i2}^n[i\in prime \quad or \quad minp(i)\geq p_j]g(n,j−1)i2∑n[i∈primeorminp(i)≥pj]
如果脑子转不过来可以强行算 [i∈primeorminp(i)≥pj][i\in prime \quad or \quad minp(i)\geq p_j][i∈primeorminp(i)≥pj]且不满足[i∈primeorminp(i)pj][i \in prime \quad or\quad minp(i)p_j][i∈primeorminp(i)pj] 即 [i∈primeorminp(i)≥pj][i\in prime \quad or \quad minp(i)\geq p_j][i∈primeorminp(i)≥pj]且[i∉prime]且[minp(i)≤pj][i \notin prime] 且 [minp(i)\leq p_j][i∈/prime]且[minp(i)≤pj] 综上 [i∉prime,minp(i)pj][i \notin prime,minp(i)p_j][i∈/prime,minp(i)pj] 所以
g(n,j)g(n,j−1)−∑i2n[i∉prime,minp(i)pj]ikg(n,j)g(n,j-1)-\sum_{i2}^n[i \notin prime,minp(i)p_j]i^kg(n,j)g(n,j−1)−i2∑n[i∈/prime,minp(i)pj]ik
可以提一个pjp_jpj出来,因为不能有质数就把pjp_jpj去掉刚好是111
g(n,j)g(n,j−1)−pjk∑i2⌊npj⌋[minp(i)≥pj]ikg(n,j)g(n,j-1)-p_j^k\sum_{i2}^{\lfloor\frac{n}{p_j}\rfloor}[minp(i)\geq p_j]i^kg(n,j)g(n,j−1)−pjki2∑⌊pjn⌋[minp(i)≥pj]ik
再次观察
g(n,j)∑i2n[i∈primeorminp(i)pj]ikg(n,j)\sum_{i2}^n[i \in prime \quad or\quad minp(i)p_j]i^kg(n,j)i2∑n[i∈primeorminp(i)pj]ik
发现可以拆成小于pjp_jpj的质数和大于等于pjp_jpj的所有数第二个就是刚才的式子
g(n,j)g(n,j−1)−pjk(g(⌊npj⌋,j−1)−∑i1j−1pik)g(n,j)g(n,j-1)-p_j^k(g(\lfloor\frac{n}{p_j}\rfloor,j-1)-\sum_{i1}^{j-1}p_i^k)g(n,j)g(n,j−1)−pjk(g(⌊pjn⌋,j−1)−i1∑j−1pik)
注意nnn由于都是一直整除所以只会有O(n)O(\sqrt n)O(n)种取值可以整除分块找出来强行离散化。在记录某个值vvv所在位置的时候如果vnv\sqrt nvn,我们另开一个数组存到⌊nv⌋\lfloor\frac{n}{v}\rfloor⌊vn⌋里面
然后后面只会用到最后一项所以第二维可以滚掉
求答案时设
S(n,j)∑i2n[minp(i)pj]f(i)S(n,j)\sum_{i2}^n[minp(i)p_j]f(i)S(n,j)i2∑n[minp(i)pj]f(i)
分质数和合数分别计算
质数部分用ggg算出来的减去前jjj个
g(n,π(n))−∑i1jpikg(n,\pi(\sqrt n))-\sum_{i1}^jp^k_ig(n,π(n))−i1∑jpik
合数部分枚举最小质因子和它的次数
∑kj1pk≤n∑e1pke≤nf(pke)(S(⌊npke⌋,k)[e1])\sum_{kj1}^{p_k\leq\sqrt n}\sum_{e1}^{p_k^e\leq n}f(p_k^e)(S(\lfloor\frac{n}{p_k^e}\rfloor,k)[e1])kj1∑pk≤ne1∑pke≤nf(pke)(S(⌊pken⌋,k)[e1])
[e1][e1][e1]指如果指数大于111它本身就是合数需要统计答案
两部分相加即可
总复杂度大概O(n23)O(n^{2\over3})O(n32)常数较小大约2s2s2s过1e101e101e10,3s3s3s过1e111e111e11
模板题
#include iostream
#include cstdio
#include cstring
#include cctype
#include cmath#define MAXN 200005
using namespace std;
typedef long long ll;
const int MOD1e97,INV6(MOD1)/6;
inline int add(const int x,const int y){return xyMOD? xy-MOD:xy;}
inline int dec(const int x,const int y){return xy? x-yMOD:x-y;}
int np[MAXN],pl[MAXN],cnt;
inline void init(const int N)
{np[1]1;for (int i2;iN;i){if (!np[i]) pl[cnt]i;int x;for (int j1;(xi*pl[j])N;j){np[x]1;if (i%pl[j]0) break;}}
}
ll val[MAXN];
int tot;
int g1[MAXN],g2[MAXN],sum1[MAXN],sum2[MAXN];
int key[MAXN],yek[MAXN];
ll n;
int m;
inline int getkey(const ll v){return vm? key[v]:yek[n/v];}
int S(ll n,int j)
{if (pl[j]n) return 0;int kgetkey(n);int ansdec(dec(g2[k],g1[k]),dec(sum2[j],sum1[j]));for (int kj1;kcnt(ll)pl[k]*pl[k]n;k)for (ll e1,pepl[k];pen;e,pe*pl[k])ansadd(ans,pe%MOD*(pe%MOD-1)%MOD*(S(n/pe,k)(e1))%MOD); return ans;
}
int main()
{scanf(%lld,n);msqrt(n);init(m);for (ll l1,r;ln;lr1){rn/(n/l);val[tot]n/l;int tval[tot]%MOD;g1[tot]((ll)t*(t1)/2)%MOD-1;g2[tot](ll)t*(t1)%MOD*(2*t1)%MOD*INV6%MOD-1;if (val[tot]m) key[val[tot]]tot;else yek[n/(n/l)]tot;}for (int j1;jcnt;j){for (int i1;itot(ll)pl[j]*pl[j]val[i];i){int kgetkey(val[i]/pl[j]);g1[i]dec(g1[i],(ll)pl[j]*dec(g1[k],sum1[j-1])%MOD);g2[i]dec(g2[i],(ll)pl[j]*pl[j]%MOD*dec(g2[k],sum2[j-1])%MOD);} sum1[j]add(sum1[j-1],pl[j]),sum2[j]add(sum2[j-1],(ll)pl[j]*pl[j]%MOD);}printf(%d\n,S(n,0)1);return 0;
}
