peltorator

Автор статьи: Александр Голованов

Задача: Дана мультипликативная функция $f(n)$, а также число $N$. Необходимо посчитать $s_f(N) = \sum_{k=1}^Nf(k)$.

Мы умеем делать это за $O(n^{2/3})$, если нам известна мультипликативная функция $g$, для которой мы умеем считать $s_g$ и $s_{f*g}$ за $O(1)$ при помощи обращения Мёбиуса. Оказывается, зачастую это можно делать быстрее, если функции обладают некоторыми дополнительными свойствами.

Пусть $g(n)$ и $h(n)$ — такие мультипликативные функции, что $f = h*g$. Тогда верна следующая лемма:

Лемма: $$s_f(n) = \sum_{k=1}^nh(k)s_g\left(\left\lfloor\frac{n}{k}\right\rfloor\right)$$

Доказательство: В сумме слева есть все возможные слагаемые вида $h(k) \cdot g(l)$, где $k \cdot l \le n$. Таким образом, для фиксированного $k$ мы должны перемножить $h(k)$ со всеми такими $g(l)$, что $l \le \left\lfloor \frac{n}{k} \right\rfloor$, то есть с префиксной суммой $g$. Что и требовалось доказать.

Пускай функция $h$ такова, что $h(p) = 0$ для любого простого $p$. Тогда заметим, что в правой части выражения из леммы те слагаемые, в которых $k$ содержит в себе некоторый простой множитель ровно в первой степени, зануляются (из свойства функции $h$). Остальные же числа содержат каждый свой простой множитель хотя бы во второй степени. Будем называть такие числа сильными (от английского powerful, которое происходит от слова power в значении «степень»).

Лемма: Всякое сильное число $n$ можно представить в виде $a^2b^3$ для некоторых натуральных $a$ и $b$.

Доказательство: Для этого нужно представить каждую степень простого $p_i^{\alpha_i}$, входящую в $n$, в таком виде, а потом перемножить их. То есть надо просто представить степень $\alpha_i$ в виде $2x + 3y$. Безусловно, это можно сделать, так как $\alpha_i \ge 2$ (для четного $n$ можно взять $x = \frac{n}{2}$ и $y = 0$, а для нечетного подойдет решение $x = \frac{n - 3}{2}$ и $y = 1$).

Тогда для вычисления $s_f$ нам нужно научиться эффективно перебирать все сильные числа $k$ и суммировать $h(k) s_g\left(\left\lfloor \frac{n}{k} \right\rfloor\right)$. Будем делать это рекурсивно, перебирая простые числа по возрастанию и выбирая их в степени либо ноль, либо хотя бы два. Обратите внимание на то, что если простое число $p$ входит в разложение сильного числа $k$, то $p \le \sqrt{k} \le \sqrt{n}$, потому что каждый простой множитель входит в $k$ хотя бы во второй степени, поэтому перебирать простые мы можем только до $\sqrt{n}$.

Рассмотрим следующий алгоритм вычисления функции $s_f$:

T s_f(long long n, long long cur, int idx, bool new_value = true) {
    T ans = new_value ? h(cur) * s_g(n / cur) : 0;
    if (n / cur / primes[idx] / primes[idx] == 0) { // can't multiply anymore
        return ans;
    }
    ans += s_f(n, cur, idx + 1, false);
    cur = cur * primes[idx] * primes[idx];
    while (cur <= n) {
        ans += s_f(n, cur, idx + 1, true);
        cur *= primes[idx];
    }
    return ans;
}

Здесь primes — массив с простыми числами хотя бы до $\sqrt{N}$; cur — текущее значение $k$, которое мы перебираем рекурсивно; idx — индекс текущего простого из списка primes; new_value — флаг, который нужен, чтобы мы не прибавили одно и то же значение к ответу несколько раз, когда мы не берем текущее простое в разложение; а тип T может быть любым типом на Ваш выбор: например, long long или Ваш класс для модульной арифметики.

Легко видеть, что этот алгоритм вычисляет $s_f$, перебирая лишь сильные числа. Оценим время его работы. Будем считать, что значение функции $s_g$ в точке $n$ мы можем вычислить за время $T_{s_g}(n)$, а вычисление $h(cur)$ работает за $O(1)$ (оно может считаться на лету по мультипликативности через значения в степенях простых). Тогда

$$T_{s_f(n)} = \sum_{k\text{ is powerful}}T_{s_g}\left(\left\lfloor\frac{N}{k}\right\rfloor\right)\leq\sum_{a, b : a^2b^3\leq N}T_{s_g}\left(\frac{N}{a^2b^3}\right)\approx\int\limits_{1}^{\sqrt{N}}\int\limits_{1}^{\sqrt[3]{N/a^2}}T_{s_g}\left(\frac{N}{a^2b^3}\right)\ \mathrm{d}b\ \mathrm{d}a$$

Для простоты предположим, что $T_{s_g}(n) = n^c$, где $c > 1/3$. Тогда

$$ \int\limits_{1}^{\sqrt{N}}\int\limits_{1}^{\sqrt[3]{N/a^2}}T_{s_g}\left(\frac{N}{a^2b^3}\right)\ \mathrm{d}b\ \mathrm{d}a = \int\limits_{1}^{\sqrt{N}}\int\limits_{1}^{\sqrt[3]{N/a^2}}\frac{N^c}{a^{2c}b^{3c}}\ \mathrm{d}b\ \mathrm{d}a = N^c \int\limits_{1}^{\sqrt{N}} a^{-2c} \int\limits_{1}^{\sqrt[3]{N / a^2}} b^{-3c} \ \mathrm{d}b\ \mathrm{d}a = $$

$$ = N^c \int\limits_{1}^{\sqrt{N}} a^{-2c} \cdot 1/(3c - 1) \cdot (1 - const_1) \ \mathrm{d}a \le $$

$$ \le 1/(3c - 1)N^c\int\limits_{1}^{\sqrt{N}} a^{-2c} \ \mathrm{d}a\leq \begin{cases} 1/((1 - 2c)(3c - 1))N^c\cdot N^{1/2 - c} = O(\sqrt{N}), & c < \frac{1}{2}, \\ 1/(3c - 1)N^c\ln\sqrt{N} = O(\sqrt{N}\log{n}), & c = \frac{1}{2}, \\ 1/((3c - 1)(2c - 1))N^c = O(T_{s_g}(N)), & c > \frac{1}{2}. \end{cases} $$

Где $const_1$ — какое-то положительное число.

Очевидно, если $c \le 1/3$, то время работы и подавно составляет $O(\sqrt{N})$.

Пример: Если $f(n)=\mathrm{rad}(n)$, то есть $f(p_1^{\alpha_1} \cdot \ldots \cdot p_k^{\alpha_k}) = p_1p_2 \cdot \ldots \cdot p_k$, то в качестве $h(n)$ можно взять мультипликативную функцию, для которой верно $h(1) = 1$, $h(p) = 0$ и $h(p^k) = p - p^2$ для всех $k > 1$, а в качестве $g$ можно взять $g(n) = n$. Легко проверить, что $f = g * h$. Тогда сумма $f(n)$ на префиксе размера $N$ будет считаться за $O(\sqrt{N})$, потому что $T_{s_g}(n) = O(1)$.

Замечание: Как видно из примера, при поиске функций $h$ и $g$ имеет смысл подобрать их значения на степенях простых.

Замечание: При оценке времени работы мы пользовались представимостью сильных чисел в виде $a^2b^3$. Казалось бы, можно было бы написать два for-а по числам $a$ и $b$, и такой алгоритм был бы проще рекурсивного перебора сильных чисел. Однако сильные числа бывают представимы в таком виде больше, чем одним способом, поэтому при таком алгоритме придётся каким-то образом проверять числа на уникальность, что может повлиять на время работы.