数学(2)——莫比乌斯反演做题记录

Reywmp · 2021-02-11 19:55:43 · 个人记录

Luogu P1447 [NOI2010] 能量采集

求

\sum_{i=1}^n\sum_{j=1}^m(2\times\gcd{(i,j)-1})

可以转换为

-(n\times m)+2\times{\color{blue}{\sum_{i=1}^n\sum_{j=1}^m\gcd(i,j)}}

考虑对蓝色部分莫比乌斯反演

设

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

这个 g(\mathrm d) 可以帮我们做数论分块。

首先可以看出前面的部分

\texttt{Part}=\sum_{\mathrm d=1}^nf(\mathrm d)\times \mathrm d

我们可以直接枚举 i 而不是 i|d

所以有

f(\mathrm d)=\sum_{i=1}^{\lfloor\frac{n}{\mathrm d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{\mathrm d}\rfloor}\sum_{k|\gcd(i,j)}\mu(k)\\ f(\mathrm d)=\sum_{k=1}^{\lfloor\frac{n}{\mathrm d}\rfloor}\mu(k)\times\lfloor\frac{n}{k\times \mathrm d}\rfloor\lfloor\frac{m}{k\times \mathrm d}\rfloor

答案就是

\texttt{ans}=[2\sum_{\mathrm d=1}^n\mathrm d\times\sum_{k=1}^{\lfloor\frac{n}{\mathrm d}\rfloor}(\mu(k)\lfloor\frac{n}{k\times \mathrm d}\rfloor\lfloor\frac{m}{k\times \mathrm d}\rfloor)]-n\times m

代码通过数论分块和线性筛实现可以就通过这题。并且首先保证 \color{red}n\leq m 。

#include<cstdio>
#include<cstring>
#include<cmath>
#include<iostream>
#include<algorithm>

#define ll long long

using namespace std;

inline int read()
{
    int x=0,w=1;char ch=getchar();
    while((ch<'0'||ch>'9')&&ch!='-')ch=getchar();
    if(ch=='-')ch=getchar(),w=-1;
    while(ch>='0'&&ch<='9')x=(x<<1)+(x<<3)+ch-48,ch=getchar();
    return x*w;
}

int n,m,tot;
const int N=100005;
ll phi[N],pri[N],sum[N];
bool vis[N];

void Phi()
{
    phi[1]=1;
    for(int i=2;i<=100000;i++)
    {
        if(!vis[i])pri[++tot]=i,phi[i]=i-1;
        for(int j=1;j<=tot&&i*pri[j]<=100000;j++)
        {
            if(i*pri[j]>100000)continue;
            vis[i*pri[j]]=1;
            if(i%pri[j])phi[i*pri[j]]=phi[i]*phi[pri[j]];
            else {phi[i*pri[j]]=phi[i]*pri[j];continue;}
        }
    }
}

int main()
{
    n=read();m=read();
    Phi();
    for(int i=1;i<=100000;i++)sum[i]=sum[i-1]+phi[i];
    ll ans=0;
    ll l=1,r=0;
    while(l<=min(n,m))
    {
        r=min(n/(n/l),m/(m/l));
        ans+=(ll)(sum[r]-sum[l-1])*(n/l)*(m/l);
        l=r+1; 
    }
    printf("%lld\n",2ll*ans-(ll)n*m);
    return 0;
}

Luogu P2257 YY的GCD

这个题有个神必的巧妙做法，参考这篇题解。

这个题是

\sum_{i=1}^n\sum_{j=1}^m[\gcd(i,j)\in \texttt{Prime}]

(下文质数集合用 \mathbb P 表示且首先同样保证 \color{red}n\leq m)

我们考虑直接枚举质数：

\sum_{p\in\mathbb P}\sum_{i=1}^n\sum_{j=1}^m[\gcd(i,j)=p]

我们考虑莫比乌斯反演技巧，设：

f(\mathrm d)=\sum_{i=1}^n\sum_{j=1}^m[\gcd(i,j)=\mathrm d]\\ g(x)=\sum_{x|\mathrm d}f(\mathrm d)=\lfloor\frac{n}{x}\rfloor\lfloor\frac{m}{x}\rfloor\\ f(x)=\sum_{x|\mathrm d}\mu(\lfloor\frac{\mathrm d}{x}\rfloor)g(\mathrm d)

那我们发现答案就是

\texttt{ans}=\sum_{p\in\mathbb P}f(p)

再对答案进行莫比乌斯反演

\begin{aligned} \texttt{ans}&=\sum_{p\in\mathbb P}\sum_{p|\mathrm d}\mu(\lfloor\frac{\mathrm d}{p}\rfloor)g(\mathrm d)\\ &=\sum_{p\in\mathbb P}\sum_{\mathrm d=1}^{\lfloor\frac{n}{p}\rfloor}\mu(\mathrm d)g(p\times \mathrm d)\\ &=\sum_{p\in\mathbb P}\sum_{\mathrm d=1}^{\lfloor\frac{n}{p}\rfloor}\mu(\mathrm d)\lfloor\frac{n}{p\times \mathrm d}\rfloor\lfloor\frac{m}{p\times \mathrm d}\rfloor\\ \end{aligned}

消去 p\times \mathrm d 可以考虑换元为 k

\begin{aligned} \texttt{ans}&=\sum_{k=1}^n\sum_{p|k,p\in\mathbb P}\mu(\lfloor\frac{k}{p}\rfloor)\lfloor\frac{n}{k}\rfloor\lfloor\frac{m}{k}\rfloor\\ &=\sum_{k=1}^n\lfloor\frac{n}{k}\rfloor\lfloor\frac{m}{k}\rfloor(\sum_{p|k}\mu(\lfloor\frac{k}{p}\rfloor)) \end{aligned}

这样数论分块，筛法，前缀和就可以实现了。

#include<cstdio>
#include<cstring>
#include<cmath>
#include<iostream>
#include<algorithm>

#define ll long long

using namespace std;

inline int read()
{
    int x=0;int w=1;char ch=getchar();
    while((ch<'0'||ch>'9')&&ch!='-')ch=getchar();
    if(ch=='-')ch=getchar(),w=-1;
    while(ch>='0'&&ch<='9')x=(x<<1)+(x<<3)+ch-48,ch=getchar();
    return x*w;
}

int n,m;
const int N=10000000;
const int M=10000005;

ll pri[M/10],phi[M],s[M],g[M];
bool vis[M];
ll mn;
ll ans,cnt;

void Phi()
{
    phi[1]=1;
    for(int i=2;i<=mn;i++)
    {
        if(!vis[i])phi[i]=-1,pri[++cnt]=i;
        for(int j=1;j<=cnt&&i*pri[j]<=mn;j++)
        {
            vis[i*pri[j]]=1;
            if(i%pri[j]==0)break;
            else phi[i*pri[j]]=-phi[i];
        }
    }
    for(int i=1;i<=cnt;i++)
    {
        for(int j=1;j*pri[i]<=mn;j++)
        g[j*pri[i]]+=phi[j];
    }
    for(int i=1;i<=mn;i++)
    {
        s[i]=s[i-1]+g[i];
    }
}

int main()
{
    int T;
    T=read();
    mn=10000000;
    Phi();ll ans=0;ll r;
    while(T--)
    {
        n=read();m=read();
        mn=(ll)min(n,m);ans=0;
        for(int i=1;i<=mn;i=r+1)
        {
            r=min(n/(n/i),m/(m/i));
            ans+=(ll)(s[r]-s[i-1])*(ll)(n/i)*(ll)(m/i);
        }
        printf("%lld\n",ans);
    }
    return 0;
}

Luogu P2261 [CQOI2007]余数求和

这个题好像很水，因为好像只是一个入门的除法分块练习题。

\because k \bmod i=k-i\times\lfloor\frac{k}{i}\rfloor\\ \therefore \sum_{i=1}^nk\bmod i=\sum_{i=1}^n(k-i\times\lfloor\frac{k}{i}\rfloor)=n\times k-\sum_{i=1}^n(i\times\lfloor\frac{k}{i}\rfloor)

这样就可以做除法分块

就是对 \displaystyle {\color{blue}{\sum_{i=1}^n(i\times\lfloor\frac{k}{i}\rfloor)}} 做除法分块。

代码如下

for(ll i=1;i<=n;i=r+1)
{
    if(k/i)
        r=min((ll)n,k/(k/i));
    else 
        r=n;
    sum+=(ll)(k/i)*(ll)(r-i+1)*(i+r)/2ll;
}

Luogu P1390 公约数的和

这个题是求

\sum_{i=1}^n\sum_{j=i+1}^n\gcd(i,j)

我们先考虑求出

\sum_{i=1}^n\sum_{j=1}^n\gcd(i,j)

采用一些神必技巧，枚举一个数 k ，然后改成对 k 求和得到

\sum_{k=1}^n(k\sum_{i=1}^n\sum_{j=1}^n[\gcd(i,j)=k])\\\rightarrow\sum_{k=1}^n(k\sum_{i=1}^{\lfloor\frac{n}{k}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{k}\rfloor}[\gcd(i,j)=1])\\ \rightarrow\sum_{k=1}^n[k\sum_{i=1}^{\lfloor\frac{n}{k}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{k}\rfloor}\epsilon(\gcd(i,j))]

我们改为 \epsilon $\epsilon$ 是为了方便我们做 \texttt{Dirichlet} 卷积。

所以我们可以转化为

\sum_{k=1}^nk\sum_{i=1}^{\lfloor\frac{n}{k}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{k}\rfloor}\sum_{\mathrm d|\gcd(i,j)}\mu(\mathrm d)\\ \rightarrow \sum_{k=1}^nk\sum_{\mathrm d=1}\mu(\mathrm d)\sum_{i=1}^{\lfloor\frac{n}{k}\rfloor}[\mathrm d|i]\sum_{j=1}^{\lfloor \frac{n}{k}\rfloor}[\mathrm d|j]=\sum_{k=1}^nk\sum_{\mathrm d=1}\mu(\mathrm d)\lfloor\frac{n}{k\times\mathrm d}\rfloor\lfloor\frac{n}{k\times \mathrm d}\rfloor

不过别忘了，这个式子只是 \displaystyle \sum_{i=1}^n\sum_{j=1}^n\gcd(i,j)

我们推导得

\begin{aligned} \texttt{ans}&=\sum_{i=1}^n\sum_{j=i+1}^{n}\gcd(i,j)\\ &=\sum_{i=1}^n\sum_{j=1}^{i-1}\gcd(i,j)\\ &=\frac{\displaystyle\sum_{i=1}^n\sum_{j=1}^n\gcd(i,j)-\sum_{i=1}^n\gcd(i,i)}{2}\\ &=\frac{\displaystyle{ \color{blue}{\sum_{i=1}^n\sum_{j=1}^n\gcd(i,j)}}-\frac{(1+n)\times n}{2}}{2}\\ \end{aligned}

用之前那个式子对蓝色的 \displaystyle{\color{blue}{\sum_{i=1}^n\sum_{j=1}^n\gcd(i,j)}} 部分做数论分块求解。

代码：

#include<cstdio>
#include<algorithm>
#include<iostream>
#include<cstring>
#include<cmath>
#include<queue>
#include<map>
#include<stack>
//#include<bits/stdc++.h>

#define ll long long
#define ull unsigned long long
#define INL inline
//Tosaka Rin Suki~
using namespace std;

const int N=2000005;

INL int read()
{
     int x=0;int w=1;
     char ch=getchar();
     while((ch<'0'||ch>'9')&&ch!='-')ch=getchar();
     if(ch=='-')w=-1,ch=getchar();
     while(ch>='0'&&ch<='9')
     {x=(x<<1)+(x<<3)+ch-48,ch=getchar();}
     return x*w;
}

INL int mx(int a,int b){return a>b?a:b;}
INL int mn(int a,int b){return a<b?a:b;}

int mu[N];
ll pri[N],cnt;
bool vis[N];
int n;

void Mu()
{
    mu[1]=1;
    for(int i=2;i<=n;i++)
    {
        if(!vis[i])pri[++cnt]=i,mu[i]=-1;
        for(int j=1;j<=cnt&&i*pri[j]<=n;j++)
        {
            vis[i*pri[j]]=1;
            if(i%pri[j])mu[i*pri[j]]=-mu[i];
            else {mu[i*pri[j]]=0;break;}
        }
    }
    for(int i=1;i<=n;i++)mu[i]+=mu[i-1];
}

int main()
{
    //freopen(".in","r",stdin);
    //freopen(".out","w",stdout);
    n=read();
    ll ans=0;
    Mu();
    for(int i=1;i<=n;i++)
    {
        ll sum=0;
        for(ll l=1,r;l<=n/i;l=r+1)
        {
            r=(n/i)/((n/i)/l);
            sum+=(ll)(mu[r]-mu[l-1])*(ll)((n/i)/l)*(ll)((n/i)/l);
        }
        ans+=sum*i*1ll;
    }
    ans-=(ll)((ll)n*(ll)(n+1))>>1ll;
    ans>>=1;
    printf("%lld\n",ans);
    return 0;
}

LibreOJ #6539. 奇妙数论题

大家先回顾一下上一题，这题的部分分就是上一题。

题如其名，一道巧妙莫反题。

给定长为 n 的排列 a 求

\sum_{i=1}^n\sum_{j=1}^n\gcd(a_i,a_j)\times \gcd(i,j)

我们有部分分是 a_i=i 的，所以一直到 70 分的部分分就是求：

\sum_{i=1}^n\sum_{j=1}^n\gcd(i,j)^2

就上一题稍微改一下就有了：

\sum_{k=1}^nk^2\sum_{d=1}^{\lfloor\frac{n}{k}\rfloor}\mu(d)\lfloor\frac{n}{k\times d}\rfloor\lfloor\frac{n}{k\times d}\rfloor

正解的瓶颈其实就是排列变成无序后无从下手。

所以我们先将式子转换成：

\begin{aligned} \texttt{ans}&=\sum_{i=1}^n\sum_{j=1}^n\gcd(a_i,a_j)\times \gcd(i,j)\\ &=\sum_{d=1}^nd\sum_{i=1}^n\sum_{j=1}^n[\gcd(i,j)=d]\times \gcd(a_i,a_j)\\ &=\sum_{d=1}^nd\sum_{x=1}^{\lfloor\frac{n}{d}\rfloor}\mu(x)\sum_{i=1}^{\lfloor\frac{n}{d\times x}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{d\times x}\rfloor}\gcd(a_{idx},a_{jdx}) \end{aligned}

这里就是根据上一题的套路，\gcd(i,j) 看成 \gcd(i,j)\times 1 于是可以在数论分块的部分直接带入对应下标 a_i 的 \gcd。

但是求 \gcd(a_i,a_j) 也就是 \gcd(a_{idx},a_{jdx}) 的复杂度降不下来，考虑怎么做。

设：

f(x)=\sum_{i=1}^{\lfloor\frac{n}{x}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{x}\rfloor}\gcd(a_{ix},a_{jx})

我们观察这个函数，其实就是要求在 a 这个排列中下标分别为 x 的倍数的元素彼此之间 \gcd 的和。

我们设 \mathbb S_x 是大小为 x 的倍数的下标的集合。

于是有：

f(x)=\sum_{i\in\mathbb S_x}\sum_{j\in\mathbb S_x}\gcd(a_i,a_j)

我们进一步抽象，我们将 i\in\mathbb S_x 映射到每一个 a_i，\forall i\in \mathbb S_x, 设 a_i 的集合为 \mathbf S_x。

现在有：

f(x)=\sum_{i\in\mathbf S_x}\sum_{j\in\mathbf S_x}\gcd(i,j)\Rightarrow \sum_d d\sum_{i\in\mathbf S_x}\sum_{j\in \mathbf S_x}[\gcd(i,j)=d]

我们有另一种方法求 $[x=d]$ 这一类的式子。我们观察 $[x=d]$ 长得很像 $\epsilon(x)=[x=1]$。我们有 $\mu*1=\epsilon$ 所以有 $\displaystyle \epsilon(x)=\sum_{d|x}\mu(d)$。那么我们的 $[x=d]$ 可以化为 $\displaystyle[\frac{T}{d}=1\space\mathrm{and}\space T|d]$。类似的我们展开得 $\displaystyle [x=d]=\sum_{d|T,T|x}\mu(\frac{T}{d})

得到（这个地方看官方题解看了很久才懂的，感觉很 tricky）：

\begin{aligned} f(x)&=\sum_d d\sum_{i\in\mathbf S_x}\sum_{j\in \mathbf S_x}\sum_{d|T,T|\gcd(i,j)}\mu(\frac{T}{d})\\ &=\sum_d d\sum_{d|T}\mu(\frac{T}{d})(\sum_{i\in \mathbf S_x}[T|i])^2\\ &=\sum_T({\color{blue}{\sum_{d|T}d\times\mu(\frac{T}{d})}})(\sum_{i\in \mathbf S_x}[T|i])^2 \end{aligned}

我们知道有 \mu * \mathrm{id}=\varphi。

观察上面蓝色部分的 \displaystyle{\color{blue}{\sum_{d|T}d\times\mu(\frac{T}{d})}} 就是 \varphi(T)

最后得到：

f(x)=\sum_T \varphi(T)(\sum_{i\in \mathbf S_x}[T|i])^2

我们回到 \texttt{ans} 的部分。

有：

\texttt{ans}=\sum_{d=1}^nd\sum_{x=1}^{\lfloor\frac{n}{d}\rfloor}\mu(x)f(dx)

具体实现方法：

首先我们筛出得到 \forall x\in[1,n] 每个 \mathbf S_x。
依次统计对于每个 \mathbf S_x 每个元素的的因子的贡献。
得到 f(x) 并筛出 \mu(x) 后直接求和即可。

复杂度分析的话，看到 loj 官方题解的评论区中有这样的结论：

我们设对于 x，x 的因子个数为 d(x)。

As Vaclav Kotesovec said Aug 30 2018, \displaystyle \sum_{k\leq n}d^2(k)is asymptotic to \Theta(\frac{1}{\pi^2}n\log^3 n+n\log^2 n).

代码（其实就是出题人的写法，好像多数采用了埃筛）：

#include<bits/stdc++.h>

#define ll long long
#define ull unsigned long long
#define INL inline
#define Re register

//Tosaka Rin Suki~

INL int read()
{
    int x=0,w=1;char ch=getchar();
    while((ch<'0'||ch>'9')&&ch!='-')ch=getchar();if(ch=='-')w=-1,ch=getchar();
    while(ch>='0'&&ch<='9')x=(x<<1)+(x<<3)+ch-48,ch=getchar();return x*w;
}

const int N=1e5+5,MOD=1e9+7;

std::vector <int> factors[N];
std::vector <int> pack;
int mu[N],g[N],f[N],d[N],n,a[N],ans;

int main()
{
    //freopen(".in","r",stdin);
    //freopen(".out","w",stdout);
    n=read();for(int i=1;i<=n;i++)a[i]=read();
    mu[1]=1;
    for(int i=1;i<=n;i++)
    {
        factors[i].push_back(i);
        for(int j=i*2;j<=n;j+=i)
            mu[j]-=mu[i],factors[j].push_back(i);
    }
    for(int i=1;i<=n;i++)
        for(int j=i;j<=n;j+=i)
            g[j]=(g[j]+mu[j/i]*i)%MOD;
    for(int i=1;i<=n;i++)
    {
        for(int j=i;j<=n;j+=i)
            for(int k=0;k<factors[a[j]].size();k++)
                ++d[factors[a[j]][k]],pack.push_back(factors[a[j]][k]);
        for(int j=0;j<pack.size();j++)
            if(d[pack[j]])
                f[i]=(f[i]+1ll*d[pack[j]]*d[pack[j]]*g[pack[j]])%MOD,d[pack[j]]=0;
        pack.clear();
    }
    for(int d=1;d<=n;d++)
    {
        int s=0;
        for(int x=1;x<=n/d;x++)s=(s+1ll*mu[x]*f[d*x])%MOD;
        ans=(ans+1ll*s*d)%MOD;
    }
    printf("%d\n",ans);
    return 0;
}

Luogu P3704 [SDOI2017]数字表格

~~在机房复盘 SDOI 2017 的模拟赛第一次场 A 黑题写题解纪念一下。~~

这个题要你求：

\prod_{i=1}^n\prod_{j=1}^m \mathrm{fib_{\gcd(i,j)}}

其中 \mathrm{fib_{i}} 表示斐波拉契数列第 i 项。

首先一般性地，我们保证 ${\color{red}{n\leq m}}$。我们枚举 $\displaystyle\sum_i\sum_j f(i,j)$ 这一类式子的时候会习惯性地枚举一个 $k$ 让其满足和式内部的条件。这题即使是 $\prod$ 也可以套路化地使用这个方法。所以原式可以转换为 $$ \prod_{k=1}^n\prod_{i=1}^n\prod_{j=1}^m \mathrm{fib}_k\times[\gcd(i,j)=k] $$ 我们发现每个 $\mathrm{fib}_d$ 都会与自己相乘很多次，于是我们可以将每个 $d$ 的贡献考虑写成幂的形式。显然每个 $\mathrm{fib}_d$ 的幂的指数就是其等于 $\gcd(i,j)$ 对应的数对 $(i,j)$ 的个数。现在我们有 $$ \prod_{k=1}^n\mathrm {fib}_k^{\displaystyle\bigg({\color{blue}{\sum_{i=1}^n\sum_{j=1}^m[\gcd(i,j)=k]}}\bigg)} $$ 我们把指数特地标记出来。因为我们发现他是一个老朋友。我们在做形似 $\displaystyle\sum_{i=1}^n\sum_{j=1}^n\gcd(i,j)$ 的一类莫反题是会将式子转换成 $\displaystyle \sum_{k=1}^n\bigg(k{\color{blue}{\sum_{i=1}^n\sum_{j=1}^n[\gcd(i,j)=k]}}\bigg)

我们考虑怎么将这个式子转化成可做的形式（以下皆是莫反题老套路，跟不上请先去做一些基础莫反题~~，不是幽灵乐团，别信 CYJian~~）。

\begin{aligned} {\color{blue}{\texttt{blue part}}}&=\sum_{i=1}^{\lfloor\frac{n}{k}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{k}\rfloor}[\gcd(i,j)=1]\\ &=\sum_{i=1}^{\lfloor\frac{n}{k}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{k}\rfloor}\epsilon(\gcd(i,j))\\ &=\sum_{i=1}^{\lfloor\frac{n}{k}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{k}\rfloor}\sum_{d|\gcd(i,j)}\mu(d)\\ &=\sum_{d=1}^{\lfloor\frac{n}{k}\rfloor}\mu(d)\sum_{i=1}^{\lfloor\frac{n}{k}\rfloor}[ d|i]\sum_{j=1}^{\lfloor \frac{m}{k}\rfloor}[d|j]\\ &=\sum_{d=1}^{\lfloor\frac{n}{k}\rfloor}\mu(d)\lfloor\frac{n}{k\times d}\rfloor\lfloor\frac{m}{k\times d}\rfloor \end{aligned}

这个东西已经可以 \Theta(\sqrt n) 级别的数论分块做了，但是我们要把其代入原式继续看。

现在原式变成了：

\prod_{k=1}^n\mathrm {fib}_k^{\bigg (\displaystyle\sum_{d=1}^{\lfloor\frac{n}{k}\rfloor}\mu(d)\lfloor\frac{n}{k\times d}\rfloor\lfloor\frac{m}{k\times d}\rfloor\bigg)}

我们设 T=k\times d。

我们发现这样方便让我们将 \displaystyle\sum_{d=1}^{\lfloor \frac{n}{k}\rfloor} 与 \mathrm{fib}_k 以幂的形式结合。所有变量用 T 表示出来后贡献到一个 \prod 中。

写出来就是：

\prod_{T=1}^n\bigg({\color{blue}{\prod_{k|T}\mathrm{fib}_k^{\mu(\frac{T}{k})}}}\bigg)^{\lfloor\frac{n}{T}\rfloor\lfloor\frac{m}{T}\rfloor}

蓝色部分我们可以通过类似筛法预处理的方法求出：

线筛 \mu 函数。
对于每一个 k 枚举其倍数 T ，乘上这个 \mathrm{fib}_k 算入这个 T 的贡献。
因为 \mu 只有三种取值，有影响的只有 \{1,-1\} 两种。
所以问题就变成了乘上 \mathrm {fib}_{k} 或除以 \mathrm {fib}_k。这里就变成 \mathrm{inv fib}（逆元）和 \mathrm{fib} 两种情况。

外层我们做数论分块的时候，每次求内部只需要一次 \Theta (\log p)（p 是模数）和 \Theta(\log(\lfloor\frac{n}{T}\rfloor\lfloor\frac{m}{T}\rfloor))的快速幂分别求逆元和与外层指数的幂。

预处理复杂度是 \Theta(n\log n+n\log p)，数论分块复杂度是 \Theta(t\sqrt n (\log p+\log \lfloor\frac{n}{T}\rfloor\lfloor\frac{m}{T}\rfloor)) 小写 t 是数据组数。

实现看代码（考试时忘了写逆元卡了好久导致 T2 LCT 都没写完）：

#include<bits/stdc++.h>

#define ll long long
#define ull unsigned long long
#define INL inline
#define Re register

//Tosaka Rin Suki~

INL int read()
{
    int x=0,w=1;char ch=getchar();
    while((ch<'0'||ch>'9')&&ch!='-')ch=getchar();if(ch=='-')w=-1,ch=getchar();
    while(ch>='0'&&ch<='9')x=(x<<1)+(x<<3)+ch-48,ch=getchar();return x*w;
}

const int N=1e6+5;
const ll MOD=1e9+7;

ll f[N],finv[N],prime[N],mu[N],cnt,prod[N];
bool vis[N];

INL ll fpow(ll x,ll p)
{
    ll res=1;
    while(p)
    {
        if(p&1)res=(1ll*res*x)%MOD;
        x=(1ll*x*x)%MOD;p>>=1;
    }
    return res;
}

INL void get_mu(ll n)
{
    f[0]=0,f[1]=1;
    for(int i=2;i<=n;i++)f[i]=(f[i-1]+f[i-2])%MOD;
    for(int i=0;i<=n;i++)finv[i]=fpow(f[i],MOD-2),prod[i]=1;
    vis[1]=1,vis[0]=1,mu[1]=1;
    for(int i=2;i<=n;i++)
    {
        if(!vis[i])prime[++cnt]=i,mu[i]=-1;
        for(int j=1;j<=cnt&&prime[j]*i<=n;j++)
        {
            vis[prime[j]*i]=1;
            if(i%prime[j])mu[i*prime[j]]=-mu[i];
            else {mu[i*prime[j]]=0;break;}
        }
    }
    for(int i=1;i<=n;i++)
    {
        for(int j=i;j<=n;j+=i)
        {
            int now=j/i;
            if(mu[now]==-1)prod[j]=1ll*prod[j]*finv[i]%MOD;
            else if(mu[now]==1)prod[j]=1ll*prod[j]*f[i]%MOD;
        }
    }
    for(int i=1;i<=n;i++)prod[i]=1ll*prod[i-1]*prod[i]%MOD;
}

int main()
{
    //freopen("product.in","r",stdin);
    //freopen("product.out","w",stdout);
    int T=read();
    get_mu(1000000);
    while(T--)
    {
        ll n=read(),m=read();
        if(n<m)std::swap(n,m);
        ll l=1,r,ans=1;
        while(l<=m)
        {
            r=std::min(n/(n/l),m/(m/l));
            ans=1ll*ans*fpow(1ll*prod[r]*fpow(prod[l-1],MOD-2)%MOD,1ll*(n/l)*(m/l)%(MOD-1))%MOD;
            l=r+1;
        }
        printf("%lld\n",(ans%MOD+MOD)%MOD);
    }
    return 0;
}