数论 莫比乌斯反演

前置需求

数论分块

概念

对于一个形如 \(\sum_{x=1}^n \lfloor{\frac{n}{x}}\rfloor\) 的式子，我们发现对于一部分的 \(x\)，它们的 \(\lfloor{\frac{n}{x}}\rfloor\) 值相同，因此我们没必要 \(\mathcal{O(n)}\) 计算，可以采用数论分块的办法将这一步的复杂度降低至 \(\mathcal{O(\sqrt{n})}\)。

对于一个已知值 \(n\)，给定一个值 \(l\)，求一个变量 \(r\) 的最大值使得：\(\lfloor{\frac{n}{l}}\rfloor=\lfloor{\frac{n}{r}}\rfloor\)。比较好得出结论，\(r=\lfloor{\frac{n}{\lfloor{\frac{n}{l}}\rfloor}}\rfloor\)。

证明

image

摘自 OI-wiki。

应用

在莫反中通常除 \(n\) 外还有一个 \(m\) 形成双层循环，对于它我们只需更改 \(r\) 的位置使得在 \(i\in \left[l,r\right]\) 中 \(\lfloor{\frac{n}{i}}\rfloor\) 与 \(\lfloor{\frac{m}{i}}\rfloor\) 的值始终不变即可。规定 \(n\le m\)。

if(n>m) swap(n,m);
int res=0;
for(int l=1,r;l<=n;l=r+1)
{
	r=min(n/(n/l),m/(m/l));
	res=/*code*/;
}

含 \(\prod\) 实现

例：求 $\prod_{x=1}^n f(x)^{\lfloor{\frac{n}{x}}\rfloor{}\ \lfloor{\frac{m}{x}}\rfloor{}} $。

思路是一样的，不过此时前缀和不再是相加而是相乘，需要用到快速幂和逆元。

由扩展欧拉定理可知：\(a^c\equiv a^{c\mod{\varphi{(p)}}}\mod{p}\quad gcd(a,p)=1\)，通常情况下模数均为质数，而欧拉函数有 $\varphi{(x)}=x-1\quad x\in prime $，所以对指数取模 \(mod-1\) 就好啦。

if(n>m) swap(n,m);
int res=1;
for(int l=1,r;l<=n;l=r+1)
{
    r=min(n/(n/l),m/(m/l));
    res=1ll*res*qpow(1ll*f[r]*qpow(f[l-1],mod-2)%mod,1ll*(n/l)*(m/l)%(mod-1))%mod;
}

线性筛

莫反中必用到的知识点之一。

首先了解莫比乌斯函数：

\[\mu(x)=\begin{cases} 1\qquad\quad n=1\\ 0 \qquad\quad n\,含有平方因子\\ (-1)^k\quad k\,为\, n\,本质不同的质因子个数\end{cases}\]

你需要知道的，一是它是一个积性函数，二是怎么线性筛出它。

从定义能看出，\(\mu\) 函数的值跟质数有关，所以在筛它的同时我们还需要线性筛出质数。

最简单的筛

\(u\) 为 \(\mu\) 函数的值，\(N\) 为值域上界，其他的都是筛素数常见的。

int u[N],pri[N],tot;
bool vis[N];
void Wprepare()
{
	u[1]=1;
	for(int i=2;i<=N;i++)
	{
		if(!vis[i]) pri[++tot]=i,u[i]=-1;
		for(int j=1;j<=tot;j++)
		{
			if(i*pri[j]>N) break;
			vis[i*pri[j]]=1;
			if(i%pri[j]==0) break;
			u[i*pri[j]]=-u[i];
		}
	}
}

莫反通常多测，大多数我们需要在准备工作中多预处理一些式子里的其他函数值，有时候还要对其做求前缀和等一系列操作，这些都是因题而异，需要我们灵活应对。

一定的推式子能力

几个显然的很常用的式子转化：

默认 \(n\le m\)。

基础一点的：

  • \[交换律\sum_x\sum\sum f(x)=\sum_x f(x)\sum\sum \]

  • \[同乘同除，变换上界\sum_{i=1}^n\sum_{j=1}^m [gcd(i,j)=d]=\sum_{i=1}^{\lfloor{\frac{n}{d}}\rfloor}\sum_{j=1}^{\lfloor{\frac{m}{d}}\rfloor}[gcd(i,j)=1] \]

  • \[等价描述，令\,T=dx：\sum_{d=1}^n\sum_{x=1}^{\lfloor{\frac{n}{d}}\rfloor}=\sum_{T=1}^n \sum_{x|T} \]

提高一点的：

  • \[综合上面\sum_{i=1}^n \sum_{j=1}^m f(gcd(i,j))=\sum_{d=1}^n f(d)\sum_{i=1}^{\lfloor{\frac{n}{d}}\rfloor} \sum_{j=1}^{\lfloor{\frac{m}{d}}\rfloor} \;[gcd(i,j)=1] \]

  • \[枚举因数，去循环\sum_{i=1}^n\sum_{j=1}^m\sum_{x|gcd(i,j)} \mu(x)=\sum_{x=1}^n \mu(x)\lfloor{\frac{n}{x}}\rfloor\lfloor{\frac{m}{x}}\rfloor \]

  • 枚举分母，方便数论分块
    令 \(T=dx\)

\[\begin{aligned} Ans &= \sum_{d=1}^n f(d) \sum_{x=1}^{\lfloor{\frac{n}{d}}\rfloor} \mu(x) \lfloor{\frac{n}{dx}}\rfloor \lfloor{\frac{m}{dx}}\rfloor \\ &= \sum_{T=1}^n \lfloor{\frac{n}{T}}\rfloor \lfloor{\frac{m}{T}}\rfloor \sum_{d|T} f(d) \mu(\frac{T}{d}) \end{aligned} \]

含 \(\prod\) 的，其实只有换位变了：

  • \[基本换位，可逆\prod_{x}\prod_{y} f(x)=\prod_{x} f(x)^{\sum_{y}} \]

具有结合律，并且适用于整除分块。

主体

由于蒟蒻不会卷积和其他的东西，所以只用到了最基础的一个反演公式

反演公式：

\[[gcd(i,j)=1]=\sum_{d|gcd(i,j)} \mu(d) \]

这是根据莫比乌斯函数的一个性质推出来的：

\[\sum_{d|n}=\begin{cases}1\quad n=1 \\ 0 \quad n\neq 1\end{cases} \]

该性质证明如下：

image

摘自 OI-wiki。

如果你和我一样对这类式子没有什么太强烈的求知欲，那么食用方法如下：

  • 根据题意，推出式子；
  • 尽可能向 \([gcd(i,j)=1]\) 的方向推，然后使用反演公式；
  • 结合上面的式子转化，推出一个可以预处理出函数值，并能运用数论分块优化复杂度的式子，然后 AC 开始调代码。

例题

做莫反的题我认为关键就在刷题，在不断地推式子中寻找共性，熟能生巧。

先推荐一篇学案，上面的推式子技巧都是循序渐进的，比较好。个人推荐一种学习方式：前几道较基础的题可以看着题解的式子一点一点理解，然后回过头来尝试自己独立完成推式子，再到不断优化筛法。

前几题式子部分会写的详细一点，可能导致篇幅冗长请见谅 qwq。

注：以下题未说明均默认 \(n\le m\)。

YY 的 GCD

这道题虽然在第一题的位置，但要考虑的其实还挺多的，综合了上面几乎所有的推式子技巧。

\[\begin{aligned} Ans&= \sum_{i=1}^n\sum_{j=1}^m \;[gcd(i,j)\in prime] \\ &= \sum_{d=1}^n\sum_{i=1}^n\sum_{j=1}^m \;[gcd(i,j)=d]\;(d\in prime)\\ &= \sum_{d=1}^n\sum_{i=1}^{\lfloor{\frac{n}{d}}\rfloor}\sum_{j=1}^{\lfloor{\frac{m}{d}}\rfloor} \;[gcd(i,j)=1]\;(d\in prime)\\ &=\sum_{d=1}^n\sum_{i=1}^{\lfloor{\frac{n}{d}}\rfloor}\sum_{j=1}^{\lfloor{\frac{m}{d}}\rfloor}\sum_{x|gcd(i,j)}\mu(x)\;(d\in prime)\\&= \sum_{d=1}^n\sum_{x=1}^{\lfloor{\frac{n}{d}}\rfloor}\;\mu(x) \lfloor{\frac{n}{dx}}\rfloor\lfloor{\frac{m}{dx}}\rfloor \;(d\in prime) \end{aligned} \]

令 \(T=dx\)，则有：

\[\begin{aligned}Ans&= \sum_{T=1}^n\sum_{d|T} \;\mu(\frac{T}{d})\lfloor{\frac{n}{T}}\rfloor\lfloor{\frac{m}{T}}\rfloor (d\in prime)\\&= \sum_{T=1}^n\;\lfloor{\frac{n}{T}}\rfloor\lfloor{\frac{m}{T}}\rfloor \sum_{d|T} \;\mu(\frac{T}{d})(d\in prime) \end{aligned} \]

发现 \(g(T)= \sum_{d|T} \mu(\frac{T}{d})\) 是可以通过枚举质数的倍数提前处理出来的，对它求前缀和，然后数论分块即可。

点击查看代码 
#include<bits/stdc++.h>
using namespace std;
const int Ratio=0;
const int N=1e7+5,M=1e7;
int u[N],pri[N],tot,yz[N],g[N];
void Wpre()
{
    u[1]=1;
    for(int i=2;i<=M;i++)
    {
        if(!yz[i]) pri[++tot]=i,u[i]=-1;
        for(int j=1;j<=tot;j++)
        {
            if(i*pri[j]>M) break;
            yz[i*pri[j]]=1;
            if(i%pri[j]==0) break;
            u[i*pri[j]]=-u[i];
        }
    }
    for(int i=1;i<=tot;i++) for(int j=1;j<=M;j++)
    {
        if(pri[i]*j>M) break;
        g[pri[i]*j]+=u[j];
    }
    for(int i=2;i<=M;i++) g[i]+=g[i-1];
}
int main()
{
    Wpre();
    int T,n,m;scanf("%d",&T);
    while(T--)
    {
        scanf("%d%d",&n,&m);
        if(n>m) swap(n,m);
        long long res=0;
        for(int l=1,r;l<=n;l=r+1)
            r=min(n/(n/l),m/(m/l)),res+=1ll*(g[r]-g[l-1])*(n/l)*(m/l);
        printf("%lld\n",res);
    }
    return Ratio;
}

约数个数和

个人感觉较简单的一道，要点在于求证 \(d(i,j)=\sum_{x|i}\sum_{y|j}[gcd(x,y)=1]\)

简证

image

摘自洛谷 Siyuan 的题解。

\[\begin{aligned}Ans&= \sum_{i=1}^n\sum_{j=1}^m\;d(ij)\\&= \sum_{i=1}^n\sum_{j=1}^m\sum_{x|i}\sum_{y|j}\;[gcd(x,y)=1]\\&=\sum_{i=1}^n\sum_{j=1}^m\sum_{x|i}\sum_{y|j}\sum_{d|gcd(x,y)}\mu(d)\\&=\sum_{i=1}^n\sum_{j=1}^m\;\lfloor{\frac{n}{i}}\rfloor\lfloor{\frac{m}{j}}\rfloor\sum_{d|gcd(i,j)}\mu(d)\\&=\sum_{d=1}^n\sum_{i=1}^{\lfloor{\frac{n}{d}}\rfloor}\sum_{j=1}^{\lfloor{\frac{m}{d}}\rfloor}\lfloor{\frac{n}{id}}\rfloor\lfloor{\frac{m}{jd}}\rfloor\mu(d)\\&=\sum_{d=1}^n\;\mu(d)\sum_{i=1}^{\lfloor{\frac{n}{d}}\rfloor}\sum_{j=1}^{\lfloor{\frac{m}{d}}\rfloor} \lfloor{\frac{n}{id}}\rfloor\lfloor{\frac{m}{jd}}\rfloor\\&=\sum_{d=1}^n\;\mu(d)\,(\sum_{i=1}^{\lfloor{\frac{n}{d}}\rfloor}\,\lfloor{\frac{n}{id}}\rfloor)\,(\sum_{j=1}^{\lfloor{\frac{m}{d}}\rfloor}\,\lfloor{\frac{m}{jd}}\rfloor) \end{aligned} \]

后面两坨看着就很熟悉是不是，没错就是数论分块的模板。这道题范围是 \(n\le 5\times10^5\)，比较小，所以直接对 \(\mu(d)\) 取前缀和，并 \(\mathcal{O(n\sqrt{n})}\) 预处理出后面 \(g(n)=\sum_{i=1}^n \lfloor{\frac{n}{i}}\rfloor\) 值，再运用数论分块即可，查询复杂度 \(\mathcal{O(T\sqrt{n})}\)。

点击查看代码 
#include<bits/stdc++.h>
typedef long long ll;
using namespace std;
const int Ratio=0;
const int N=5e4+5,M=5e4;
int u[N],pri[N],tot,yz[N],sum[N];
ll g[N];
void Wpre()
{
    u[1]=1;
    for(int i=2;i<=M;i++)
    {
        if(!yz[i]) pri[++tot]=i,u[i]=-1;
        for(int j=1;j<=tot;j++)
        {
            if(i*pri[j]>M) break;
            yz[i*pri[j]]=1;
            if(i%pri[j]==0) break;
            u[i*pri[j]]=-u[i];
        }
    }
    for(int i=1;i<=M;i++) sum[i]=sum[i-1]+u[i];
    for(int i=1;i<=M;i++)
    {
        ll res=0;
        for(int l=1,r;l<=i;l=r+1)
            r=i/(i/l),res+=1ll*(r-l+1)*(i/l);
        g[i]=res;
    }
}
int main()
{
    Wpre();
    int T,n,m;scanf("%d",&T);
    while(T--)
    {
        scanf("%d%d",&n,&m);
        if(n>m) swap(n,m);
        ll res=0;
        for(int l=1,r;l<=n;l=r+1)
            r=min(n/(n/l),m/(m/l)),res+=(sum[r]-sum[l-1])*g[n/l]*g[m/l];
        printf("%lld\n",res);
    }
    return Ratio;
}

数表

都说是板子题，但我打了三天。

按照能把题做下去的思路，我们先不考虑 \(a\) 的限制，设 \(\sigma (x)\) 为 \(x\) 的因数和，那么显然有：

\[\begin{aligned}Ans&=\sum_{i=1}^n\sum_{j=1}^m\;\sigma(gcd(i,j))\\&= \sum_{d=1}^n\sum_{i=1}^n\sum_{j=1}^m \;[gcd(i,j)=d]\;\sigma(d)\\&=\sum_{d=1}^n\sum_{i=1}^{\lfloor{\frac{n}{d}}\rfloor{}}\sum_{j=1}^{\lfloor{\frac{m}{d}}\rfloor{}}\;[gcd(i,j)=1]\;\sigma(d)\\&=\sum_{d=1}^n\sum_{i=1}^{\lfloor{\frac{n}{d}}\rfloor{}}\sum_{j=1}^{\lfloor{\frac{m}{d}}\rfloor{}}\sum_{x|gcd(i,j)}\mu(x)\,\sigma(d)\\&= \sum_{d=1}^n\;\sigma(d)\sum_{x=1}^{\lfloor{\frac{n}{d}}\rfloor{}}\;\mu(x)\lfloor{\frac{n}{dx}}\rfloor{}\lfloor{\frac{m}{dx}}\rfloor{} \end{aligned} \]

令 \(T=dx\)，则有：

\[\begin{aligned}Ans&= \sum_{T=1}^n \sum_{d|T} \;\sigma(d)\,\mu(\frac{T}{d})\,\lfloor{\frac{n}{T}}\rfloor{}\lfloor{\frac{m}{T}}\rfloor{}\\&=\sum_{T=1}^n\;\lfloor{\frac{n}{T}}\rfloor{}\lfloor{\frac{m}{T}}\rfloor{}\sum_{d|T}\;\sigma(d)\,\mu(\frac{T}{d}) \end{aligned} \]

不会求 \(\sigma\) 函数的看这里。

化到这一步，后面的部分已经可以提前处理了，此时考虑加入 \(a\) 的限制条件。

令 $g(T)=\sum_{d|T} \sigma(d),\mu(\frac{T}{d}) $，当 \(\sigma(d)\le a\) 时会对 \(g(T)\) 产生贡献，而对于所有的多测，改变的只是矩形的规格，并不是数值。。。

那么我们完全可以离线下来做！讲所有询问按 \(a\) 升序排序，每当 \(a\) 的值变大时，对有变动的地方插入相应的值。恰好我们所需的是前缀和且带修，所以考虑用树状数组去维护 \(g(T)\) 函数的前缀和。

插入值的复杂度近似于调和级数，树状数组更新复杂度为 \(\mathcal{O(\log n)}\)，故修改总复杂度为 \(\mathcal{O(n\log^2 n)}\)；采用数论分块优化查询 \(\mathcal{O(\sqrt{n})}\)，树状数组查询 \(\mathcal{O(\log n)}\)，故查询总复杂度为 \(\mathcal{O(q\sqrt{n}\log n)}\)。程序总复杂度 \(\mathcal{O(n\log^2 n+q\sqrt{n}\log n)}\)。

点击查看代码 
#include<bits/stdc++.h>
#define _(a,b) make_pair(a,b)
#define fi first
#define se second
typedef long long ll;
using namespace std;
const int Ratio=0;
const int N=1e5+5,M=1e5;
int u[N],pri[N],tot,yz[N],low[N];
int ans[N],v[N];
pair<int,int> o[N];
struct rmm{int n,m,a,id;}q[N];
bool cmp(rmm A,rmm B){return A.a<B.a;}
void Wpre()
{
    u[1]=1,o[1]=_(1,1);
    for(int i(2);i<=M;i++)
    {
        if(!yz[i]) pri[++tot]=i,u[i]=-1,low[i]=i+1,o[i]=_(i+1,i);
        for(int j(1);j<=tot;j++)
        {
            if(i*pri[j]>M) break;
            yz[i*pri[j]]=1;
            if(i%pri[j]==0)
            {
                low[i*pri[j]]=low[i]*pri[j]+1;
                o[i*pri[j]]=_(o[i].fi/low[i]*low[i*pri[j]],i*pri[j]);
                break;
            }
            u[i*pri[j]]=-u[i];
            low[i*pri[j]]=pri[j]+1;
            o[i*pri[j]]=_(o[i].fi*o[pri[j]].fi,i*pri[j]);
        }
    }
    sort(o+1,o+M+1);
}
void Wupd(int x,int k){for(;x<=M;x+=(x&-x)) v[x]+=k;}
int Wq(int x)
{
    int res=0;
    while(x) res+=v[x],x-=(x&-x);
    return res;
}
int main()
{
    Wpre();
    int T;scanf("%d",&T);
    for(int i(1);i<=T;i++) scanf("%d%d%d",&q[i].n,&q[i].m,&q[i].a),q[i].id=i;
    sort(q+1,q+1+T,cmp);
    int j=1;
    for(int i(1);i<=T;i++)
    {
        while(j<=M&&o[j].fi<=q[i].a)
        {
            for(int k(o[j].se);k<=M;k+=o[j].se)
                Wupd(k,o[j].fi*u[k/o[j].se]);
            j++;
        }
        int res=0,n=q[i].n,m=q[i].m;
        if(n>m) swap(n,m);
        for(int l(1),r;l<=n;l=r+1)
            r=min(n/(n/l),m/(m/l)),res+=(Wq(r)-Wq(l-1))*(n/l)*(m/l);
        ans[q[i].id]=res;
    }
    for(int i(1);i<=T;i++) printf("%d\n",(ans[i]&(~(1<<31))));
    return Ratio;
}

DZY Loves Math

有新函数，先不管它，推式子。

\[\begin{aligned}Ans&=\sum_{i=1}^n\sum_{j=1}^m \;f(gcd(i,j))\\ &=\sum_{d=1}^n\sum_{i=1}^n\sum_{j=1}^m\;[gcd(i,j)=d] \,f(d)\\&=\sum_{d=1}^n\sum_{i=1}^{\lfloor{\frac{n}{d}}\rfloor{}}\sum_{j=1}^{\lfloor{\frac{m}{d}}\rfloor{}}\;[gcd(i,j)=1]\,f(d)\\&=\sum_{d=1}^n\;f(d)\sum_{i=1}^{\lfloor{\frac{n}{d}}\rfloor{}}\sum_{j=1}^{\lfloor{\frac{m}{d}}\rfloor{}}\sum_{x|gcd(i,j)}\mu(x)\\&=\sum_{d=1}^n\;f(d)\sum_{x=1}^{\lfloor{\frac{n}{d}}\rfloor{}}\;\mu(x)\lfloor{\frac{n}{dx}}\rfloor{}\lfloor{\frac{m}{dx}}\rfloor{} \end{aligned} \]

令 \(T=dx\)，则有：

\[\begin{aligned}Ans&=\sum_{T=1}^n\;\lfloor{\frac{n}{T}}\rfloor{}\lfloor{\frac{m}{T}}\rfloor{}\sum_{d|T}\;f(d)\,\mu(\frac{T}{d}) \end{aligned} \]

到这一步似乎就可以做了，不过后面的那坨 $g(T)=\sum_{d|T}\ f(d)\ \mu(\frac{T}{d}) $ 需要调和级数预处理，\(n\le 10^7\)，在一个月前的题库上是能跑过的，不过现在不知道为什么速度慢了一半，所以需要我们继续优化，用线性复杂度处理出。

我们着眼于 \(g(T)\) 这个函数。设 \(T\) 共有 \(k\) 个不同的质因子，将其分解为 $T=\prod_{i=1}^k\ p_i^{a_i} $ 。令 \(T\) 的一个因子 \(d\) 分解为 \(d=\prod_{i=1}^k \;p_i^{b_i}\)，发现 \(d\) 有贡献当且仅当 \(\mu(\frac{T}{d})\neq 0\)，即 \(\forall i\in[1,k], \;a_i-b_i\le 1\)，换句话说 \(b_i\) 的取值只能是 \(a_i\) 或 \(a_i-1\)。

然后再分两种情况：

一种是 \(a_i\) 不全相等的情况，此时必有 $存在 \ i,j\in\ [1,k],a_i\lt a_j $，考虑 \(b_i\) 无论取 \(a_i\) 或 \(a_i-1\)，不会对 \(f(T)\) 的值产生影响，且两种取法的 \(\mu\) 值之和为 \(0\)，即：此时无论如何分配 \(b_i\) 的值，最终结果总能两两分配出相加为 \(0\) 的情况，故 \(g(T)=0\)。

另外就是 \(a_i\) 全部相等的情况。此时有贡献的 \(d\) 共有 \(2^k\) 种取值，已知 \(\binom{k}{0}+\binom{k}{2}+\cdots+\binom{k}{2\times \lfloor{\frac{k}{2}}\rfloor{}}=\binom{k}{1}+\binom{k}{3}+\cdots+\binom{k}{2\times \lfloor{\frac{k-1}{2}}\rfloor{}+1}\)，即杨辉三角中每行的奇数项和等于偶数项和。括号中上面的数恰好可对应 \(a_i-b_i\) 中值为 \(1\) 的数量，即在 \(f(d)\) 均相同的情况下，\(g(T)\) 的值是 \(0\)。但当括号中上面的数等于 \(k\)，即 \(\forall \;i\in[1,k],b_i=a_i-1\) 时，\(f(d)=a_i-1\)。

稍微理解下，然后再来分情况。若 \(k\) 为奇数，则当 \(\forall\;i\in[1,k],b_i=a_i-1\) 时，\(\mu(\frac{T}{d})=(-1)^k=-1\)，也就是说，在对结果产生贡献的那一步中本应该减去 \(a_i\)，但实际上减去了 \(a_i-1\)，那么最终 \(g(T)=1\)；若 \(k\) 为偶数同理，此时 \(\mu(\frac{T}{d})=1\)，那么本应该加上 \(a_i\)，但实际上只加了 \(a_i-1\)，最终 \(g(T)=-1\)。

整合一下，发现 \(g\) 函数的取值如下：

\[g(T)=\begin{cases}0\quad\quad\quad\quad 存在 \; i,j\in[1,k],a_i\lt a_j\\ (-1)^{k+1}\quad\, \forall\;i,j\in[1,k],a_i=a_j \end{cases} \]

其中 \(k\) 为 \(T\) 的质因子个数。

实现中，记录最小质因子的指数和它的指数次幂，转移时判断除去最小质因子后次小质因子的指数是否与之相同即可。

预处理 \(g\) 函数及前缀和的复杂度为 \(\mathcal{O(n)}\)，数论分块优化后查询复杂度为 \(\mathcal{O(T\sqrt{n})}\)，总复杂度 \(\mathcal{O(n+T\sqrt{n})}\)。

点击查看代码 
#include<bits/stdc++.h>
typedef long long ll;
using namespace std;
const int Ratio=0;
const int N=1e7+5,M=1e7;
int n,m;
int u[N],pri[N],tot,yz[N],Index[N],low[N],g[N];
void Wpre()
{
    u[1]=1,g[1]=0;
    for(int i(2);i<=M;i++)
    {
        if(!yz[i]) low[i]=pri[++tot]=i,u[i]=-1,Index[i]=g[i]=1;
        for(int j(1);j<=tot;j++)
        {
            if(i*pri[j]>M) break;
            yz[i*pri[j]]=1;
            if(i%pri[j]==0)
            {
                Index[i*pri[j]]=Index[i]+1;
                low[i*pri[j]]=low[i]*pri[j];
                if(low[i]==i) g[i*pri[j]]=1;
                else g[i*pri[j]]=Index[i/low[i]]==Index[i*pri[j]]?-g[i/low[i]]:0;
                break;
            }
            u[i*pri[j]]=-u[i];
            Index[i*pri[j]]=1;
            low[i*pri[j]]=pri[j];
            g[i*pri[j]]=Index[i]==1?-g[i]:0;
        }
    }
    for(int i(2);i<=M;i++) g[i]+=g[i-1];
}
int main()
{
    Wpre();
    int T;scanf("%d",&T);
    while(T--)
    {
        scanf("%d%d",&n,&m);
        if(n>m) swap(n,m);
        ll res=0;
        for(int l(1),r;l<=n;l=r+1)
            r=min(n/(n/l),m/(m/l)),res+=1ll*(g[r]-g[l-1])*(n/l)*(m/l);
        printf("%lld\n",res);
    }
    return Ratio;
}

数字表格

从这道题开始就要接触含 \(\prod\) 的式子了。

令 \(Fib_i\) 表示 Fibonacci 数列的第 \(i\) 项。

\[\begin{aligned}Ans&= \prod_{i=1}^n\prod_{j=1}^m Fib_{gcd(i,j)}\\&=\prod_{d=1}^n\prod_{i=1}^n\prod_{j=1}^m\ [gcd(i,j)=d]\ Fib_d\\&=\prod_{d=1}^n\prod_{i=1}^{\lfloor{\frac{n}{d}}\rfloor{}}\prod_{j=1}^{\lfloor{\frac{m}{d}}\rfloor{}}\ [gcd(i,j)=1]\ Fib_d\\&=\prod_{d=1}^n \ Fib_d^{\sum_{i=1}^{\lfloor{\frac{n}{d}}\rfloor{}}\sum_{j=1}^{\lfloor{\frac{m}{d}}\rfloor{}}\ [gcd(i,j)=1]} \end{aligned} \]

上面那一坨已经不陌生了，拿下来化简。

\[\begin{aligned}\sum_{i=1}^{\lfloor{\frac{n}{d}}\rfloor{}}\sum_{j=1}^{\lfloor{\frac{m}{d}}\rfloor{}}\sum_{x|gcd(i,j)} \mu{(x)}&=\sum_{x=1}^{\lfloor{\frac{n}{d}}\rfloor{}}\ \mu{(x)}\lfloor{\frac{n}{dx}}\rfloor{}\lfloor{\frac{m}{dx}}\rfloor{} \end{aligned} \]

令 \(T=dx\)，放上去化简有：

\[\begin{aligned} Ans&=\prod_{d=1}^n\ Fib_d^{\sum_{x=1}^{\lfloor{\frac{n}{d}}\rfloor{}}\mu{(x)}\lfloor{\frac{n}{T}}\rfloor{}\lfloor{\frac{m}{T}}\rfloor{}}\\&=\prod_{T=1}^n\prod_{d|T}\ Fib_d^{\mu{(\frac{T}{d})}\lfloor{\frac{n}{T}}\rfloor{}\lfloor{\frac{m}{T}}\rfloor{}}\\ &=\prod_{T=1}^n(\prod_{d|T}\ Fib_d^{\mu{(\frac{T}{d})}})^{\lfloor{\frac{n}{T}}\rfloor{}\lfloor{\frac{m}{T}}\rfloor{}} \end{aligned} \]

括号外层的指数部分已经可以用数论分块解决了，考虑如何得到内层的 $g(T)=\prod_{d|T}Fib_d^{\mu{(\frac{T}{d})}} $ 函数。发现 \(n\le 10^6\)，不是很大，可调和级数预处理，枚举每一个因数 \(d\)，枚举它的倍数 \(T\)，可以实现。

预处理部分调和级数复杂度约为 \(\mathcal{O(n\ln n)}\)，查询数论分块复杂度为 \(\mathcal{O(T\sqrt{n})}\)，总复杂度 \(\mathcal{O(n\ln n+T\sqrt{n})}\)。

点击查看代码 
#include<bits/stdc++.h>
typedef long long ll;
using namespace std;
const int Ratio=0;
const int N=1e6+5,M=1e6;
const int mod=1e9+7;
int n,m;
int u[N],pri[N],tot,fib[N],invf[N],g[N];
bool yz[N];
int Wqp(ll x,int y)
{
    ll res=1;
    while(y){if(y&1) res=res*x%mod;x=x*x%mod;y>>=1;}
    return res;
}
void Wpre()
{
    u[1]=1,g[0]=g[1]=1,fib[1]=1,invf[1]=1;
    for(int i(2);i<=M;i++)
    {
        fib[i]=(fib[i-1]+fib[i-2])%mod;
        invf[i]=Wqp(fib[i],mod-2);
        g[i]=1;
        if(!yz[i]) pri[++tot]=i,u[i]=-1;
        for(int j(1);j<=tot;j++)
        {
            if(i*pri[j]>M) break;
            yz[i*pri[j]]=1;
            if(i%pri[j]==0) break;
            u[i*pri[j]]=-u[i];
        }
    }
    for(int i(1);i<=M;i++)
    {
        if(!u[i]) continue;
        for(int j(1);j<=M/i;j++)
            g[i*j]=1ll*g[i*j]*(u[i]==1?fib[j]:invf[j])%mod;
    }
    for(int i(2);i<=M;i++) g[i]=1ll*g[i]*g[i-1]%mod;
}
int main()
{
    Wpre();
    int T;scanf("%d",&T);
    while(T--)
    {
        scanf("%d%d",&n,&m);
        if(n>m) swap(n,m);
        ll res=1;
        for(int l(1),r;l<=n;l=r+1)
            r=min(n/(n/l),m/(m/l)),
            res=1ll*res*Wqp(1ll*g[r]*Wqp(g[l-1],mod-2)%mod,1ll*(n/l)*(m/l)%(mod-1))%mod;
        printf("%lld\n",(res+mod)%mod);
    }
    return Ratio;
}

于神之怒加强版

这题难点同样不在推式子，而是在后续函数的预处理。

\[\begin{aligned}Ans&=\sum_{i=1}^n\sum_{j=1}^m\ gcd(i,j)^k\\&=\sum_{d=1}^n\sum_{i=1}^{\lfloor{\frac{n}{d}}\rfloor{}}\sum_{j=1}^{\lfloor{\frac{m}{d}}\rfloor{}}\ [gcd(i,j)=1]\ d^k\\&=\sum_{d=1}^n\ d^k\sum_{x=1}^{\lfloor{\frac{n}{d}}\rfloor{}}\ \mu{(x)}\lfloor{\frac{n}{dx}}\rfloor{}\lfloor{\frac{m}{dx}}\rfloor{} \end{aligned} \]

令 \(T=dx\)：

\[Ans=\sum_{T=1}^n\ \lfloor{\frac{n}{T}}\rfloor{}\lfloor{\frac{m}{T}}\rfloor{}\sum_{d|T}\ d^k\mu{(\frac{T}{d})} \]

考虑 $g(T)=\sum_{d|T}d^k\mu{(\frac{T}{d})} $ 这个函数。因为该函数为两个积性函数的积，所以它也是个积性函数。考虑一个质数 \(p\)，有 $g(pa)=\sum_{i=0}a\ (pi)k\mu{(p^{a-i})} $。发现 \(a-i\ge 2\) 时 \(\mu{(p^{a-i})}\) 均为 \(0\)，所以只需考虑 \(i=a\) 和 \(a-1\) 的贡献，此时可以简化为 $g(pa)=(pa)k-(p)^k $。发现和 \(\mu\) 没关系，震惊，莫反题不需要筛莫比乌斯函数。

考虑由积性函数的性质推出转移式：取任意一数 \(x\)，设其本质不同质因子数量为 \(s\)，即 $x=\prod_{i=1}^s p_i^{a_i} $，则 $g(x)=\prod_{i=1}^s g(p_i^{a_i}) $。也就是说，只要保证转移来源的两个因数互质即可。实现中仍然记录一个 \(low_i\) 表示 \(i\) 的最小质因子的幂，转移 $g_{i\times pri_j}= g_{\frac{i}{low_i}}\times g_{low_i\times pri_j} $。

预处理是接近线性，数论分块处理询问，总复杂度 \(\mathcal{O(n+T\sqrt{n})}\)。

点击查看代码 
#include<bits/stdc++.h>
using namespace std;
const int Ratio=0;
const int N=5e6+5,M=5e6;
const int mod=1e9+7;
int k;
int pri[N],tot,low[N],g[N];
bool yz[N];
int Wqp(int x,int y)
{
    int res=1;
    while(y){if(y&1) res=1ll*res*x%mod;x=1ll*x*x%mod;y>>=1;}
    return res;
}
void Wpre()
{
    g[1]=1;
    for(int i(2);i<=M;i++)
    {
        if(!yz[i]) low[i]=pri[++tot]=i,g[i]=Wqp(i,k)-1;
        for(int j(1);j<=tot;j++)
        {
            if(i*pri[j]>M) break;
            yz[i*pri[j]]=1;
            if(i%pri[j]==0)
            {
                low[i*pri[j]]=low[i]*pri[j];
                if(low[i]==i) g[i*pri[j]]=1ll*(Wqp(i*pri[j],k)-Wqp(i,k)+mod)%mod;
                else g[i*pri[j]]=1ll*g[i/low[i]]*g[low[i]*pri[j]]%mod;
                break;
            }
            low[i*pri[j]]=pri[j];
            g[i*pri[j]]=1ll*g[i]*g[pri[j]]%mod;
        }
    }
    for(int i(2);i<=M;i++) g[i]=(g[i]+g[i-1])%mod;
}
int main()
{
    int T;scanf("%d%d",&T,&k);
    Wpre();
    while(T--)
    {
        int n,m;scanf("%d%d",&n,&m);
        if(n>m) swap(n,m);
        int res=0;
        for(int l(1),r;l<=n;l=r+1)
            r=min(n/(n/l),m/(m/l)),
            res=1ll*(res+1ll*(g[r]-g[l-1]+mod)%mod*(n/l)%mod*(m/l)%mod)%mod;
        printf("%d\n",(res+mod)%mod);
    }
    return Ratio;
}

jzptab

第一道自己推出来的，%%% 5k 的引导。

\[\begin{aligned}Ans&= \sum_{i=1}^n\sum_{j=1}^m\ lcm(i,j)\\&=\sum_{i=1}^n\sum_{j=1}^m\ \frac{i\times j}{gcd(i,j)}\\&=\sum_{d=1}^n\ \frac{1}{d}\ \sum_{i=1}^n\sum_{j=1}^m\ i\times j\ [gcd(i,j)=d]\\&=\sum_{d=1}^n\ \frac{1}{d}\ \sum_{i=1}^{\lfloor{\frac{n}{d}}\rfloor{}}\sum_{j=1}^{\lfloor{\frac{m}{d}}\rfloor{}}\ d^2\times i\times j\ \sum_{x|d}\ \mu{(x)}\\&=\sum_{d=1}^n\ d\ \sum_{x=1}^{\lfloor{\frac{n}{d}}\rfloor{}}\ \mu{(x)}\ \sum_{i=1}^{\lfloor{\frac{n}{dx}}\rfloor{}}\ i\times x\ \sum_{j=1}^{\lfloor{\frac{m}{dx}}\rfloor{}}\ j\times x\\&=\sum_{d=1}^n\ d\ \sum_{x=1}^{\lfloor{\frac{n}{d}}\rfloor{}}\ \mu{(x)}\times x^2\ \sum_{i=1}^{\lfloor{\frac{n}{dx}}\rfloor{}}\ i\ \sum_{j=1}^{\lfloor{\frac{m}{dx}}\rfloor{}}\ j \end{aligned} \]

发现 $\sum_{i=1}^x\ i=\frac{x(x+1)}{2} $，记 \(s_i=\sum_{j=1}^i\ j\)，并令 \(T=dx\)，则有：

\[\begin{aligned}Ans&=\sum_{d=1}^n\ d\ \sum_{x=1}^{\lfloor{\frac{n}{d}}\rfloor{}}\ \mu{(x)}\times x^2\times s_{\lfloor{\frac{n}{dx}}\rfloor{}}\times s_{\lfloor{\frac{m}{dx}}\rfloor{}}\\&=\sum_{T=1}^n\ s_{\lfloor{\frac{n}{T}}\rfloor{}}\times s_{\lfloor{\frac{m}{T}}\rfloor{}}\ \sum_{x|T}\ \mu{(x)}\times x^2\times \frac{T}{x}\\&=\sum_{T=1}^n\ s_{\lfloor{\frac{n}{T}}\rfloor{}}\times s_{\lfloor{\frac{m}{T}}\rfloor{}}\times T\ \sum_{x|T}\ \mu{(x)}\times x \end{aligned} \]

推到这一步我们完全可以记函数 $g(T)=T\times \sum_{x|T}\ \mu{(x)}\times x $，筛出 \(\mu\) 后调和级数预处理即可题库上通过本题，不过时间复杂度很不优秀，我们若想在 3000ms 内通过同题 Crash 的数字表格，必须进一步优化。

拿出 $g(T) $ 这个函数，首先由定义可知其为奇函数，那么进一步探索。设 \(x=p^a\)，则 $g(x)=x\times [\mu{(p^0)}\times p0+\mu{(p1)}\times p^1+\cdots +\mu{(p^a)}\times p^a] $，发现中括号内只有前两项是有贡献的，化简可写为 \(g(x)=x\times (1-p)=p^a\times (1-p)\)。那么很显然，若设 \(y=p^{a-1}\)，则 $g(y)=p^{a-1}\times(1-p) $，得到 \(g(p^a)=g(p^{a-1})\times p\)。这样，我们就得到了函数 \(g\) 在所有条件下的转移式，可以在线性时间内求得 \(g\) 函数了！

最后注意坑点：模数是 1e8+9！

预处理复杂度 \(\mathcal{O(n)}\)，询问数论分块复杂度 \(\mathcal{O(T\sqrt{n})}\)，总复杂度 \(\mathcal{O(n+T\sqrt{n})}\)，可以通过本题。

点击查看代码 
#include<bits/stdc++.h>
typedef long long ll;
using namespace std;
const int Ratio=0;
const int N=1e7+5,M=1e7;
const int mod=1e8+9;
int k;
int pri[N],tot;
ll s[N],g[N];
bool yz[N];
void Wpre()
{
    g[1]=1,s[1]=1;
    for(int i(2);i<=M;i++)
    {
        s[i]=1ll*i*(i+1)/2%mod;
        if(!yz[i]) pri[++tot]=i,g[i]=(i-1ll*i*i%mod+mod)%mod;
        for(int j(1);j<=tot;j++)
        {
            if(i*pri[j]>M) break;
            yz[i*pri[j]]=1;
            if(i%pri[j]==0)
            {
                g[i*pri[j]]=g[i]*pri[j]%mod;
                break;
            }
            g[i*pri[j]]=g[i]*g[pri[j]]%mod;
        }
    }
    for(int i(2);i<=M;i++) g[i]=(g[i]+g[i-1])%mod;
}
int main()
{
    Wpre();
    int T;scanf("%d",&T);
    while(T--)
    {
        int n,m;scanf("%d%d",&n,&m);
        if(n>m) swap(n,m);
        ll res=0;
        for(int l(1),r;l<=n;l=r+1)
            r=min(n/(n/l),m/(m/l)),
            res=(res+s[n/l]*s[m/l]%mod*(g[r]-g[l-1]+mod)%mod)%mod;
        printf("%lld\n",res);
    }
    return Ratio;
}

Steps to One

CF1139D 用到一点莫反也是莫反。

题目大意：每次从 \(\left[1,n\right]\) 随机取一个数加入数组 \(a_i\)，当 \(gcd_{i=1}^{len}\ a_i = 1\) 时停止，问 \(len\) 的期望。

直接用期望式子推：

\[\begin{aligned}ans&=\sum_{i=1}^{\infty}\ P_{(len=i)}\times i \\ &=\sum_{i = 1}^{\infty}\sum_{j=1}^i\ P_{(len=i)}\\&=\sum_{j=1}^{\infty}\sum_{i\ge j}\ P_{(len=i)}\\&=\sum_{i=1}^{\infty}\ P_{(len\ge i)}\\&=1+\sum_{i=1}^{\infty}\ P_{(len\gt i)}\\&=1+\sum_{i=1}^{\infty}\ P_{(gcd_{j=1}^i\ a_j\gt 1)}\\&= 1+\sum_{i=1}^{\infty}1-P_{(gcd_{j=1}^i\ a_j=1)}\\&=1+\sum_{i=1}^{\infty}1-\frac{\sum_{a_1=1}^n\sum_{a_2=1}^n\cdots\sum_{a_i=1}^n\ \left[gcd_{j=1}^i\ a_j=1 \right]}{n^i}\\&=1+\sum_{i=1}^{\infty}1-\frac{\sum_{a_1=1}^n\sum_{a_2=1}^n\cdots\sum_{a_i=1}^n\ \sum_{d|gcd_{j=1}^i\ a_j}\ \mu{(d)}}{n^i}\\&=1+\sum_{i=1}^{\infty}1-\frac{\sum_{d=1}^n\ \mu{(d)}\ (\lfloor{\frac{n}{d}}\rfloor{})^i}{n^i}\\&=1-\sum_{i=1}^{\infty}\frac{\sum_{d=2}^n\ \mu{(d)}\ (\lfloor{\frac{n}{d}}\rfloor{})^i}{n^i}\\&=1-\sum_{i=1}^{\infty}\frac{1}{n^i}\sum_{d=2}^n\ \mu{(d)}\ (\lfloor{\frac{n}{d}}\rfloor{})^i\\&=1-\sum_{d=2}^n\ \mu{(d)}\sum_{i=1}^{\infty}\left(\frac{\lfloor{\frac{n}{d}}\rfloor{}}{n}\right)^i \end{aligned} \]

发现最后是个无穷等比数列，可以如下求值：

\[S=\sum_{i=1}^{\infty}\ x^i \]

\[Sx=\sum_{i=2}^{\infty}\ x^i \]

\[S-Sx=x \]

\[S=\frac{x}{1-x} \]

然后代入即可：

\[\begin{aligned}ans&=1-\sum_{d=2}^n\ \mu{(d)}\frac{\lfloor{\frac{n}{d}}\rfloor{}}{n-\lfloor{\frac{n}{d}}\rfloor{}} \end{aligned} \]

这道题范围很小，只有 \(10^5\)，\(\mathcal{O(n)}\) 跑就行了，整除分块优化能到 \(\mathcal{O(\sqrt{n})}\)，但筛还是 \(\mathcal{O(n)}\) 的。

点击查看代码 
#include<bits/stdc++.h>
#define fo(x, y, z) for(register int (x) = (y); (x) <= (z); (x)++)
#define fu(x, y, z) for(register int (x) = (y); (x) >= (z); (x)--)
using namespace std;
typedef long long ll;
#define lx ll
inline lx qr()
{
	char ch = getchar();lx x = 0 , f = 1;
	for(;ch<'0'||ch>'9';ch = getchar()) if(ch == '-') f = -1;
	for(;ch>= '0' && ch<= '9';ch = getchar()) x = (x<<3) + (x<<1) + (ch^48);
	return x*f;
}
#undef lx
#define qr qr()
#define pii pair<int, int>
#define M_P(a, b) make_pair(a, b)
#define fi first
#define se second
#define P_B(x) push_back(x)
const int Ratio = 0;
const int N = 1e5 + 5, M = 1e5;
const int mod = 1e9 + 7;
int n;
int pri[N], tot, u[N];
bool yz[N];
ll ans = 1;
namespace Wisadel
{
    ll Wqp(ll x, int y)
    {
        ll res = 1;
        while(y){if(y & 1) res = res * x % mod; x = x * x % mod; y >>= 1;}
        return res;
    }
    void Wpre()
    {
        u[1] = 1;
        fo(i, 2, M)
        {
            if(!yz[i]) pri[++tot] = i, u[i] = -1;
            fo(j, 1, tot)
            {
                if(i * pri[j] > M) break;
                yz[i * pri[j]] = 1;
                if(i % pri[j] == 0) break;
                u[i * pri[j]] = -u[i];
            }
        }
    }
    short main()
    {
        // freopen(".in", "r", stdin) , freopen(".out", "w", stdout);
        // freopen(".err", "w", stderr);
        n = qr;
        Wpre();
        ll zc;
        fo(i, 2, n)
            zc = n / i, ans = (ans - u[i] * zc % mod * Wqp(n - zc, mod - 2) % mod + mod) % mod;
        printf("%lld\n", ans);
        return Ratio;
    }
}
signed main(){return Wisadel::main();}

其他例题会逐渐更新。

有道理

数论既然上不了首页，那就推首歌吧

  • 有道理——李荣浩

太理智的消费 对不对
成夜不睡
没一句说的对 也很对
这个年岁
有爱的人不追 还斗嘴
显得般配
有时也为自己感到惭愧
穿的也不算贵 黑白灰
尽量自费
被自己往里推 动动嘴
也能后退
总会被动吃亏 不防备
深知有罪
我现在要开始表现伤悲和承认理亏
人文和法规
其实都算是有理想啊
都是废寝食忘啊
有深度有主张啊
难免会累会慌
难免也会迷茫
太理智的消费 对不对
成夜不睡
没一句说的对 也很对
这个年岁
有爱的人不追 还斗嘴
显得般配
有时也为自己感到惭愧
穿的也不算贵 黑白灰
尽量自费
被自己往里推 动动嘴
也能后退
总会被动吃亏 不防备
深知有罪
我现在要开始表现伤悲和承认理亏
人文和法规
其实都算是有理想啊
都是废寝食忘啊
有深度有主张啊
难免会累会慌
难免也会迷茫
其实都算是有理想啊
都是废寝食忘啊
有深度有主张啊
难免会累会慌
难免也会迷茫
其实都算是有理想啊
都是废寝食忘啊
有深度有主张啊
难免会累会慌
难免也会迷茫
其实都算是有理想啊
都是废寝食忘啊
有深度有主张啊
难免会累会慌
难免也会迷茫
其实都算是有理想啊

制作不易，如有帮助，请推荐支持，感谢！

posted @ 2024-09-13 10:56  DrRatio  阅读(218)  评论(16)  编辑  收藏  举报