企业官网建设流程全解析

1. 先把题意翻译成机器能算的语言

1.1 这道题在问什么，数据长什么样

题单里这道"对数难题"的原始描述一般很朴素：给一个长度为 n 的正整数序列，问有多少对下标 (i, j)，满足 i < j 且 log 以 a_i 为底、a_j 的对数值是一个整数。样本和数据范围通常给到 n ≤ 2×10^5，a_i ≤ 10^9。这个规模很有意思，它把暴力解按在地上摩擦，又没大到必须动用 Pollard Rho 那种重型武器，正好卡在"逼你想清楚数论性质"的区间上。

我第一次读题时的直觉是：底数 a_i 和真数 a_j 都是给定的，那我只要算出 log 值，判断它与最近的整数差得够不够小就行了。这个想法三分钟就能写完，然后在第三个测试点直接崩盘。原因很朴素——浮点数没有能力精确表达"这个对数值恰好是整数"这件事，尤其是在 10^9 这个量级上。后来我把整个思路推翻重来，才意识到出题人真正想考的不是对数运算，而是质因数分解之后的指数向量结构。所谓"对数"，只是把这个结构包装了一层外衣。

换个角度看，题目问的是"a_j 能不能写成 a_i 的某次幂"。当 a_i ≥ 2 且 a_j ≥ 1 时，log 值为整数和不成立，取决于 a_j 是否为 a_i 的整数次幂（包括 0 次幂，也就是 a_j = 1 的情况）。这个问题一旦翻译成"幂"的判定，浮点数就彻底没用了，剩下的全是整数运算。

1.2 浮点解为什么一定会翻车

很多人（包括我）第一反应是double v = log((double)a[j]) / log((double)a[i]);，然后判断fabs(v - round(v)) < 1e-9。这个写法在 a_i = 2、a_j = 8 这种小数据上永远是对的，于是产生了一种"我写对了"的错觉。

真正的问题在于误差是双向的、不可预测的。a_i = 3、a_j = 3^18 = 387420489 时，正确的 log 值应该是 18，但log(387420489)/log(3)在很多 libm 实现上会算出一个像 17.999999999999996 这样的结果。反过来，a_i = 7、a_j 取一个接近 7^10 但不等于它的数，算出来的比值可能落在 9.999999999999998 或 10.000000000000002，方向完全随机。你根本没法用一个固定的 eps 同时兜住这两种情况，因为误差量级和数值大小并不成简单的线性关系，它还跟分解方式、编译器的优化级别、甚至 CPU 的指令集有关。

这里有个很坑的细节：同一份代码在本地-O2下跑出错误答案，换到评测机上可能又是对的。别以为这是玄学，只是浮点运算路径变了。真正的解法是彻底不碰浮点。

我后来用 long double 试过一轮，把小数据全部对拍通过，但在 a_i = 2、a_j 接近 2^29 的边界上还是挂了。这不是精度不够的问题，是思路方向错了。判断两个大整数之间的幂关系，应该用整数的方式做。

1.3 把对数问题拆成指数向量

正整数的质因数分解是唯一的，这一点是整道题的基石。把 a_i 写成 ∏ p_k^{e_k}，把 a_j 写成 ∏ p_k^{f_k}，那么 a_j = a_i^m 等价于对每一个质因子 p_k，都有 f_k = m · e_k。也就是说，指数向量之间必须成统一的整数倍关系，倍数 m 必须对所有质因子都一样，不能一个因子是 2 倍、另一个是 3 倍。

这一步翻译带来的直接好处是：原本需要算对数的东西，变成了逐质因子的整数比较；原本需要担心精度的地方，全部消失了。剩下的难点只有两个——怎么快速分解 10^9 量级的数，以及怎么在 n = 2×10^5 的规模下统计出所有满足关系的数对。这两个问题正好对应标题里的"质因数分解"和"数论优化"，整篇文章其实就是在解决这两件事。

顺便说一个容易忽略的边界：a_j = 1 时，对任意 a_i ≥ 2 都有 a_i^0 = 1，所以 log 值等于 0，是整数，这一对要算进去。而 a_i = 1 时，log 以 1 为底本身没有定义，任何一对含 a_i = 1 的都不该计入。这个不对称性非常容易漏，我第一次提交就漏在这里，第二次又反向修过头，把 a_j = 1 的情况也一起删了。

2. 方案选型：三条路为什么最后只剩一条

2.1 暴力快速幂校验：想法对，复杂度不对

在决定用分解之前，我认真评估过一个更直接的做法：对每一对 (i, j) 暴力判断。判断逻辑是维护一个cur = a_i，然后不断乘 a_i，看看能不能正好等于 a_j。因为 a_i ≥ 2，最多乘 30 次就会超过 10^9，所以单次判断是 O(log a_j) 的常数操作。这个做法的正确性没有任何问题，整数乘法加上溢出剪枝，非常稳。

问题在于对数数量。n = 2×10^5，总对数大约是 2×10^10，哪怕每次判断只花 5 纳秒，也需要 100 秒。这条路在数据规模面前完全走不通，只能拿来做对拍程序。我在整场调试里，用这个暴力版本作为真理来源，随机生成几百组小数据去验证主程序，效果非常好，后面会详细讲对拍怎么搭。

还有一个中间方案：先把所有 a_i 去重，只对不同值之间做判定。如果数据重复率高，这个优化能把实际计算量砍掉一大截。但它会被构造数据针对——出题人只要让所有 a_i 互不相同，优化就失效了。所以它可以作为锦上添花的加速手段（我在代码里也确实加了去重缓存），但不能作为主方案的基础。

2.2 本原根：把"长得不同的底数"归纳成同一个 key

真正让这道题变简单的，是一个数学上的规范化操作。看两个数：a = 12 = 2^2 · 3^1，b = 144 = 2^4 · 3^2。它们的指数向量分别是 (2, 1) 和 (4, 2)，能不能配对？显然能，144 = 12^2。再来看 a = 36 = 2^2 · 3^2，b = 6 = 2^1 · 3^1，36 = 6^2，所以 log_6 36 = 2。这两个数的指数向量是 (2,2) 和 (1,1)。

把这两种情况放一起观察，会发现关键的规律：把指数向量除以它们的最大公约数，得到的"最简指数向量"是同一个。12 的指数 (2,1) 的 gcd 是 1，除以 1 还是 (2,1)；144 的指数 (4,2) 的 gcd 是 2，除以 2 得到 (2,1)，完全一致。同理，36 的 (2,2) 除以 2 得 (1,1)，6 的 (1,1) 除以 1 还是 (1,1)。

这个最简指数向量对应的那个整数，就是所谓的本原根。12 和 144 的本原根都是 12，6 和 36 的本原根都是 6。任何一个正整数 a 都能唯一地写成 c^h，其中 c 是它的本原根，h 是它的所有质因子指数的最大公约数。有了这个分解，判定就变得极其干净：两个数的本原根相同，才可能构成幂关系；如果本原根不同，直接判否。

再看指数那一层。设 a = c^{h_a}，b = c^{h_b}，那么 b = a^m 等价于 c^{h_b} = c^{m·h_a}，也就是 h_b = m · h_a，即 h_a 整除 h_b。整道题的判定条件被压缩成了一句话：本原根相同，且前者的 h 整除后者的 h。这个结论我推导完之后对着纸验算了十几个例子，包括 1、质数、平方数、完全立方数，全部成立。

2.3 一个让实现变简单的事实：本原根不会超过 10^9

推导到这里，我一直打算用哈希来表示本原根向量，因为直觉上 c 可能很大。但仔细想了一下：c^h = a ≤ 10^9 且 h ≥ 1，所以 c ≤ a ≤ 10^9，它本身就是一个能装进 int 的数。这个发现让整个实现成本骤降——我可以直接把 c 算出来当 key 用，根本不需要什么双哈希、随机权重、冲突检测。

把 c 算出来也不难，分解完之后把每个质因子的指数除以 g，再乘回去就行。中途乘积最大也就是 a 本身，不会溢出 32 位有符号整数，用 long long 更是绰绰有余。哈希方案我保留在代码注释里，作为值域放大到 10^18 时的备选，那时候 c 可能超出 64 位，才需要退回哈希。

选型这件事我踩过一个思维惯性：看到"向量"两个字就想哈希，实际上先算一下这个向量能不能塌缩成一个标量，很多时候能。省下来的不只是代码量，还有调试时间。

3. 分解环节：把 10^9 的数拆干净

3.1 筛法上界到底取多少

试除法的上界必须覆盖到 √(最大值)，这是基本要求。10^9 开方是 31622.776…，向上取整得到 31623。这个范围内的质数一共 3401 个，这个数字在很多题解里都出现过，我实际跑筛法验证过，确实是 3401。

筛法本身没什么争议，线性筛写起来最舒服，每个合数只被它的最小质因子筛掉一次，复杂度是严格 O(LIM)。代码只有十来行，我通常直接背下来默写。有一点要注意：内层循环里1LL * i * p > LIM这个判断必须在访问comp[i * p]之前做，而且要用 long long 乘，否则 i 和 p 都到三万级别的时候乘积会溢出 int，导致越界写。这个坑我在别的题上栽过，数组被写坏之后报的是段错误，查了半天才发现是类型问题。

static const int LIM = 31623; static vector<int> primes; static bool comp[LIM + 1]; void build_primes() { for (int i = 2; i <= LIM; ++i) { if (!comp[i]) primes.push_back(i); for (int p : primes) { if (1LL * i * p > LIM) break; comp[i * p] = true; if (i % p == 0) break; } } }

值域如果放宽到 10^12，上界就要取到 10^6，质数个数变成 78498 个。这时候对每个数都完整试除一遍就危险了，2×10^5 乘以 78498 是 1.5×10^10 次取模，绝对超时。这种情况需要换策略，要么上 Pollard Rho，要么先只筛到 10^6 以内的质数做部分分解，剩下的残余部分最多是"一个大质数"或"两个大质数的积"，再用完全平方检测和 Miller-Rabin 素性测试收尾。这道题的数据范围还没到那一步，但把这条退路想清楚，遇到加强版就不会慌。

3.2 试除循环的剪枝写法

分解的标准写法是枚举质数表，对每个质数 p 反复除尽。但直接写完整循环会做很多无用功，真正有效的剪枝是把循环条件改成p * p > rem，其中 rem 是当前还没分解完的剩余部分。理由很直接：如果 rem 的最小质因子都比 √rem 大，那 rem 只能是 1 或者一个质数，不可能再有更小的因子了。

这个剪枝的收益非常明显。考虑 a_i 是一个接近 10^9 的大质数，完整循环要走完 3401 个质数；但如果 a_i 是 2 乘以某个大质数，第一次除完 2 之后 rem 变成 5×10^8，下一轮 p = 3 时 9 < 5×10^8 继续，一直到 p 接近 22360 才退出，实际只走了 2500 个左右的质数。平均下来，整个序列的试除次数远小于 n × 3401 这个上界。

还有一个更狠的优化是针对重复值做缓存。用unordered_map<long long, Info>记住每个数分解后的结果，遇到相同的数直接查表。这在随机数据上收益有限，但在"所有数都相同"这类构造数据上能把分解次数从 2×10^5 降到 1。我一般都会加上，代价只是一次哈希查找。

注意：缓存一定要在分解之前查，不能分解完再存，否则该省的时间一分没省。这个低级错误我犯过不止一次。

3.3 残余部分、1 和特殊值的处理

循环退出时，rem 有三种可能：等于 1，说明所有质因子都被除干净了；大于 1 且是质数，说明它是最后剩下的那个大质因子；大于 1 且是合数——这种情况不可能发生，因为如果 rem 是合数，它必然有一个不超过 √rem 的质因子，而这个质因子早就该在循环里被枚举到了。

最后那个大质因子的指数必然是 1，这一点也需要单独验证。假设 rem = p^2 且 p > 31623，那么 rem > 10^9，超出值域，矛盾。所以残余部分一定是单次幂。

数字 1 需要单独拎出来。1 的质因数分解是空向量，gcd 是 0，本原根算出来也是 1，h 是 0，整个规范化流程会崩。我的处理方式是在主循环最开始就判断 x == 1，直接跳过分解，走单独的分支统计。这个分支的逻辑后面讲统计的时候展开。

还有 a_i = 2 这种最小的情况，分解出来 fac 里只有 (2, 1)，g = 1，本原根 c = 2，h = 1，一切正常。质数都是这个模式：c 等于它自己，h 等于 1。

3.4 指数向量的 gcd 和本原根的还原

分解出来的结果存成vector<pair<long long,int>>，每个元素是"质因子, 指数"。g 的初值设成 0，然后对每个指数做g = __gcd(g, e)。这里有个小技巧：__gcd(0, e)的结果就是 e，所以初值设 0 完全正确，不需要特判第一个元素。

g 算出来之后，再遍历一遍 fac，把每个指数除以 g 得到最简指数，然后累乘还原成本原根。累乘的时候可以加一个断言root <= 1000000000LL做自检，理论上它一定成立，但如果分解逻辑写错了，这个断言能第一时间抓住问题，比让错误数据流到统计阶段好得多。

struct Info { long long root; int h; }; Info normalize(long long x) { long long rem = x, root = 1; int g = 0; vector<pair<long long,int>> fac; for (int i = 0; i < (int)primes.size(); ++i) { long long p = primes[i]; if (p * p > rem) break; if (rem % p == 0) { int e = 0; do { rem /= p; ++e; } while (rem % p == 0); fac.emplace_back(p, e); g = __gcd(g, e); } } if (rem > 1) { fac.emplace_back(rem, 1); g = __gcd(g, 1); } for (auto &pr : fac) { int e = pr.second / g; while (e--) root *= pr.first; } return {root, g}; }

这段代码是整个方案的骨架，正确性直接决定后面所有统计。我当时用了大概二十组手工构造的数据验证它：完全平方数、完全立方数、质数的幂、多个质因子混合、指数互质导致 g = 1 的情况，全部对得上。

4. 统计环节：把 O(n²) 压到 O(n log)

4.1 一个被忽略的常数：h 最大只有 29

规范化之后，每个数变成了一个二元组 (root, h)，判定条件是 root 相同且 h_i 整除 h_j。第一反应是分组之后组内两两比较，最坏情况所有数本原根都一样，退化成 O(n²)。那就得想办法优化。

关键在于 h 的取值范围。h 是本原根的指数，满足 c^h ≤ 10^9 且 c ≥ 2，所以 h ≤ log2(10^9) ≈ 29.897，取整就是 29。2^29 = 536870912 在范围内，2^30 = 1073741824 超了。所以h 只可能是 1 到 29 之间的整数，一共 29 种取值。

这个观察把问题彻底改变了。组内不再需要两两比较，只需要统计每个 h 值出现了多少次。对于当前元素 (root, h)，它作为"真数"能配对的前面元素，是那些 h' 能整除 h 的。因为 h ≤ 29，直接枚举 h 的所有约数就够了，最多 8 个（24 有 8 个约数），单次查询是常数时间。

4.2 计数数组加枚举约数

具体做法是给每个 root 维护一个长度为 32 的计数数组cnt[32]，cnt[h]表示前面已经出现过多少个本原根为 root、指数为 h 的数。处理当前元素时，枚举 h 的所有约数 d，把cnt[d]累加到答案里，然后把cnt[h]加一。

枚举约数的写法我用的是标准的一半枚举：

for (int d = 1; d * d <= h; ++d) { if (h % d == 0) { ans += arr[d]; if (d * d != h) ans += arr[h / d]; } }

这里不需要按顺序遍历，也不需要排序，因为计数数组本身已经把所有信息聚合了。复杂度是 O(n × 8)，总共不到 200 万次操作，完全不构成瓶颈。

整个算法的复杂度因此变成：分解阶段 O(n × π(√V)) 的最坏情况，实际远小于这个上界；统计阶段 O(n × d(h))。两部分加起来，在 n = 2×10^5、V = 10^9 的数据上，我实测跑完大概 0.3 到 0.5 秒，非常安全。

4.3 方向性、自反对和 1 的特殊处理

这里有几个非常容易搞错的方向性问题，我逐个说清楚。

题目要求 i < j 且 log 以 a_i 为底，方向是固定的：前面的数当底，后面的数当真数。所以遍历时是"先查后插"，查询的是已经插入的前缀，而不是全量统计。顺序一旦搞反，答案就会多算一倍。

自反对的问题。按上面的写法，当前元素查询时还没有插入自己，所以不会自己配自己，正确。但如果先插入再查询，就会出现 h 整除自己的情况，多算 n 个。这个顺序必须严格执行。

1 的处理。a_j = 1 时，对前面任意 a_i ≥ 2 都有 a_i^0 = 1，log 值是 0，是整数，这一对要算。a_i = 1 时，log 以 1 为底没有定义，这一对不能算。所以代码里对 x = 1 单独走一个分支：把当前下标减去前面出现过的 1 的个数，就是前面非 1 元素的个数，全部累加进答案，然后 1 的计数加一。如果 x > 1，则它和前面的 1 不能配对，正常走分解流程即可。

还有一种变体题型要提一下：有的版本问的是"存在整数 k 使得 a_i^k = a_j 或 a_j^k = a_i"，这就没有方向性了，两个方向都要统计。处理办法是对每个 root 维护完整的计数，遍历完之后再做一次全局的整除关系统计，或者干脆跑两遍把方向拆开。看清题面里的"以谁为底"这几个字，能省掉很多返工。

4.4 完整实现

把分解和统计拼起来，主体代码大概是这样：

int main() { build_primes(); int n; scanf("%d", &n); unordered_map<long long, array<int,32>> cnt; cnt.reserve(n * 2); long long ans = 0, ones = 0; for (int i = 0; i < n; ++i) { long long x; scanf("%lld", &x); if (x == 1) { ans += i - ones; // 前面所有非 1 的数都能以它为真数配对 ++ones; continue; } Info in = normalize(x); auto &arr = cnt[in.root]; for (int d = 1; d * d <= in.h; ++d) { if (in.h % d == 0) { ans += arr[d]; if (d * d != in.h) ans += arr[in.h / d]; } } ++arr[in.h]; } printf("%lld\n", ans); return 0; }

unordered_map的operator[]对std::array<int,32>会做值初始化，所有元素清零，所以不需要手动memset。这一点我特意确认过，因为如果换成裸数组指针就会留下垃圾值。答案的类型必须是 long long，最坏情况 n = 2×10^5 且所有数都相等，答案接近 2×10^10，超过 int 范围。

5. 同一套骨架的三种变体

5.1 乘积为完全平方数对

题目换个问法：统计有多少对 (i, j) 使得 a_i × a_j 是一个完全平方数。这道题的解法思路和上面完全一致，只是规范化的方式从"除 gcd"变成了"指数模 2"。

a × b 是完全平方，等价于每个质因子的指数之和是偶数，也就是 e_a + e_b ≡ 0 (mod 2)。而在模 2 的意义下，-1 就等于