【学习笔记】杜教筛

2023-10-12 10:28
文章标签 学习 笔记 杜教

本文主要是介绍【学习笔记】杜教筛,希望对大家解决编程问题提供一定的参考价值,需要的开发者们随着小编来一起学习吧!

补充几个前置知识

莫比乌斯反演
φ ∗ 1 = i d \varphi*1=id φ1=id

杜教筛

杜教筛能够在 Θ ( n 2 3 ) \Theta(n^\frac{2}{3}) Θ(n32)求出一个积性函数的前缀和。

给出一个积性函数 f ( i ) f(i) f(i),求 ∑ i = 1 n f ( i ) \sum_{i=1}^{n}f(i) i=1nf(i)

构造一个前缀和容易求的积性函数 g g g,令 h = f ∗ g h=f*g h=fg(其中 ∗ * 是狄利克雷卷积),使得 h h h的前缀和也容易求,那么

∑ i = 1 n h ( i ) = ∑ i = 1 n ∑ d ∣ i g ( d ) f ( i d ) = ∑ d = 1 n g ( d ) ∑ i = 1 ⌊ n d ⌋ f ( i ) = ∑ d = 1 n g ( d ) S ( ⌊ n d ⌋ ) = g ( 1 ) S ( n ) + ∑ d = 2 n g ( d ) S ( ⌊ n d ⌋ ) \sum_{i=1}^{n}h(i) \\ =\sum_{i=1}^{n}\sum_{d|i}g(d)f(\frac{i}{d}) \\ =\sum_{d=1}^{n}g(d)\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}f(i) \\ =\sum_{d=1}^{n}g(d)S(\lfloor\frac{n}{d}\rfloor) \\ =g(1)S(n)+\sum_{d=2}^{n}g(d)S(\lfloor\frac{n}{d}\rfloor) i=1nh(i)=i=1ndig(d)f(di)=d=1ng(d)i=1dnf(i)=d=1ng(d)S(⌊dn⌋)=g(1)S(n)+d=2ng(d)S(⌊dn⌋)

所以

g ( 1 ) S ( n ) = ∑ i = 1 n h ( i ) − ∑ d = 2 n g ( d ) S ( ⌊ n d ⌋ ) g(1)S(n)=\sum_{i=1}^{n}h(i)-\sum_{d=2}^{n}g(d)S(\lfloor\frac{n}{d}\rfloor) g(1)S(n)=i=1nh(i)d=2ng(d)S(⌊dn⌋)
不难发现左边部分就是 h ( i ) h(i) h(i)的前缀和,右边可以通过整除分块做。
实现的时候需要先预处理前 Θ ( n 2 3 ) \Theta(n^\frac{2}{3}) Θ(n32)项,然后用map存下求出的已知项。时间复杂度 Θ ( n 2 3 ) \Theta(n^\frac{2}{3}) Θ(n32)(不会证qwq)

int S(int n)
{if(n<=5000000) return sum[n];if(mp[n]!=0) return mp[n];int ans=g(n);for(int l=2,r; l<=n; l=r+1){r=n/(n/l);(ans-=S(n/l)*(t(r)-t(l-1))%mod)%=mod;}return mp[n]=ans;
}

如何构造g?

我们的目的是让 h = f ∗ g h=f*g h=fg的前缀和也容易求,所以可以先把 ( f ∗ g ) ( n ) (f*g)(n) (fg)(n)写出来,看看 h ( n ) h(n) h(n)长啥样,努力凑出一个可以求前缀和的形式。

例题

【模板】杜教筛

题目描述

给定一个正整数,求

a n s 1 = ∑ i = 1 n φ ( i ) ans_1=\sum_{i=1}^n\varphi(i) ans1=i=1nφ(i)

a n s 2 = ∑ i = 1 n μ ( i ) ans_2=\sum_{i=1}^n \mu(i) ans2=i=1nμ(i)

输入格式

本题单测试点内有多组数据

输入的第一行为一个整数,表示数据组数 T T T

接下来 T T T 行,每行一个整数 n n n,表示一组询问。

输出格式

对于每组询问,输出一行两个整数,分别代表 a n s 1 ans_1 ans1 a n s 2 ans_2 ans2

样例 #1

样例输入 #1
6
1
2
8
13
30
2333
样例输出 #1
1 1
2 0
22 -2
58 -3
278 -3
1655470 2

提示

数据规模与约定

对于全部的测试点,保证 1 ≤ T ≤ 10 1 \leq T \leq 10 1T10 1 ≤ n < 2 31 1 \leq n \lt 2^{31} 1n<231

题解

莫比乌斯函数前缀和

μ ∗ 1 = ε μ∗1=ε μ1=ε
不难想到令 g = 1 g=1 g=1,那么 h = ε h=ε h=ε
所以式子就变成了
1 × S ( n ) = 1 − ∑ d = 2 n S ( ⌊ n d ⌋ ) 1\times S(n)=1-\sum_{d=2}^{n}S(\lfloor\frac{n}{d}\rfloor) 1×S(n)=1d=2nS(⌊dn⌋)

欧拉函数前缀和

方法一:根据定义用莫比乌斯反演做。
∑ i = 1 n φ ( i ) = ∑ i = 1 n ∑ j = i + 1 n [ g c d ( i , j ) = 1 ] \sum_{i=1}^{n}\varphi(i)=\sum_{i=1}^{n}\sum_{j=i+1}^{n}[gcd(i,j)=1] i=1nφ(i)=i=1nj=i+1n[gcd(i,j)=1]
这不就是莫比乌斯反演板子题吗

方法二:杜教筛
φ ∗ 1 = i d \varphi*1=id φ1=id
g = 1 g=1 g=1
S ( n ) = ∑ i = 1 n i − ∑ i = 2 n S ( ⌊ n d ⌋ ) S(n)=\sum_{i=1}^{n}i-\sum_{i=2}^nS(\lfloor\frac{n}{d}\rfloor) S(n)=i=1nii=2nS(⌊dn⌋)

代码

#include<bits/stdc++.h>
#define int long long
using namespace std;
const int N=1e6+7;
int mu[N],phi[N],s1[N],s2[N];
bool bz[N];
vector<int> p;
map<int,int> mp1,mp2;
void init(int n)
{mu[1]=1; phi[1]=1;for(int i=2; i<=n; i++){if(!bz[i]){p.push_back(i);mu[i]=-1; phi[i]=i-1;}for(int j:p){if(i*j>n) break;bz[i*j]=1;if(i%j==0){phi[i*j]=phi[i]*j;break;}mu[i*j]=-mu[i];phi[i*j]=phi[i]*(j-1);}}for(int i=1; i<=n; i++){s1[i]=mu[i]+s1[i-1];s2[i]=phi[i]+s2[i-1];}
}
int S1(int n)
{if(n<2) return n;if(n<=N-7) return s1[n];if(mp1[n]) return mp1[n];int ans=1;for(int l=2,r; l<=n; l=r+1){r=n/(n/l);ans-=S1(n/l)*(r-l+1);}return mp1[n]=ans;
}
int S2(int n)
{if(n<2) return n;if(n<=N-7) return s2[n];if(mp2[n]) return mp2[n];int ans=n*(n+1)/2;for(int l=2,r; l<=n; l=r+1){r=n/(n/l);ans-=S2(n/l)*(r-l+1);}return mp2[n]=ans;
}
void O_o()
{int n;cin>>n;cout<<S2(n)<<" "<<S1(n)<<"\n";
}
signed main()
{ios::sync_with_stdio(false);cin.tie(0); cout.tie(0);	int T=1;init(N-7);cin>>T;while(T--){O_o();}
}

简单的数学题

题目描述

由于出题人懒得写背景了,题目还是简单一点好。

输入一个整数 n n n 和一个整数 p p p,你需要求出:

( ∑ i = 1 n ∑ j = 1 n i j gcd ⁡ ( i , j ) ) m o d p \left(\sum_{i=1}^n\sum_{j=1}^n ij \gcd(i,j)\right) \bmod p (i=1nj=1nijgcd(i,j))modp

其中 gcd ⁡ ( a , b ) \gcd(a,b) gcd(a,b) 表示 a a a b b b 的最大公约数。

输入格式

一行两个整数 p , n p,n p,n

输出格式

一行一个整数表示答案。

样例 #1

样例输入 #1
998244353 2000
样例输出 #1
883968974

提示

对于20%的数据, n ≤ 1000 n \leq 1000 n1000

对于30%的数据, n ≤ 5000 n \leq 5000 n5000

对于60%的数据, n ≤ 1 0 6 n \leq 10^6 n106,时限1s。

对于另外20%的数据, n ≤ 1 0 9 n \leq 10^9 n109,时限3s。

对于最后20%的数据, n ≤ 1 0 10 n \leq 10^{10} n1010,时限4s。

对于100%的数据, 5 × 1 0 8 ≤ p ≤ 1.1 × 1 0 9 5 \times 10^8 \leq p \leq 1.1 \times 10^9 5×108p1.1×109 p p p 为质数。

题解

引入欧拉反演(用莫比乌斯反演能做,推出来是一样的,但是很麻烦)

∑ i = 1 n ∑ j = 1 n g c d ( i , j ) i j = ∑ i = 1 n ∑ j = 1 n i d ( g c d ( i , j ) ) i j = ∑ i = 1 n ∑ j = 1 n ∑ d ∣ g c d ( i , j ) φ ( d ) i j = ∑ d = 1 n φ ( d ) ∗ d 2 ∑ i = 1 ⌊ n d ⌋ ∑ j = 1 ⌊ n d ⌋ i j \sum_{i=1}^{n}\sum_{j=1}^{n}gcd(i,j)ij \\ =\sum_{i=1}^{n}\sum_{j=1}^{n}id(gcd(i,j))ij \\ =\sum_{i=1}^{n}\sum_{j=1}^{n}\sum_{d|gcd(i,j)}\varphi(d)ij \\ =\sum_{d=1}^{n}\varphi(d)*d^2\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{d}\rfloor}ij i=1nj=1ngcd(i,j)ij=i=1nj=1nid(gcd(i,j))ij=i=1nj=1ndgcd(i,j)φ(d)ij=d=1nφ(d)d2i=1dnj=1dnij

s u m ( n ) = ∑ i = 1 n ∑ j = 1 n i j sum(n)=\sum_{i=1}^{n}\sum_{j=1}^{n}ij sum(n)=i=1nj=1nij
显然 s u m ( n ) = n 2 ( n + 1 ) 2 4 sum(n)=\frac{n^2(n+1)^2}{4} sum(n)=4n2(n+1)2

f ( i ) = φ ( i ) ∗ i 2 f(i)=\varphi(i)*i^2 f(i)=φ(i)i2,现在我们如果能求出它的前缀和,我们就可以整除分块解决问题。这个时候就要杜教筛了:
( f ∗ g ) ( n ) = ∑ d ∣ n φ ( d ) ∗ d 2 ∗ g ( n d ) (f*g)(n)\\=\sum_{d|n}\varphi(d)*d^2*g(\frac{n}{d}) (fg)(n)=dnφ(d)d2g(dn)
d 2 d^2 d2特别烦人怎么办?消掉!
g ( i ) = i 2 g(i)=i^2 g(i)=i2
( f ∗ g ) ( n ) = ∑ d ∣ n φ ( d ) ∗ d 2 ∗ g ( n d ) = n 2 ∑ d ∣ n φ ( d ) = n 3 (f*g)(n)\\=\sum_{d|n}\varphi(d)*d^2*g(\frac{n}{d})\\=n^2\sum_{d|n}\varphi(d)\\=n^3 (fg)(n)=dnφ(d)d2g(dn)=n2dnφ(d)=n3

用我们小学二年级就学过的立方和公式得
∑ i = 1 n i 3 = n 2 ( n + 1 ) 2 4 = s u m ( n ) \sum_{i=1}^{n}i^3=\frac{n^2(n+1)^2}{4}=sum(n) i=1ni3=4n2(n+1)2=sum(n)
于是乎,整除分块即可

代码

#include<bits/stdc++.h>
#define int long long
using namespace std;
const int N=5e6+7;
int mod,inv6,phi[N],sum[N];
bool bz[N];
vector<int> p;
map<int,int> mp;
int power(int x,int t)
{int b=1;while(t){if(t&1) b=b*x%mod;x=x*x%mod; t>>=1;}return b;
}
void init(int n)
{phi[1]=1;for(int i=2; i<=n; i++){if(!bz[i]){p.push_back(i); phi[i]=i-1;}for(int j:p){if(i*j>n) break;bz[i*j]=1;if(i%j==0){phi[i*j]=phi[i]*j%mod;break;}phi[i*j]=phi[i]*phi[j]%mod;}}for(int i=1; i<=n; i++) sum[i]=(sum[i-1]+phi[i]*i%mod*i%mod)%mod;
}
int t(int n)
{n%=mod;return n*(n+1)%mod*(2*n%mod+1)%mod*inv6%mod;
}
int g(int n)
{n%=mod;return (n*(n+1)/2%mod)*(n*(n+1)/2%mod)%mod;
}
int S(int n)
{if(n<=5000000) return sum[n];if(mp[n]!=0) return mp[n];int ans=g(n);for(int l=2,r; l<=n; l=r+1){r=n/(n/l);(ans-=S(n/l)*(t(r)-t(l-1))%mod)%=mod;}return mp[n]=ans;
}
void O_o()
{int n,ans=0;cin>>mod>>n;inv6=power(6,mod-2);init(5000000);for(int l=1,r; l<=n; l=r+1){r=n/(n/l);(ans+=(S(r)-S(l-1))*g(n/l)%mod)%=mod;}cout<<(ans+mod)%mod;
}
signed main()
{ios::sync_with_stdio(false); cin.tie(0); cout.tie(0);cout<<fixed<<setprecision(2);int T=1;//	cin>>T;while(T--){O_o();}
}

这篇关于【学习笔记】杜教筛的文章就介绍到这儿,希望我们推荐的文章对编程师们有所帮助!



http://www.chinasem.cn/article/195203

相关文章

Java进阶学习之如何开启远程调式

《Java进阶学习之如何开启远程调式》Java开发中的远程调试是一项至关重要的技能,特别是在处理生产环境的问题或者协作开发时,:本文主要介绍Java进阶学习之如何开启远程调式的相关资料,需要的朋友... 目录概述Java远程调试的开启与底层原理开启Java远程调试底层原理JVM参数总结&nbsMbKKXJx

Java深度学习库DJL实现Python的NumPy方式

《Java深度学习库DJL实现Python的NumPy方式》本文介绍了DJL库的背景和基本功能,包括NDArray的创建、数学运算、数据获取和设置等,同时,还展示了如何使用NDArray进行数据预处理... 目录1 NDArray 的背景介绍1.1 架构2 JavaDJL使用2.1 安装DJL2.2 基本操

HarmonyOS学习(七)——UI(五)常用布局总结

自适应布局 1.1、线性布局(LinearLayout) 通过线性容器Row和Column实现线性布局。Column容器内的子组件按照垂直方向排列,Row组件中的子组件按照水平方向排列。 属性说明space通过space参数设置主轴上子组件的间距,达到各子组件在排列上的等间距效果alignItems设置子组件在交叉轴上的对齐方式,且在各类尺寸屏幕上表现一致,其中交叉轴为垂直时,取值为Vert

Ilya-AI分享的他在OpenAI学习到的15个提示工程技巧

Ilya(不是本人,claude AI)在社交媒体上分享了他在OpenAI学习到的15个Prompt撰写技巧。 以下是详细的内容: 提示精确化:在编写提示时,力求表达清晰准确。清楚地阐述任务需求和概念定义至关重要。例:不用"分析文本",而用"判断这段话的情感倾向:积极、消极还是中性"。 快速迭代:善于快速连续调整提示。熟练的提示工程师能够灵活地进行多轮优化。例:从"总结文章"到"用

【前端学习】AntV G6-08 深入图形与图形分组、自定义节点、节点动画(下)

【课程链接】 AntV G6:深入图形与图形分组、自定义节点、节点动画(下)_哔哩哔哩_bilibili 本章十吾老师讲解了一个复杂的自定义节点中,应该怎样去计算和绘制图形,如何给一个图形制作不间断的动画,以及在鼠标事件之后产生动画。(有点难,需要好好理解) <!DOCTYPE html><html><head><meta charset="UTF-8"><title>06

学习hash总结

2014/1/29/   最近刚开始学hash,名字很陌生,但是hash的思想却很熟悉,以前早就做过此类的题,但是不知道这就是hash思想而已,说白了hash就是一个映射,往往灵活利用数组的下标来实现算法,hash的作用:1、判重;2、统计次数;

零基础学习Redis(10) -- zset类型命令使用

zset是有序集合,内部除了存储元素外,还会存储一个score,存储在zset中的元素会按照score的大小升序排列,不同元素的score可以重复,score相同的元素会按照元素的字典序排列。 1. zset常用命令 1.1 zadd  zadd key [NX | XX] [GT | LT]   [CH] [INCR] score member [score member ...]

【机器学习】高斯过程的基本概念和应用领域以及在python中的实例

引言 高斯过程(Gaussian Process,简称GP)是一种概率模型,用于描述一组随机变量的联合概率分布,其中任何一个有限维度的子集都具有高斯分布 文章目录 引言一、高斯过程1.1 基本定义1.1.1 随机过程1.1.2 高斯分布 1.2 高斯过程的特性1.2.1 联合高斯性1.2.2 均值函数1.2.3 协方差函数(或核函数) 1.3 核函数1.4 高斯过程回归(Gauss

【学习笔记】 陈强-机器学习-Python-Ch15 人工神经网络(1)sklearn

系列文章目录 监督学习:参数方法 【学习笔记】 陈强-机器学习-Python-Ch4 线性回归 【学习笔记】 陈强-机器学习-Python-Ch5 逻辑回归 【课后题练习】 陈强-机器学习-Python-Ch5 逻辑回归(SAheart.csv) 【学习笔记】 陈强-机器学习-Python-Ch6 多项逻辑回归 【学习笔记 及 课后题练习】 陈强-机器学习-Python-Ch7 判别分析 【学

系统架构师考试学习笔记第三篇——架构设计高级知识(20)通信系统架构设计理论与实践

本章知识考点:         第20课时主要学习通信系统架构设计的理论和工作中的实践。根据新版考试大纲,本课时知识点会涉及案例分析题(25分),而在历年考试中,案例题对该部分内容的考查并不多,虽在综合知识选择题目中经常考查,但分值也不高。本课时内容侧重于对知识点的记忆和理解,按照以往的出题规律,通信系统架构设计基础知识点多来源于教材内的基础网络设备、网络架构和教材外最新时事热点技术。本课时知识