P3768 简单的数学题 [狄利克雷卷积,杜教筛,莫比乌斯反演]

简单的数学题

题目连接

https://www.luogu.org/problemnew/show/P3768

题目描述

输入一个正整数n,n≤1010n,n\le 10^{10}n,n1010p,p≤1.1×109p,p \le 1.1 \times 10^9p,p1.1×109.且ppp为质数.

计算∑i=1n∑j=1nijgcd(i,j)\sum_{i=1}^n\sum_{j=1}^nijgcd(i,j)i=1nj=1nijgcd(i,j)ppp取模.

题解

方法一.暴力推公式

  1. 考虑把枚举gcd(i,j)gcd(i,j)gcd(i,j)
    原式=∑d=1nd3∑i=1n/d∑j=1n/dij×[gcd(i,j)=1]\sum_{d=1}^nd^3\sum_{i=1}^{n/d}\sum_{j=1}^{n/d}ij\times [gcd(i,j)=1]d=1nd3i=1n/dj=1n/dij×[gcd(i,j)=1]
  2. 莫比乌斯反演处理后边的式子
    f(x)=∑i=1n/d∑j=1n/dij×[gcd(i,j)=x]f(x)=\sum_{i=1}^{n/d}\sum_{j=1}^{n/d}ij\times [gcd(i,j)=x]f(x)=i=1n/dj=1n/dij×[gcd(i,j)=x]
    g(x)=∑x∣df(d)=x2∑i=1n/dx∑j=1n/dxijg(x) = \sum_{x|d}f(d)=x^2\sum_{i=1}^{n/dx}\sum_{j=1}^{n/dx}ijg(x)=xdf(d)=x2i=1n/dxj=1n/dxij
    f(x)=∑x∣pμ(p/x)g(p)=∑x∣pμ(p/x)p2∑i=1n/dp∑j=1n/dpijf(x)=\sum_{x|p}\mu(p/x)g(p)=\sum_{x|p}\mu(p/x)p^2\sum_{i=1}^{n/dp}\sum_{j=1}^{n/dp}ijf(x)=xpμ(p/x)g(p)=xpμ(p/x)p2i=1n/dpj=1n/dpij
    因此f(1)=∑p=1nμ(p)p2∑i=1n/dp∑j=1n/dpijf(1)=\sum_{p=1}^n\mu(p)p^2\sum_{i=1}^{n/dp}\sum_{j=1}^{n/dp}ijf(1)=p=1nμ(p)p2i=1n/dpj=1n/dpij
  3. 合并两部分
    原式=∑d=1nd3∑p=1nμ(p)p2∑i=1n/dp∑j=1n/dpij\sum_{d=1}^nd^3\sum_{p=1}^n\mu(p)p^2\sum_{i=1}^{n/dp}\sum_{j=1}^{n/dp}ijd=1nd3p=1nμ(p)p2i=1n/dpj=1n/dpij
    转而枚举a=dpa=dpa=dp.我们得到:
    原式=∑a=1n∑d∣aa2dμ(a/d)∑i=1n/ai∑j=1n/aj=\sum_{a=1}^n\sum_{d|a}a^2d\mu(a/d)\sum_{i=1}^{n/a}i\sum_{j=1}^{n/a}j=a=1ndaa2dμ(a/d)i=1n/aij=1n/aj
    =∑a=1n(⌊na⌋(⌊na⌋+1)2)2a2∑d∣adμ(a/d)=\sum_{a=1}^n (\frac{\lfloor \frac{n}{a} \rfloor(\lfloor \frac{n}{a} \rfloor+1)}{2})^2a^2\sum_{d|a}d\mu(a/d)=a=1n(2an(an+1))2a2dadμ(a/d)
    我们发现Id∗μ=ϕId*\mu=\phiIdμ=ϕ是常见的狄利克雷卷积.
    因此原式=∑a=1n(⌊na⌋(⌊na⌋+1)2)2a2ϕ(a)=\sum_{a=1}^n (\frac{\lfloor \frac{n}{a} \rfloor(\lfloor \frac{n}{a} \rfloor+1)}{2})^2a^2\phi(a)=a=1n(2an(an+1))2a2ϕ(a)
    化简到这一步的时候,前面部分可以分块计算,后面的部分a2ϕ(a)a^2\phi(a)a2ϕ(a)要快速计算前缀和.
    于是我们想到了杜教筛:
    f(x)=x2ϕ(x)f(x)=x^2\phi(x)f(x)=x2ϕ(x),找到积性函数g(x)=x2g(x)=x^2g(x)=x2与它做卷积使得g(x)g(x)g(x)的前缀和可以快速求出,而(f∗g)(x)=∑d∣xxd2×d2ϕ(d)=x2∑d∣xϕ(x)=x3(f*g)(x)=\sum_{d|x}\frac{x}{d}^2\times d^2\phi(d)=x^2\sum_{d|x}\phi(x)=x^3(fg)(x)=dxdx2×d2ϕ(d)=x2dxϕ(x)=x3的前缀和也可以快速求出.
    S(x)=∑i=1xf(x)S(x)=\sum_{i=1}^xf(x)S(x)=i=1xf(x)
    根据杜教筛公式g(1)S(n)=∑x=1n(f∗g)(x)−∑x=2ng(x)S(nx)g(1)S(n)=\sum_{x=1}^n(f*g)(x) - \sum_{x=2}^ng(x)S(\frac{n}{x})g(1)S(n)=x=1n(fg)(x)x=2ng(x)S(xn).
    S(n)=∑x=1nx3−∑x=2nx2S(nx)S(n)=\sum_{x=1}^nx^3-\sum_{x=2}^nx^2S(\frac{n}{x})S(n)=x=1nx3x=2nx2S(xn)
    写到这就完了.

方法二.ϕ\phiϕ卷积

根据公式:

gcd(i,j)=∑d∣gcd(i,j)ϕ(d)=∑d∣i,d∣jϕ(d)gcd(i,j) = \sum_{d|gcd(i,j)}\phi(d)=\sum_{d|i,d|j}\phi(d)gcd(i,j)=dgcd(i,j)ϕ(d)=di,djϕ(d)

因此

∑i=1n∑j=1nijgcd(i,j)=∑i=1n∑j=1nij∑d∣i,d∣jϕ(d)=∑d=1nϕ(d)d2∑i=1n/d∑j=1n/dij=∑d=1nϕ(d)d2(⌊nd⌋(⌊nd⌋+1)2)2\sum_{i=1}^n\sum_{j=1}^nijgcd(i,j)=\sum_{i=1}^n\sum_{j=1}^nij\sum_{d|i,d|j}\phi(d)=\sum_{d=1}^n\phi(d)d^2\sum_{i=1}^{n/d}\sum_{j=1}^{n/d}ij=\sum_{d=1}^n\phi(d)d^2(\frac{\lfloor \frac{n}{d} \rfloor(\lfloor \frac{n}{d} \rfloor+1)}{2})^2i=1nj=1nijgcd(i,j)=i=1nj=1nijdi,djϕ(d)=d=1nϕ(d)d2i=1n/dj=1n/dij=d=1nϕ(d)d2(2dn(dn+1))2

同样也得到了我们的式子,是不是推导方法简单多了?

遇见gcd(i,j)gcd(i,j)gcd(i,j)的时候,多想想是不是可以用ϕ\phiϕ卷积来做?

可能会减少很多不必要的推导过程.

代码

#include <iostream>
#include <algorithm>
#include <cstring>
#include <unordered_map>
#define pr(x) std::cout << #x << ':' << x << std::endl
#define rep(i,a,b) for(LL i = a;i <= b;++i)
const int N = 1e7;
typedef long long LL;
LL n,p; 
LL phi[N+10];
int prime[N+10],zhi[N+10],low[N+10],pcnt;
LL mod_pow(LL x,LL n) {LL res = 1;while(n) {if(n&1) res = res * x % p;x = x * x % p;n >>= 1;}return res;
}
LL inv6,inv2,inv4;
void sieve() {pcnt = 0;low[1] = phi[1] = zhi[1] = 1;rep(i,2,N) {if(!zhi[i]) {phi[i] = i-1;prime[pcnt++] = i;low[i] = i;}for(LL j = 0;j < pcnt && prime[j]*i <= N;++j) {zhi[i*prime[j]] = 1;if(i % prime[j] == 0) {low[i*prime[j]] = low[i] * prime[j];if(i == low[i]) {phi[i*prime[j]] = phi[i]*prime[j];}else {phi[i*prime[j]] = phi[i/low[i]] * phi[low[i]*prime[j]];}break;}else{low[i*prime[j]] = prime[j];phi[i*prime[j]] = phi[i] * phi[prime[j]];}}}
}
LL sum3(LL n) {n %= p;LL res = n*(n+1)%p*(2*n+1)%p*inv6%p;return res;
}
LL sum4(LL n) {n %= p;LL x = n*(n+1)%p;return x*x%p*inv4%p;
}
std::unordered_map<LL,LL> vis,rec;
LL F(LL n) {if(n <= N) return phi[n];if(vis[n]) return rec[n];LL res = sum4(n);for(LL x = 2,last;x <= n;x = last) {last = n/(n/x)+1;res = (res - ((sum3(last-1)-sum3(x-1)+p)%p*F(n/x)%p) + p) % p;}vis[n] = 1;return rec[n] = res;
}signed main() {sieve();std::ios::sync_with_stdio(false);std::cin >> p >> n;inv6 = mod_pow(6,p-2);inv4 = mod_pow(4,p-2);inv2 = mod_pow(2,p-2);rep(i,1,N) phi[i] = ((phi[i]*i%p*i%p) + phi[i-1]) % p;LL last,ans = 0;for(LL x = 1;x <= n;x = last+1) {last = n/(n/x);LL y = n/x%p;LL z = (1+y)*(y)/2%p;ans = (ans + (z*z%p*((F(last)-F(x-1)+p)%p))) % p;}std::cout << ans << std::endl;return 0;
}

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处:http://www.mzph.cn/news/320985.shtml

如若内容造成侵权/违法违规/事实不符,请联系多彩编程网进行投诉反馈email:809451989@qq.com,一经查实,立即删除!

相关文章

.NET Core微服务之基于Exceptionless实现分布式日志记录

一、Exceptionless极简介绍Exceptionless 是一个开源的实时的日志收集框架&#xff0c;它可以应用在基于 ASP.NET&#xff0c;ASP.NET Core&#xff0c;Web API&#xff0c;Web Forms&#xff0c;WPF&#xff0c;Console&#xff0c;ASP.NET MVC 等技术开发的应用程序中&#x…

P2801-教主的魔法【分块,二分】

正题 题目链接:https://www.luogu.com.cn/problem/P2801 题目大意 nnn个数字&#xff0c;要求支持 区间加上一个数字www询问一个区间内不小于www的数的个数 解题思路 考虑分块&#xff0c;对于块内我们维护一个排序后的数组&#xff0c;查询时直接在整个块中二分答案即可。修…

【多重背包】太空电梯(jzoj 1286)

太空电梯 题目大意&#xff1a; 有n&#xff08;1<n<400&#xff09;种石头&#xff0c;每种石头有它的数量c&#xff08;1<c<10&#xff09;&#xff0c;高度h&#xff08;1<h<100&#xff09;&#xff0c;可搭到的最高高度a&#xff08;1<a<40000&…

P2522 HAOI2011 Problem b [莫比乌斯反演,数论分块]

P2522 HAOI2011 题意 对于给出的n个询问&#xff0c;每次求有多少个数对(x,y)(x,y)(x,y)&#xff0c;满足a≤x≤ba≤x≤ba≤x≤b&#xff0c;c≤y≤dc≤y≤dc≤y≤d&#xff0c;且gcd(x,y)kgcd(x,y) kgcd(x,y)k&#xff0c;gcd(x,y)gcd(x,y)gcd(x,y)函数为xxx和yyy的最大公约…

.netcore 整合 log4net

1.背景前两天&#xff0c;曾经的一个同事咨询我&#xff0c;怎样将log4net以中间件的形式整合到core里边去。我不假思索的回答&#xff0c;这种问题应该有人做过吧&#xff0c;他说没有。于是&#xff0c;我去博客园搜了下&#xff0c;发现还真没有&#xff0c;全部都是传统.NE…

P5459-[BJOI2016]回转寿司【树状数组】

正题 题目链接:https://www.luogu.com.cn/problem/P5459 题目大意 nnn个数&#xff0c;求有多少个区间和在[L,R][L,R][L,R]范围内。 解题思路 显然我们做了前缀和之后&#xff0c;枚举右端点就只需要找到有多少个左端点满足在[x−R,x−L][x-R,x-L][x−R,x−L]这个范围内就好了…

【暴力】穹妹的求助

穹妹的求助 题目大意&#xff1a; 输入两个数&#xff0c;输出这两个数之间因数最多的数&#xff0c;和这个数的的因数个数 原题&#xff1a; 题目描述 由于穹妹很聪明&#xff0c;她的数学老师给她布置了一个作业&#xff0c;让她求出L到R之间不同因子数最多的那个数和这…

Asp.Net Core中利用Seq组件展示结构化日志功能

在一次.Net Core小项目的开发中&#xff0c;掌握的不够深入&#xff0c;对日志记录并没有好好利用&#xff0c;以至于一出现异常问题&#xff0c;都得跑动服务器上查看&#xff0c;那时一度怀疑自己肯定没学好&#xff0c;不然这一块日志不可能需要自己扒服务器日志来查看&…

P2633-Count on a tree【主席树,LCA】

正题 题目链接:https://www.luogu.com.cn/problem/P2633 题目大意 nnn个点的树&#xff0c;每个点有点权&#xff0c;求u∼vu\sim vu∼v的路径上第kkk小的权值&#xff0c;强制在线。 解题思路 考虑在树上维护主席树&#xff0c;我们不难发现如果sizxsiz_xsizx​表示根节点到…

P3327 约数的个数和 [约数函数性质,数论分块]

P3327 约数的个数和 题意 d(x)d(x)d(x)为约数的个数,对于每个询问,回答∑i1n∑j1md(ij)\sum_{i1}^n\sum_{j1}^md(ij)∑i1n​∑j1m​d(ij). 题解 这个题推得我头皮发麻,然后还没推出来,后来发现要做这题的先知道一个性质: d(ij)∑x∣i∑y∣j[gcd(x,y)1]d(ij)\sum_{x|i}\sum_{…

【结论】小X的矩阵

小X的矩阵 题目大意&#xff1a; 有一个nn的矩阵&#xff0c;要你执行g此操作&#xff0c;然后根据操作输出&#xff08;详情见原题&#xff09; 原题&#xff1a; 题目描述 小X最近迷上了矩阵&#xff0c;他定义了一个对于一种特殊矩阵的特征函数G。对于N∗NN*NN∗N的矩阵…

Apache SkyWalking的架构设计【译文】

Apache SkyWalking提供了一个功能强大并且很轻量级的后端。在此&#xff0c;将介绍为什么采用以下方式来设计它&#xff0c;以及它又是如何工作的。架构图对于APM而言&#xff0c;agent或SDKs仅是如何使用libs的技术细节。手动或自动的形式与架构无关&#xff0c;因此在本文中&…

P4735-最大异或和【可持久化Trie】

正题 题目链接:https://www.luogu.com.cn/problem/P4735 题目大意 nnn个数字&#xff0c;有操作 在末尾加入一个数字xxx询问[l,r][l,r][l,r]范围内的一个ppp使得ap⊕ap1⊕ap2...⊕an⊕xa_p\oplus a_{p1}\oplus a_{p2}...\oplus a_{n}\oplus xap​⊕ap1​⊕ap2​...⊕an​⊕x的…

【SPFA】桐人的约会

桐人的约会 题目大意&#xff1a; 删掉一条边&#xff0c;让一个图中的最短路最长 原题&#xff1a; 题目描述 这是一个风和日丽的日子&#xff0c;桐人和诗乃在约会。他们所在的城市共有N个街区&#xff0c;和M条道路&#xff0c;每条道路连接两个不同的街区&#xff0c;…

北京区域赛I题,Uva7676,A Boring Problem,前缀和差分

A Boring Problem 题解 其实这题不难,只要想到了前缀和差分就基本OK了. 我们要求的是第iii项的式子: F(i)(a1a2...ai)k(a2...ai)k...(ai)kF(i)(a_1a_2...a_i)^k(a2...a_i)^k...(a_i)^kF(i)(a1​a2​...ai​)k(a2...ai​)k...(ai​)k 记Sia1a2...ai,S00S_i a_1a_2...a_i,S_…

P4331-[BalticOI2004]Sequence数字序列【左偏树】

正题 题目链接:https://www.luogu.com.cn/problem/P4331 题目大意 给出一个序列aaa&#xff0c;求一个单调上升的序列bbb使得∑i1n∣ai−bi∣\sum_{i1}^n|a_i-b_i|∑i1n​∣ai​−bi​∣最小。 解题思路 巧妙的解法 首先我们让所有的ai−ia_i-iai​−i这样我们求的bbb序列就…

通俗易懂,什么是.NET?什么是.NET Framework?什么是.NET Core?

什么是.NET&#xff1f;什么是.NET Framework?本文将从上往下&#xff0c;循序渐进的介绍一系列相关.NET的概念&#xff0c;先从类型系统开始讲起&#xff0c;我将通过跨语言操作这个例子来逐渐引入一系列.NET的相关概念&#xff0c;这主要包括&#xff1a;CLS、CTS(CLI)、FCL…

【离散化】【DP】命运石之门的选择

命运石之门的选择 题目大意&#xff1a; 有n个盒子&#xff0c;高度为ai&#xff0c;可以数值刷盒子&#xff0c;也可以横着刷&#xff0c;但如果前面没盒子了&#xff0c;就要停下&#xff0c;问刷完这些盒子最少要刷多少次 原题&#xff1a; 题目描述 在某一条不知名世界…

P2153 晨跑,费用流裸题

晨跑 题目连接 https://www.luogu.org/problemnew/show/P2153 题解 求最大不相交路径数,并在路径数最大前提下,求总路程最短. 太裸了. 求不相交路径数:将除1,n1,n1,n两点外的所有点拆分,中间连一条容量为111,费用为000的边.然后所有的原边u→vu \rightarrow vu→v视作从u…

P3812-[模板]线性基

正题 题目链接:https://www.luogu.com.cn/problem/P3812 题目大意 给出nnn个数&#xff0c;求在其中选出若干个数使得它们的异或和最大。 解题思路 序列aaa的线性基bbb满足以下性质 aaa中的任何一个数都可以由bbb中的若干个数异或得到bbb中的任何一个数都不可由bbb中的若干个…