【UR #27】红场阅兵

题意

给定一个积性函数 \(S\)，保证 \(S(p^k)\) 是关于 \(p^k\) 的二次函数，求 \(\sum_{i=1}^n \sum_{j=1}^n S(ij)\)。

数据范围：\(n\le 10^9\)。

算法一

为了解决上面的问题，我们需要把积性函数推广到二维。

首先二维数论函数 \(g(x,y)\) 还是没啥性质的普通函数。

我们称二维数论函数 \(f(x,y)\) 是二维积性函数，当且仅当 \(\forall a,b,c,d\) 使得 \(\gcd(ab,cd)=1\)，都有 \(f(ac,bd) = f(a,b)f(c,d)\)。

上面这个定义还是有点太抽象了（这是按照一维积性函数的定义改的）。实际上可以更简单地理解这件事情。假设 \(p_i\) 是从小到大第 \(i\) 个质数，\(x=\prod_i p_i^{a_i}\)，\(y=\prod_i p_i^{b_i}\)，那么 \(f(x,y) = \prod_{i} f(p_i^{a_i},p_i^{b_i})\)。

容易发现，我们要求的是一个二维积性函数 \(f(x,y) = S(xy)\) 的二维前缀和。

类似地，我们定义二维 Dirichlet 卷积。设 \(a,b\) 的 Dirichlet 卷积为 \(c\)，那么：

\[c(x,y) = \sum_{u|x,v|y} a(u,v) b(x/u,y/v) \]

类似一维地，\(c(1,1) \neq 0\) 的数论函数存在逆元；两个二维积性函数的 Dirichlet 卷积和 Dirichlet 除法都是二维积性函数。

类似 Powerful Number 地，我们考虑用 Powerful Number 的思想来解决问题。

考虑设计积性函数 \(g(x,y)\) 满足 \(g(p^{a},p^b) = f(p^a,1) f(1,p^b)\)。对于这个函数，其前缀和 \(\sum_{i\le n_1} \sum_{j\le n_2} f(ij)\) 是相当容易计算的。

观察 \(h=f/g\)。对于单个 \(p\) 来说，我们发现这个 Dirichlet 除法相当于是二维多项式的除法。于是我们发现 \(h(p^k,1) = h(1,p^k) = 0\)。

\[\begin{aligned} \sum_{i=1}^n \sum_{j=1}^n f(i,j) &= \sum_{i=1}^n \sum_{j=1}^n \sum_{u|i} \sum_{v|j} h(u,v) g(i/u,j/v)\\&= \sum_{u,v} h(u,v) (\sum_{i\le \frac{n}{u},j\le \frac{n}{v}} g(i,j)) \end{aligned} \]

筛出 \(f(p^a,1),f(1,p^b)\) 的前缀和以及一些预处理之后，考虑暴力枚举 \(h\) 有值的位置。复杂度就是 \(x,y\in [1,n]\) 且 $ h(x,y)\neq 0$ 的 \((x,y)\) 个数。

实际上，该算法的时间复杂度是 \(\mathcal O(n)\)。。期望得分 \(70\) 分。

接下来给出复杂度证明。

考虑将条件 \(x,y\in [1,n]\) 弱化成 \(xy\le n^2\)。接下来设 \(s(xy)\) 表示 \(xy=s\) 的合法 \((x,y)\) 对数。\((x,y)\) 是合法的当且仅当，且对于 \(x=\prod_i p_i^{a_i},y=\prod_i p_i^{b_i}\)，都有 \([a_i=0] = [b_i=0]\)。可以发现 \(s\) 是积性函数。

对于一个质数 \(p\) 和正整数 \(k\)，\(s\) 满足 \(s(p^k)=k-1\)。

实际上，\(s\) 前 \(n\) 项的和是 \(\mathcal O(\sqrt n)\) 级别的。

考虑设计一些类似 Powerful Number 的复杂度分析。

考虑找出一个积性函数 \(g\) 使得 \(\forall i, g(i) \ge s(i)\) 且 \(g\) 的前缀和好算。

定义 \(T_k(x)=1\) 当且仅当存在 \(y\) 使得 \(x=y^k\)。

那么，一个合法的选择是 \(g = T_2*T_3^3*T_4^3*T_5^3\) 是一个合法的选择（随便打的一串东西）。

考虑分析 \(g\) 的前缀和的上界。

将生成 \(g\) 的过程看成：开始 \(g=T_2\)。接下来每次乘 \(T_k\)（\(k\ge 3\)）。

刚开始的时候，\(g\) 的前缀和是 \(\mathcal O(\sqrt n)\) 级别的。

乘一次（比如 \(T_3\)）后，在 \(n\) 处的前缀和是：

\[\sum_{i\le n^{1/3}} \sqrt {\frac{n}{i^3}} \]

类似 PN 的分析，这个值还是是 \(\mathcal O(\sqrt n)\) 级别的。当然，对于 \(k\ge 4\)，这个值不会更大，也只能是 \(\mathcal O(\sqrt n)\) 级别的了。

所以，多乘一个 \(T_k\)（\(k\ge 3\)）后并不会改变复杂度量级，那多乘常数个也是一样。因此，\(g\) 函数在 \(n\) 处的前缀和是 \(\mathcal O(\sqrt n)\) 的。

因此，原问题的复杂度是 \(\mathcal O(n)\) 的。

在对满足存在序列 \(a_k\) 使得 \(f(p^k)=a_k\)（\(k\ge 1\)） 的一维积性函数求和时，我们都可以利用这一类方法。在这里需要满足 \(a\) 序列非负而且 \(a\) 序列是常数次多项式级别的。

此时我们我们其实只用关心 \(a\) 第一个非零的位置 \(p\)（\(p\ge 1\)）。\(a\) 的前缀和就不会超过 \(\mathcal O(n^{1/p}\log^{a_p-1} n)\) 级别。

举例：\(d(x)^2\) 的前缀和就是 \(\mathcal O(n\log^3n)\)级别的。

尽管这类分析只给出了一个上界，但是我觉得下界在很大程度上也是对的。

算法二

考虑优化之前的算法，设计二维积性函数 \(g_3\) 满足 \(g_3\) 只在 \(g_3(x,x)\) 处有值且值为 \(h(x,x)\)，接下来令 \(h' = h/g_3\)。

此时有值的 \(h'(p^a,p^b)\)（\(p\in \mathbb{P}, a,b\ge 0\)，\(a+b\neq 0\)）就满足 \(a,b\ge 1\) 且 \(a\neq b\) 了。

接下来暴力枚举 \(h'\) 中的有值元素 \((x,y)\)，计算 \(\sum_{u \le n/x,v\le n/y} (g*g_3) (u,v)\)，把这些答案带上系数加起来就能得到最终答案了。

设 \(A=n/x,B=n/y\)，那么：

\[\begin{aligned} \sum_{i=1}^{A}\sum_{j=1}^{B} (g*g_3)&= \sum_{k} g_3(k,k) \sum_{i=1}^{A/k}\sum_{j=1}^{B/k} g(i,j) \end{aligned}\]

用整除分块优化这个过程，直接算的复杂度便是 \(\mathcal O(\min(\sqrt A + \sqrt B, A, B))\)。

实际上，除去筛 \(g_3\) 和 \(f(x,1)\) 的前缀和的复杂度，上述做法的时间复杂度是 \(\mathcal O(n^{2/3}\log n)\)。接下来我们给出证明。

首先观察 \(h'\) 有值位置个数。如果用之前的做法分析，我们只能分析出 \(\mathcal O(n^{2/3}\log n)\) 的上界；但实际上可以分析到 \(\mathcal O(n^{2/3})\) 的。

我们要分析的是满足 \(v(p^a,p^b) = [a\neq b,a,b\ge 1]\)（\(a+b>0\)）的二元积性函数前缀和。不妨考虑拆成二维积性函数 \(v_1\) 和 \(v_2\)，其中 \(v_1(p^a,p^b)\)只在 \(a=1,b=2\) 或 \(a=2,b=1\) 处有值且值为 \(1\)，\(v_2(p^a,p^b)=v(p^a,p^b)-v_1(p^a,p^b)\)（这里都是 \(a+b>0\)）。

\(v\) 不会超过 \(v_1\times v_2\)。

首先分析 \(\sum_{i=1}^n \sum_{j=1}^m v_1(i,j)\)。

这个不会超过 \(\sum_{a,b} [a^2b \le n] [b^2a\le m]\)，即 \(\sum_{a^2\le n, a\le m} \min(\frac{n}{a^2}, \sqrt {\frac{m}{a}})\)。

这个 \(\min\) 在 \(\frac{n}{a^2} = \sqrt {\frac{m}{a}}\) 时取等，此时 \(a = \sqrt[3]{\frac{n^2}{m}}\)，设其为 \(B\)。

所以复杂度不会超过：

\[\sum_{a=1}^{B} \sqrt {\frac{m}{a} } + \sum_{a>B} \frac{n}{a^2} \]

左边是 \(\mathcal O(\sqrt {mB})\)，右边是 \(\mathcal O(\frac{n}{B})\)，带入 \(B\) 发现这两者都是 \(\mathcal O((nm)^{1/3})\)，所以复杂度是 \(\mathcal O((nm)^{1/3})\)。

根据前面的分析，满足 \(xy\le L\) 的 \(v_2(x,y)\) 是不超过 \(\mathcal O(L^{1/4}\log L)\) 的。

复杂度不会超过 \(v_1\times v_2\) 的前缀和，也不会超过：

\[\sum_{x,y,v_2(x,y)>0,xy \le n^2} \sqrt [3] {(n/x) * (n/y)} \]

考虑将所有 \(x,y\) 分组。更为具体地，把 \(\lfloor\log_2 xy \rfloor\) 相同的 \(x,y\) 分到一组。我们可以将复杂度近似地分析为：

\[\sum_{k=1} k \sqrt[4] {2^{k}} \sqrt[3] {n^2/2^k} \]

计算得复杂度是 \(\mathcal O(n^{2/3})\)。这个分析对于 \(\sum_{i=1}^n\sum_{j=1}^m v(i,j)\) 的分析也有效，能分析出 \(\mathcal O((nm)^{1/3})\)。

我们需要分析的是 \(\sum_{x,y} v_2(x,y) \min(n/x,n/y,\sqrt {n/x} + \sqrt{n/y})\)。

首先，\(\min(A,B,\sqrt A + \sqrt B) \le \min(A,\sqrt B) + \min(B,\sqrt A)\)。两者是对称的，我们只拿 \(\min(A,\sqrt B)\) 分析。

那么复杂度可以分析为

\[\begin{aligned} &\sum_{x,y} v_2(x,y) \min(n/x,\lfloor\sqrt{n/y} \rfloor) \\&= \sum_{x,y} v_2(x,y) \sum_{tx\le n,t^2y\le n} 1 \\&= \sum_{t^2\le n}\sum_{x\le n/t,y\le n/t^2} v_2(x,y) \end{aligned} \]

这个是 \(\sum_{t \le \sqrt n} (n^2/t^3)^{1/3}\) 级别的，即 \(\mathcal O(n^{2/3}\log n)\)。期望得分 \(95\sim 100\) 分。

算法三

在 \(n,m\) 同阶的时候，时间复杂度实际上可以做到 \(\mathcal O(n^{2/3})\)。

对于计算 \(\sum_{i} g_3(i,i) preS(\lfloor \frac{x}{i} \rfloor) preS(\lfloor \frac{y}{i} \rfloor)\)（\(x\ge y\)），如果我们能做到 \(\mathcal O(n^{2/3})\) 预处理和 \(\mathcal O(y^{0.6})\) 查询一组 \((x,y)\)，那么总复杂度就可以做到 \(\mathcal O(n^{2/3})\)。

原因是，设 \(w(y) = y^{0.6}-(y-1)^{0.6}\)，是 \(\mathcal O(y^{-0.4})\) 级别的。我们的算法的复杂度是：

\[\sum_{x,y} v_2(x,y) \sum_{j=1}^{\lfloor \frac{n}{\max(x,y)} \rfloor} w(j)= \sum_{t} w(t) \sum_{x,y\le n/t} v_2(x,y) \]

也就是 \(\mathcal O(\sum_{t} t^{-0.4} (n/t)^{2/3}) = \mathcal O(n^{2/3}\sum_{t} t^{-16/15}) = \mathcal O(n^{2/3})\)。

接下来的问题就是怎么做到 \(\mathcal O(n^{2/3})\) 预处理和 \(\mathcal O(y^{0.6})\) 查询。

对于 \(x\ge n^{2/3}\)，我们可以对所有可能的 \(x\) 预处理出来他们的整除分块，即处理出 \(\sum_{i=1}^{\lfloor \frac{x}{i} \rfloor} g(i,i) preS(\lfloor \frac{x}{i} \rfloor)\)。查询直接对 \(y\) 整除分块即可做到 \(\mathcal O(\sqrt y\log y)\) 或 \(\mathcal O(\sqrt y)\)。这部分的总复杂度和杜教筛一样，为 \(\mathcal O(n^{2/3})\)。

对于 \(x\le n^{2/3}\)，一个简单的想法是从小到大枚举 \(i\)，并直接维护 \(a_i = g(i,i) preS(\lfloor \frac{x}{i} \rfloor)\)。使用 \(\mathcal O(1)\) 修改 \(\mathcal O(p^{o(1)})\) 查询 \(p\) 位置的前缀和的数据结构，即可做到 \(\mathcal O(y^{0.5+o(1)})\) 做单次查询。但是这样做的问题在于我们需要 \(\mathcal O(n^{2/3}\log n)\) 预处理。即使调整目前是 \(n^{2/3}\) 的阈值，我们也只能做到 \(\mathcal O(n^{2/3}\log^{1/3} n)\)。

首先，我们其实只需要维护 \(i\le \sqrt x\) 的 \(g(i,i) preS(\lfloor \frac{x}{i} \rfloor)\) 和 \(S(i) (preg(\lfloor \frac{x}{i} \rfloor) - preg(\lfloor \frac{x}{i+1} \rfloor))\)。只要能维护这两者的前缀和就够了。

考虑按从小到大的顺序维护所有可能的 \(x\)（即满足存在整数 \(p\) 使得 \(x=\lfloor \frac{n}{p} \rfloor\) 的 \(x\)）对应的数组 \([\lfloor \frac{x}{1} \rfloor, \lfloor \frac{x}{2} \rfloor, ..., \lfloor \frac{x}{\sqrt x} \rfloor]\)。

结论：这个数组的总变化是 \(\mathcal O(n^{2/3})\) 的。这也意味着我们可以 \(\mathcal O(n^{2/3})\) 解决原问题。

我们需要在意的是，对于 \(p \in [n^{1/3},n^{1/2}]\)，\(\lfloor \frac{n}{p} \rfloor\) 和 \(\lfloor \frac{n}{p+1} \rfloor\) 对应数组的差别。

设 \(r=\frac{n}{p^2}\)。那么对于数组 \(\le r\) 的部分，变化量显然不会超过 \(r\)；对于数组 \(>r\) 的部分，变化量大约为

\[\begin{aligned} &\sum_{i=r}^{\sqrt {n/p}} \frac{n}{ip} - \frac{n}{i (p+1)} \\&= \mathcal O(\frac{n}{p^2} (\ln(\sqrt {\frac{n}{p} }) - \ln(r))) \\&= \mathcal O(\frac{n}{p^2} (\ln(p^3) - \ln n)) \end{aligned}\]

直接用调和级数分析并不严谨，但是可以用倍增等方式让证明变得严谨。

所以总变化量是 \(\mathcal O(\sum_{p \in [n^{1/3},n^{1/2}]} \frac{n}{p^2} (\ln(p^3) - \ln n))\) 的。

为了分析这个，考虑将 \(p \in [n^{1/3}2^k,n^{1/3}2^k)\) 的 \(p\) 放在一起分析。

所以复杂度为 \(\mathcal O(\sum_{k=0} n^{1/3}2^k\times \frac{n^{1/3}}{4^k}\times k) = \mathcal O(n^{2/3}\sum_{k=0} \frac{k}{2^k})\)。

\(\sum_{k=0} \frac{k}{2^k}\) 是 \(\mathcal O(1)\) 的，所以复杂度就是 \(\mathcal O(n^{2/3})\)。期望得分 \(100\) 分。