loj #509. 「LibreOJ NOI Round #1」动态几何问题「反演」

#509. 「LibreOJ NOI Round #1」动态几何问题

题意：

给了一个初中几何题，使用三角函数、相似、勾股定理等技巧，最后变成了：

给定 (n,m leq 1.5 imes 10^{16})，求多少对有序二元组 ((x,y)) 满足 (xin[1,n],yin[1,m],xy) 是完全平方数。

题解：

官方题解

省略亿点点内容，直接跳到题解的算法 5

算法 5

不妨设 (nleq m)。

容易发现令 (x = da,y=d'b)，其中 (a,b) 分别是 (x,y) 最大的平方因子，也就意味着 (d,d') 是形如 (prod p_i) 的数。由于 (xy=dd'ab) 是完全平方数，能推出 (d=d')。于是我们想到可以枚举无平方因子数 (d)，然后计算多少个 (x,y) 是 (d) 的倍数，且 (frac{x}{d},frac{y}{d}) 是完全平方数。答案为：

[sum_{d = 1}^{n} mu^2(i) leftlfloor sqrt{leftlfloordfrac{n}{i} ight floor} ight floor leftlfloor sqrt{leftlfloordfrac{m}{i} ight floor} ight floor ]

稍微解释一下，(mu^2(i)) 表示的是 (i) 是否含平方因子；(1) 到 (x) 中的完全平方数有 (lfloorsqrt x floor) 个。

有一个重要的式子不得不提，即 (mu^2) 的前缀和：

[sum_{i = 1}^n mu^2(i) = sum_{i = 1}^{lfloorsqrt n floor} mu(i)leftlfloor {dfrac{n}{i^2}} ight floor ]

证明其实就是容斥，等价于证明 (mu^2(n) = sum_{d^2 mid n} mu(d))，若 (n) 有 (k) 个素因子的次数至少为 (2)，(n) 被计算到的次数为 (sum_{i = 0}^k {kchoose i}(-1)^k=[k=0])。

于是我们可以 (mathcal O(sqrt n)) 计算 (mu^2) 的前缀和了。

注意到 (leftlfloor sqrt{dfrac{n}{i}} ight floor,leftlfloor sqrt{dfrac{m}{i}} ight floor) 分别有 (mathcal O( n^{frac{1}{3}})) 和 (mathcal O(min(n, m^{frac{1}{3}}))) 种取值，所以我们可以分段求：对于 (ileq m^{frac{1}{3}}) 求，再对于 (leftlfloor sqrt{dfrac{m}{i}} ight floor leq m^{frac{1}{3}}) 求。时间复杂度是多少？我们可以通过积分求，最后大概得到 (n = m) 时是 (mathcal O(n^{0.5} log n))。

算法 6

如果我们处理了 (sqrt n) 以内 (mu,mu^2) 的前缀和，那么 (mu) 的前缀和 (S(n)) 就可以 (mathcal O(n^{frac{1}{3}})) 求了。

用积分精细分析一下复杂度大约是 (mathcal O(sqrt n))。

算法 7

如果我们预处理 (S = min(n,m^{frac{3}{7}})) 内的 (mu,mu^2) 前缀和，使用杜教筛求 (mu) 就可以在 (mathcal O(max(n,m)^{frac37})) 时间内计算，具体推导见题解（太毒瘤了）

求解有很多技巧。

对于外层有个结论：

[令leftlfloor sqrt{leftlfloordfrac{n}{i} ight floor} ight floor=a，则这一段右端点是 leftlfloordfrac{n}{a^2} ight floor ]

对于筛 (mu^2)，也有个技巧。前面说我们对于 (ileq m^{frac{1}{3}}) 求，再对于 (leftlfloor sqrt{dfrac{m}{i}} ight floor leq m^{frac{1}{3}}) 求。具体如何实现？先记录最大的 (i,i^3 leq m^{frac{1}{3}})，求出 (sum_{jleq i} mu(j)leftlfloor {dfrac{n}{j^2}} ight floor)，并记录 (t=sum_{jleq i} mu(j)) 。然后依次求 (leftlfloor {dfrac{n}{j^2}} ight floor = leftlfloor {dfrac{n}{(i+1)^2}} ight floor...1) 的答案，根据上面的结论右端点可以求出来记为 (r)，我们只需给答案加上 (-t+sum_{jleq r} mu(j))。你会发现每个 (mu(i)) 恰好被前缀和覆盖了 (leftlfloor {dfrac{n}{i^2}} ight floor) 次，也就神仙地完成了计算。

#include <unordered_map>
#include <algorithm>
#include <cstdio>
#include <cmath>
using namespace std;
typedef long long ll;
const int N = 2e7 + 30;
int p[N >> 3], mu[N], sn, tot;
ll mu1[N], mu2[N];
bool tag[N];
void sieve() {
   mu[1] = 1;
   for(int i = 2; i <= sn; i++) {
      if(!tag[i]) p[++ tot] = i, mu[i] = -1;
      for(int j = 1; j <= tot && i * p[j] <= sn; j++) {
         tag[i * p[j]] = 1;
         if (i % p[j] == 0) break ;
         mu[i * p[j]] = -mu[i];
      }
   }
   for(int i = 1; i <= sn; i ++) {
      mu1[i] = mu1[i - 1] + mu[i];
      mu2[i] = mu2[i - 1] + (mu[i] != 0);
   }
}

ll n, m;
unordered_map<int, int> Map;
ll Mu(int n) {
   if(n <= sn) return mu1[n];
   if(Map.count(n)) return Map[n];
   int res = 1;
   for(int i = 2, j; i <= n; i = j + 1) {
      j = n / (n / i); res -= Mu(n / i) * (j - i + 1);
   }
   Map[n] = res;
   return res;
}
ll Mu2(ll n) {
   if (n <= sn) return mu2[n];
   ll res = 0, i;
   for(i = 1; i * i * i <= n; i ++) {
      res += mu[i] * (n / (i * i));
   }
   ll tmp = Mu(i - 1);
   for(ll j = n / (i * i); j >= 1; j --) {
      res += Mu(sqrt(n / j)) - tmp;
   }
   return res;
}
int main() {
   scanf("%lld%lld", &n, &m);
   if(n > m) swap(n, m);
   sn = min((ll) (pow(m, 3 / 7.0)), n); sieve();
   ll ans = 0, la = 0, cur = 0;
   for(ll l = 1, r; l <= n; l = r + 1) {
      ll a = sqrt(n / l), b = sqrt(m / l);
      r = min(n / (a * a), m / (b * b));
      ans += ((cur = Mu2(r)) - la) * a * b;
      la = cur;
   }
   printf("%lld
", ans);
   return 0;
}