danihao123's Blog

求有多少正整数序列，满足所有数的最大公约数为\(x\)，所有数的和为\(y\)。

\(1\leq x, y\leq 10^9\)。

沿用约数个数和一题的思路……可以列出此式：

\[\sigma_1(ij) = \sum_{a | i}\sum_{b | j} ab[\gcd(\frac{i}{a}, b) = 1]\]

然后我们考虑去化简原式：

\[
\begin{aligned}
\sum_{i = 1}^n\sum_{j = 1}^n\sum_{a | i}\sum_{b | j} ab\sum_{d | \frac{i}{a}, d | j}\mu(d)
\end{aligned}
\]

接下来我们考虑枚举\(d\)，但是\(a\)和\(b\)的贡献又有点麻烦了……

\(b\)的贡献还好说，直接枚举\(b\)本身是\(d\)的多少倍就是\(\sum_{b = 1}^{\lfloor\frac{n}{d}\rfloor} bd\lfloor\frac{n}{bd}\rfloor\)。那\(a\)的贡献如何考虑？

我们考虑直接去枚举\(a\)本身……可以注意到一定有\(a\le\lfloor\frac{n}{d}\rfloor\)（因为\(\frac{i}{a}\ge d\)），所以直接枚举\(a\)本身之后再考虑\(\lfloor\frac{n}{a}\rfloor\)范围内\(d\)的倍数的数量即可（相当于找\(\frac{i}{a}\)），因此贡献为\(\sum_{a = 1}^{\lfloor\frac{n}{d}\rfloor} a\lfloor\frac{n}{ad}\rfloor\)。

那么继续化简原式：

\[
\begin{aligned}
\quad&\sum_{d = 1}^n\mu(d)\sum_{b = 1}^{\lfloor\frac{n}{d}\rfloor} bd\lfloor\frac{n}{bd}\rfloor\sum_{a = 1}^{\lfloor\frac{n}{d}\rfloor} a\lfloor\frac{n}{ad}\rfloor\\
=&\sum_{d = 1}^n\mu(d)d(\sum_{i = 1}^{\lfloor\frac{n}{d}\rfloor} i\lfloor\frac{n}{id}\rfloor)^2\\
=&\sum_{d = 1}^n\mu(d)d S_1^2(\lfloor\frac{n}{d}\rfloor)
\end{aligned}
\]

此处\(S_1\)表示\(\sigma_1\)的前缀和。

接下来预处理\(\mu(d)d\)的前缀和，考虑杜教筛。我们发现该函数和\(\mathrm{id}\)卷出来就是\(\epsilon\)（证明可以考虑贝尔级数……这种方式非常有效，并且还能帮我们找到需要卷的函数，我有时间会专门撰文写一下），所以杜教筛一波即可。这部分总复杂度为\(O(n^{\frac{2}{3}})\)。

至于\(S_1\)，我们考虑不大于\(n^{\frac{2}{3}}\)可以直接预处理，剩下的用时就用\(O(\sqrt{n})\)的方法求（考虑到反演不会用到\(S_1\)的重复状态所以不需要记忆化）。用类似于杜教筛复杂度证明的方法可以证明该部分总复杂度为\(O(n^{\frac{2}{3}})\)。

代码：

#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cmath>
#include <algorithm>
#include <functional>
#include <utility>
#include <unordered_map>
const int N = 1000000;
using ll = long long;
const ll ha = 1000000007LL;
const ll i2 = 500000004LL;
ll mud[N + 5], sig_1[N + 5];
int prm[N + 5]; bool vis[N + 5];
void sieve() {
  int cnt = 0;
  mud[1] = 1;
  for(int i = 2; i <= N; i ++) {
    if(!vis[i]) {
      mud[i] = (ha - i);
      prm[cnt ++] = i;
    }
    for(int j = 0; j < cnt; j ++) {
      int v = i * prm[j];
      if(v > N) break;
      vis[v] = true;
      if(i % prm[j] == 0) {
        mud[v] = 0; break;
      } else {
        mud[v] = (mud[i] * mud[prm[j]]) % ha;
      }
    }
  }
  for(int i = 1; i <= N; i ++) {
    for(int j = i; j <= N; j += i) {
      sig_1[j] = (sig_1[j] + (ll(i))) % ha;
    }
  }
  for(int i = 1; i <= N; i ++) {
    mud[i] = (mud[i - 1] + mud[i]) % ha;
    sig_1[i] = (sig_1[i - 1] + sig_1[i]) % ha;
  }
}

std::unordered_map<ll, ll> ma;
ll calc_id(ll x) {
  ll v1 = x, v2 = x + 1LL;
  if(x & 1LL) {
    v2 /= 2LL;
  } else {
    v1 /= 2LL;
  }
  return ((v1 * v2) % ha);
}
ll calc_mud(ll n) {
  if(n <= (ll(N))) return mud[n];
  if(ma.count(n)) return ma[n];
  ll ans = 1, las = 1;
  for(ll i = 2; i <= n;) {
    ll next = n / (n / i);
    ll nv = calc_id(next);
    ll ns = (nv - las + ha) % ha;
    ans = (ans - (ns * calc_mud(n / i)) % ha + ha) % ha;
    las = nv; i = next + 1LL;
  }
  ma[n] = ans; return ans;
}
ll calc_s1(ll n) {
  if(n <= (ll(N))) return sig_1[n];
  ll ans = 0, las = 0;
  for(ll i = 1; i <= n;) {
    ll next = n / (n / i);
    ll nv = calc_id(next);
    ll ns = (nv - las + ha) % ha;
    ans = (ans + (ns * (n / i)) % ha) % ha;
    las = nv; i = next + 1LL;
  }
  return ans;
}
ll calc(ll n) {
  ll ans = 0, las = 0;
  for(ll i = 1; i <= n;) {
    ll next = n / (n / i);
    ll nv = calc_mud(next); ll ns = (nv - las + ha) % ha;
    ll pv = calc_s1(n / i); pv = (pv * pv) % ha;
    ans = (ans + (ns * pv) % ha) % ha;
    las = nv; i = next + 1LL;
  }
  return ans;
}

int main() {
  sieve();
  ll n; scanf("%lld", &n);
  printf("%lld\n", calc(n));
  return 0;
}

首先我们有这样一个式子：

\[d(ij) = \sum_{a | i}\sum_{b | j} [\gcd(\frac{i}{a}, b) = 1]\]

怎样理解这一结论呢？我们知道\(ij\)的约数\(d\)一定可以被分解为\(ab\)使得\(a | i, b | j\)，但这种划分可能会非常的多，我们要保证每个约数\(d\)只有一种划分被统计过。

那么考虑上述划分方法。对于\(i\)和\(j\)中一个有一个没有的质因子，是一定会划分到\(a\)和\(b\)中固定的一个的。那么考虑一个在\(i\)和\(j\)中都出现的质因子\(p\)，假设\(p\)在\(d\)中出现的次数超过了\(i\)和\(j\)中任何一个的话，那么\(p\)在\(a\)中一定会尽可能多的出现（否则\(\frac{i}{a}\)中会出现\(p\)，而这样无论如何一定有\(p | b\)，因此会不满足\(\gcd(\frac{i}{a}, b) = 1\)）。如果没超过呢？那么一定会全部划分到\(a\)中（其他情况下有\(p | b\)并且会出现\(p | \frac{i}{a}\)）。

下面考虑应用这个结论（为了方便，使用等价的\(d(ij) = \sum_{a | i}\sum_{b | j} [\gcd(a, b) = 1]\)，假设\(n\le m\)）：

\[
\begin{aligned}
\quad&\sum_{i = 1}^n\sum_{j = 1}^m \sum_{a | i}\sum_{b | j}\sum_{d | i, d | j}\mu(d)\\
=&\sum_{d = 1}^n\mu(d)\sum_{a = 1}^{\lfloor\frac{n}{d}\rfloor}\lfloor\frac{n}{ad}\rfloor\sum_{b = 1}^{\lfloor\frac{m}{d}\rfloor}\lfloor\frac{m}{bd}\rfloor\\
=&\sum_{d = 1}^n\mu(d)\sum_{a = 1}^{\lfloor\frac{n}{d}\rfloor}d(a)\sum_{b = 1}^{\lfloor\frac{m}{d}\rfloor}d(b)
\end{aligned}
\]

后面两个和式可以随便\(O(n)\)欲处理一下（然后我当初写了个\(O(n\sqrt{n})\)的神必预处理……），然后这个题做完力，，，

代码：

#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cctype>
#include <algorithm>
#include <queue>
#include <utility>
typedef long long int_t;
const int maxn = 50005;
int_t miu[maxn], miu_s[maxn], f[maxn];
void process() {
  static bool vis[maxn];
  static int prm[maxn]; int cnt = 0;
  const int N = 50000;
  miu[1] = 1LL;
  for(int i = 2; i <= N; i ++) {
    if(!vis[i]) {
      miu[i] = -1LL; prm[cnt ++] = i;
    }
    for(int j = 0; j < cnt && prm[j] * i <= N; j ++) {
      int v = prm[j] * i;
      vis[v] = true;
      if(i % prm[j] == 0) {
        miu[v] = 0; break;
      } else {
        miu[v] = miu[i] * -1;
      }
    }
  }
  for(int i = 1; i <= N; i ++) {
    miu_s[i] = miu_s[i - 1] + miu[i];
  }
  for(int i = 1; i <= N; i ++) {
    for(int j = 1; j <= i;) {
      int nx = i / (i / j);
      f[i] += (int_t(nx - j + 1)) * (int_t(i / j));
      j = nx + 1;
    }
  }
}

int_t solve(int n, int m) {
  int_t ans = 0;
  for(int i = 1; i <= std::min(n, m);) {
    int nx = std::min(m / (m / i), n / (n / i));
    ans += (miu_s[nx] - miu_s[i - 1]) * f[m / i] * f[n / i];
    i = nx + 1;
  }
  return ans;
}
#ifdef LOCAL
#define LO "%I64d"
#else
#define LO "%lld"
#endif
int main() {
  process(); int T;
  scanf("%d", &T);
  while(T --) {
    int n, m; scanf("%d%d", &n, &m);
    printf(LO, solve(n, m)); puts("");
  }
  return 0;
}

首先这个题和能量采集几乎完全一致……稍微有点反演常识的人都知道答案事\(\sum_{i = 1}^n\varphi(i)\lfloor\frac{n}{i}\rfloor^2\)，然后你需要用\(\varphi\)的前缀和，这个东西直接杜教筛搞一波就行了。这样总复杂度事\(O(n^{\frac{2}{3}})\)的，下面给出证明。

先考虑杜教筛的部分。杜教筛可以看做一个转移复杂度为\(O(\sqrt{n})\)的DP，他的状态一定是通过总的\(n\)整除一个数得到的。不大于\(n^{\frac{2}{3}}\)的状态我们都事先预处理，那么我们假设所有大于\(n^{\frac{2}{3}}\)的状态都被计算了一遍，这些状态一定事通过\(n\)整除一个不大于\(n^{\frac{1}{3}}\)的数得到的。因此复杂度为：

\[\sum_{i = 1}^{n^{\frac{1}{3}}} \sqrt{\lfloor\frac{n}{i}\rfloor}\leq \int_0^{n^{\frac{1}{3}}} \sqrt{\frac{n}{x}}\mathrm{d}x = 2n^{\frac{2}{3}} = O(n^{\frac{2}{3}})\]

至于莫比乌斯反演，不考虑杜教筛的复杂度的话就只是\(O(\sqrt{n})\)而已了。

代码：

#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cctype>
#include <cmath>
#include <algorithm>
#include <functional>
#include <unordered_map>
#include <utility>
using ll = long long;
constexpr ll ha = 1000000007LL;
ll pow_mod(ll a, ll b) {
  ll ans = 1, res = a;
  while(b) {
    if(1LL & b) ans = (ans * res) % ha;
    res = (res * res) % ha; b >>= 1;
  }
  return ans;
}
ll inv(ll x) {
  return pow_mod(x, ha - 2LL);
}

const int N = 10000000;
ll phi[N + 1]; int prm[N + 1]; bool vis[N + 1];
void sieve() {
  phi[1] = 1; vis[1] = true;
  int cnt = 0;
  for(int i = 2; i <= N; i ++) {
    if(!vis[i]) {
      phi[i] = i - 1; prm[cnt ++] = i;
    }
    for(int j = 0; j < cnt; j ++) {
      int v = i * prm[j]; if(v > N) break;
      vis[v] = true;
      if(i % prm[j] == 0) {
        phi[v] = phi[i] * prm[j]; break;
      } else {
        phi[v] = phi[i] * phi[prm[j]];
      }
    }
  }
  for(int i = 1; i <= N; i ++) {
    phi[i] = (phi[i] + phi[i - 1]) % ha;
  }
}

std::unordered_map<ll, ll> ma;
const ll i2 = inv(2LL);
ll phi_sum(ll n) {
  if(n <= (ll(N))) return phi[n];
  if(ma.count(n)) return ma[n];
  ll fg = (n % ha) * ((n + 1LL) % ha) % ha * i2 % ha;
  for(ll i = 2; i <= n;) {
    ll next = n / (n / i);
    ll len = (next - i + 1LL) % ha;
    fg = (fg - (len * phi_sum(n / i)) % ha + ha) % ha;
    i = next + 1LL;
  }
  ma[n] = fg; return fg;
}
ll calc(ll n) {
  ll ans = 0, las = 0;
  for(ll i = 1; i <= n;) {
    ll next = n / (n / i);
    ll b = ((n / i) % ha) * ((n / i) % ha) % ha;
    ll ns = phi_sum(next);
    ans = (ans + (((ns - las + ha) % ha) * b % ha)) % ha; las = ns;
    i = next + 1LL;
  }
  return ans;
}

int main() {
  sieve();
  ll n; scanf("%lld", &n);
  printf("%lld\n", calc(n));
  return 0;
}

卡常的题见过，这么变态的卡常题……可能出题人过于文明，，，下面是卡常记录的一部分：

下面开始颓式子：

本来想做JZPKIL的……但卡常太恶心了

上来先颓式子：

\[
\DeclareMathOperator{\id}{id}
\begin{aligned}
f_d(n) &= \sum_{i = 1}^n [(i, n) = 1]i^d\\
&= \sum_{i = 1}^n i^d\sum_{k | n, k | i}\mu(k)\\
&= \sum_{k | n}\mu(k)\sum_{i = 1}^{\frac{n}{k}} (ik)^d\\
&= \sum_{k | n}\mu(k)k^d h_d(\frac{n}{d})\\
&= ((\mu\cdot\id^k)\ast h_d)
\end{aligned}
\]

其中\(h_d\)表示指数为\(d\)时的等幂求和函数。然后我们可能会想，如果这个函数是积性函数，那么我们可以对质因数分解的每一种质因数都求一遍答案，而由于\(\mu\)的存在（对于质数的幂\(p^k\)的因数，只有\(1\)和\(p\)的莫比乌斯函数值不为0），所以需要处理的情况只有两种。

不过很可惜，\(h_d\)显然不是积性函数，所以最后的整个卷积也不是积性函数。但是我们注意到\(h_d\)事一个\(d + 1\)次多项式，所以我们可以把多项式每一项都当成单独的函数，然后每个单独函数卷出来都是一个积性函数。至于求多项式每一项的系数，用伯努利数即可。

代码：

/**************************************************************
    Problem: 3601
    User: danihao123
    Language: C++
    Result: Accepted
    Time:220 ms
    Memory:948 kb
****************************************************************/
 
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <algorithm>
#include <utility>
typedef long long ll;
ll ha = 1000000007LL;
ll pow_mod(ll a, ll b) {
  ll ans = 1LL, res = a;
  while(b) {
    if(1LL & b) ans = (ans * res) % ha;
    res = (res * res) % ha;
    b >>= 1;
  }
  return ans;
}
ll inv(ll x) {
  return pow_mod(x, ha - 2LL);
}
 
const int maxn = 105;
ll C[maxn][maxn];
void process_C() {
  C[0][0] = 1LL;
  for(int i = 1; i < maxn; i ++) {
    C[i][0] = C[i][i] = 1LL;
    for(int j = 1; j < i; j ++) {
      C[i][j] = (C[i - 1][j] + C[i - 1][j - 1]) % ha;
    }
  }
}
ll B[maxn];
void process_B() {
  B[0] = 1LL;
  for(int i = 1; i <= 101; i ++) {
    B[i] = 1LL;
    for(int j = 0; j < i; j ++) {
      ll d = (((C[i][j] * B[j]) % ha) * inv(i - j + 1)) % ha;
      B[i] = (B[i] - d + ha) % ha;
    }
  }
}
 
const int maxw = 1005;
typedef std::pair<ll, ll> pii;
pii P[maxw]; ll dv[maxw][3];
int d, w;
void process_dv() {
  for(int i = 1; i <= w; i ++) {
    ll p = P[i].first, t = P[i].second;
    p %= ha;
    dv[i][0] = pow_mod(p, t);
    dv[i][1] = pow_mod(p, t - 1LL);
    dv[i][2] = pow_mod(p, d);
  }
}
ll calc(ll k, int z) {
  ll ans = k;
  for(int i = 1; i <= w; i ++) {
    ll p = P[i].first, t = P[i].second;
    p %= ha;
    ll f1 = pow_mod(dv[i][0], z), f2 = pow_mod(dv[i][1], z);
    ll v1 = (f1 - (dv[i][2] * f2) % ha + ha) % ha;
    ans = (ans * v1) % ha;
  }
#ifdef LOCAL
  printf("Case (%lld, %d) : %lld\n", k, z, ans);
#endif
  return ans;
}
 
int main() {
  process_C(); process_B();
  scanf("%d%d", &d, &w);
  for(int i = 1; i <= w; i ++) {
    scanf("%lld%lld", &P[i].first, &P[i].second);
  }
  process_dv();
  ll ans = 0LL;
  for(int i = 0; i <= d; i ++) {
    ll g = calc((((C[d + 1][i] * B[i]) % ha) * inv(d + 1)) % ha, d + 1 - i);
    ans = (ans + g) % ha;
  }
  printf("%lld\n", ans);
  return 0;
}

一年前看的反演神题……（然后现在才A）

MIN-MAX容斥第一题……

众所周知斐波那契数有个非常好的性质：

\[\gcd(\{f_i | i\in T\}) = f_{\gcd\{T\}}\]

但很不幸的事情是，lcm没有类似的性质。那我们就考虑怎么把lcm改成gcd。

根据MIN-MAX容斥，可以得到这个式子（可以理解为对指数的容斥。这里设\(S\)为全集）：

\[\prod_{T\subseteq S, T\neq\emptyset} f_{\gcd\{T\}}^{(-1)^{|T| + 1}}\]

gcd套在上面非常烦人，考虑用狄利克雷卷积去化解它。考虑有个函数\(g\)满足\(f_n = \prod_{d|n} g(d)\)，然后xjb化简一下之后柿子变成：

\[\prod_{d} g(d)^{\sum_{T\subseteq S, T\neq\emptyset, \forall x\in T, d | x}(-1)^{|T| + 1}}\]

再去考虑那个指数。当且仅当\(d\in S\)时，\(g(d)\)才会对答案做贡献。并且无论\(S\)中\(d\)的约数有多少个（只要大于），最后指数肯定都是1。因为其实大于零的情况那个指数也是一个所有最小值都是1的MIN-MAX容斥……这样求出来最大值肯定是1了，那个指数自然就是1。

于是乎柿子变成了：

\[\prod_{\exists x\in S, d | x} g(d)\]

然后搞出来\(S\)中所有数的约数很容易。当务之急就是求\(g\)了（这里其实就是莫比乌斯反演了）。不过球莫比乌斯函数反而让问题复杂化了……移下项可以发现：

\[g(n) = f_n\cdot (\prod_{d | n, d\neq n} g(d))^{-1}\]

然后就很愉快的搞出来了（逃

代码：

#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cmath>
#include <cctype>
#include <algorithm>
#include <utility>
#include <queue>
const int maxn = 1000005;
using ll = long long;
const ll ha = 1000000007LL;
ll pow_mod(ll a, ll b) {
  ll ans = 1LL, res = a;
  while(b) {
    if(1LL & b) {
      ans = (ans * res) % ha;
    }
    res = (res * res) % ha;
    b /= 2LL;
  }
  return ans;
}
ll inv(ll x) {
  return pow_mod(x, ha - 2LL);
}

ll f[maxn], g[maxn];
void process_f() {
  int N = 1000000;
  f[0] = 0LL; f[1] = 1LL;
  for(int i = 2; i <= N; i ++) {
    f[i] = (f[i - 1] + f[i - 2]) % ha;
  }
  std::fill(g + 1, g + 1 + N, 1LL);
  for(int i = 1; i <= N; i ++) {
    g[i] = (f[i] * inv(g[i])) % ha;
    for(int j = i * 2; j <= N; j += i) {
      g[j] = (g[j] * g[i]) % ha;
    }
  }
#ifdef LOCAL
  puts("g has been processed!");
#endif
}

bool vis[maxn];
void process_vis() {
  int N = 1000000;
  for(int i = 1; i <= N; i ++) {
    for(int j = i * 2; j <= N; j += i) {
      vis[i] = (vis[i] || vis[j]);
    }
  }
}

int main() {
  process_f();
  int n; scanf("%d", &n);
  for(int i = 1; i <= n; i ++) {
    int a; scanf("%d", &a);
    vis[a] = true;
  }
  process_vis();
  ll ans = 1LL;
  for(int i = 1; i <= 1000000; i ++) {
    if(vis[i]) {
      ans = (ans * g[i]) % ha;
    }
  }
  printf("%lld\n", ans);
  return 0;
}

当年在考场上一点反演都不会啊啊啊啊啊啊

这题虽然牵扯到了对积的反演，但其实也不是很难。先让我们大力推导一波（逃

考虑枚举柿子中的\(f[(i, j)]\)中的最大公约数（设为\(k\)），然后式子变成了（这里不妨偷个懒，钦定\(n\leq m\)）：

\[\prod_{k = 1}^n f[k]^{\sum_{i = 1}^n \sum_{j = 1}^m\:[(i, j) = k]}\]

然后发现上面那个指数就是Problem b的式子，显然可化为\(\sum_{d = 1}^n \mu(d) \lfloor\frac{n}{kd}\rfloor \lfloor\frac{m}{kd}\rfloor\)。然后似乎化简没有余地了，但是根据直觉，我们可以钦定\(T = kd\)，然后枚举\(T\)。

然后柿子变成了：

\[\prod_{T = 1}^n (\prod_{k\mid T} f[k]^{\mu(\tfrac{T}{k})})^{\lfloor\tfrac{n}{T}\rfloor \lfloor\tfrac{m}{T}\rfloor}\]

柿子中的\(\prod_{k\mid T} f[k]^{\mu(\tfrac{T}{k})}\)是一个类似狄利克雷卷积的东西（取完对数就是了），可以枚举约数预处理。并且这个东西很不错的一点就是上面的指数是莫比乌斯函数，可以取到的值只有那三种，不需要快速幂。

然后剩下的部分就和通常的反演题一般无二了……

代码（卡线过的蒟蒻……求轻喷……但我在LOJ上跑的还挺快的）：

/**************************************************************
    Problem: 4816
    User: danihao123
    Language: C++
    Result: Accepted
    Time:50840 ms
    Memory:33048 kb
****************************************************************/
 
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <algorithm>
#include <utility>
#include <cctype>
#include <cassert>
const int N = 1000000;
const int maxn = N + 5;
typedef long long ll;
const ll ha = 1000000007LL;
bool p_flag = false;
ll pow_mod(ll a, ll b) {
  if(!b) return 1LL;
  ll ans = pow_mod(a, b / 2LL);
  ans = (ans * ans) % ha;
  if(1LL & b) ans = (ans * a) % ha;
#ifdef LOCAL
  if(p_flag) printf("%lld^%lld : %lld\n", a, b, ans);
  if(b == ha - 2LL) p_flag = false;
#endif
  return ans;
}
ll inv(ll x) {
  // if(x == 1LL) p_flag = true;
  return pow_mod(x, ha - 2LL);
}
int miu[maxn];
void sieve() {
  static int prm[maxn]; int cnt = 0;
  static bool vis[maxn];
  miu[1] = 1;
  for(int i = 2; i <= N; i ++) {
    if(!vis[i]) {
      miu[i] = -1;
      prm[cnt ++] = i;
    }
    for(int j = 0; j < cnt && prm[j] * i <= N; j ++) {
      int v = prm[j] * i;
      vis[v] = true;
      if(i % prm[j] == 0) {
        miu[v] = 0;
        break;
      } else {
        miu[v] = miu[i] * -1;
      }
    }
  }
}
ll f[maxn], inv_d[maxn], F[maxn];
void process() {
  sieve();
  f[1] = 1LL; inv_d[1] = 1LL;
  for(int i = 2; i <= N; i ++) {
    f[i] = (f[i - 1] + f[i - 2]) % ha;
    inv_d[i] = inv(f[i]);
    assert(f[i] != 0LL); assert(inv_d[i] != 0LL);
  }
  for(int i = 1; i <= N; i ++) F[i] = 1LL;
  for(int i = 1; i <= N; i ++) {
    for(int j = i; j <= N; j += i) {
      int res = j / i;
      if(miu[res] == 1) {
        F[j] = (F[j] * f[i]) % ha;
#ifdef LOCAL
        if(j <= 3) printf("f[%d](%lld) -> F[%d]\n", i, f[i], j);
#endif
        continue;
      }
      if(miu[res] == -1) {
        F[j] = (F[j] * inv_d[i]) % ha;
#ifdef LOCAL
        if(j <= 3) printf("inv_f[%d](%lld) -> F[%d]\n", i, inv_d[i], j);
#endif
      }
    }
  }
  F[0] = 1LL;
  bool flag = true;
  for(int i = 1; i <= N; i ++) {
    F[i] = (F[i] * F[i - 1]) % ha;
    if(flag && F[i] == 0LL) {
      printf("SF[%d] has been 0.\n", i);
      flag = false;
    }
  }
}
 
ll calc(ll n, ll m) {
  if(n > m) std::swap(n, m);
  ll ans = 1LL;
  for(ll i = 1; i <= n;) {
    ll nx = std::min(n / (n / i), m / (m / i));
#ifdef LOCAL
    // printf("nx of %lld : %lld\n", i, nx);
#endif
    ll mulv = pow_mod((F[nx] * inv(F[i - 1])) % ha, (n / i) * (m / i));
    ans = (ans * mulv) % ha;
    i = nx + 1LL;
  }
  return ans;
}
 
int main() {
  process();
  int T; scanf("%d", &T);
  while(T --) {
    int n, m; scanf("%d%d", &n, &m);
    printf("%lld\n", calc(n, m));
  }
  return 0;
}

开个坑，记录一些和反演以及狄利克雷卷积的东西。

首先积性函数、狄利克雷卷积等基本概念不写了，就讲讲性质吧。

有几个一定要记住的东西：

\[\mu \ast 1 = e\]

\[\varphi \ast 1 = id\]

\[\mu \ast id = \varphi\]

这几个在推式子的过程中都有很大的作用，务必要记住。

所谓莫比乌斯反演，其实就是：

\[F = f \ast 1\Leftrightarrow f = F \ast \mu\]

（谜之音：其实很多所谓“反演题”都没用到这俩性质啊……）

关于莫比乌斯函数本身，还有一个好康的性质：

\[(\mu\ast 1)(k) = \sum_{i = 0}^k (-1)^i C_k^i\]