Занятие 1. Параллель A — Математика

1. Асимптотики и практические оценки 5:49

Лекция начинается с базового языка оценок сложности. Формально запись означает следующее.

Именно поэтому про бинарный поиск говорят «работает за $O(\log n)$»: нотация $O$ допускает и любую меньшую сложность, тогда как $o(\log n)$ уже означало бы строго меньше логарифма.

Это теоретический факт: на практике размеры входа ограничены, и константы у «страшной» на вид функции могут быть крошечными. Поэтому программисты держат в голове ещё и прикладные ориентиры: сортировка — $O(n\log n)$, поиск в unordered_set — между $O(1)$ и $O(\log)$, обход DFS — $O(n+m)$, непересекающиеся множества (DSU) — почти линейно (обратная функция Аккермана), алгоритм Куна — порядка $O(m)$ на реальных данных, потоки — около $O(m)$ при малых весах.

2. НОД, расширенный Евклид и обратный по модулю 11:07

Наибольший общий делитель устроен просто и быстро. Оценку времени работы удаётся уточнить до логарифма от отношения $\min$ к самому НОД.

Расширенный алгоритм Евклида

Он не только возвращает $g=\gcd(a,b)$, но и коэффициенты $x,y$ с $a x+b y=g$. Рекурсия повторяет обычный Евклид, только теперь на каждом шаге нужно переписать коэффициенты.

Пусть рекурсивный вызов для пары $(b,\ a\bmod b)$ нашёл $x_1,y_1$: $$ b\,x_1+(a\bmod b)\,y_1=g. $$ Подставим $a\bmod b=a-b\lfloor a/b\rfloor$: $$ b\,x_1+\bigl(a-b\lfloor a/b\rfloor\bigr)y_1 = a\,y_1 + b\bigl(x_1-\lfloor a/b\rfloor\,y_1\bigr)=g. $$ Значит, для исходной пары надо взять $$ x=y_1,\qquad y=x_1-\left\lfloor \frac{a}{b}\right\rfloor y_1. $$ База: если $b=0$, то $g=a$, $x=1$, $y=0$.

// ax + by = gcd(a, b)
long long extgcd(long long a, long long b, long long &x, long long &y) {
    if (b == 0) { x = 1; y = 0; return a; }
    long long x1, y1;
    long long g = extgcd(b, a % b, x1, y1);
    x = y1;
    y = x1 - (a / b) * y1;
    return g;
}

Обратный элемент через расширенный Евклид

Найти $x$ такое, что $a x\equiv 1\pmod n$, — это в точности решить диофантово уравнение: $$ a x\equiv 1\pmod n \iff a x - n y = 1. $$ Такое уравнение разрешимо тогда и только тогда, когда $\gcd(a,n)=1$. Решение даёт расширенный Евклид, а ответ берётся по модулю $n$: $$ a^{-1}\equiv x \pmod n. $$ Если $\gcd(a,n)>1$, обратного не существует: умножение на такое $a$ «склеивает» классы и может дать нуль.

3. Функция Эйлера и теорема Эйлера 17:03

Этого свойства достаточно, чтобы получить все формулы. Для простого $p$ взаимно просты все числа, кроме кратных $p$: $$ \varphi(p)=p-1. $$ Для степени простого выкидываем числа $p,2p,\dots,p^{\alpha-1}p$ — их ровно $p^{\alpha-1}$: $$ \varphi(p^{\alpha})=p^{\alpha}-p^{\alpha-1}=p^{\alpha-1}(p-1). $$ Разложив $n=\prod p_i^{\alpha_i}$ и воспользовавшись мультипликативностью, получаем удобную вычислительную формулу $$ \varphi(n)=n\prod_{p\mid n}\left(1-\frac1p\right). $$

Теорема Эйлера

Отсюда получается способ обращения: $$ x^{-1}\equiv x^{\varphi(m)-1}\pmod m, $$ поскольку $x\cdot x^{\varphi(m)-1}=x^{\varphi(m)}\equiv 1$. Этот способ работает, но требует знать $\varphi(m)$, то есть разложить $m$. Поэтому на практике для произвольного модуля надёжнее расширенный Евклид; формула Эйлера хороша как готовый приём при известном простом $p$: $x^{-1}\equiv x^{p-2}\pmod p$.

4. Быстрое возведение в степень, modint и все обратные за $O(n)$ 20:23

Стандартный бинарный алгоритм возведения в степень обрабатывает двоичную запись показателя: квадрат основания на каждом шаге, умножение результата — по младшему биту.

long long mulmod(long long a, long long b, long long mod) {
    return (long long)((__int128)a * b % mod);
}

long long power(long long a, long long e, long long mod) {
    long long res = 1 % mod;
    a %= mod;
    while (e > 0) {
        if (e & 1) res = mulmod(res, a, mod);   // младший бит показателя
        a = mulmod(a, a, mod);                  // переходим к следующему биту
        e >>= 1;
    }
    return res;
}

Структура modint

Чтобы не ловить баги при вычислениях по модулю, всю арифметику заворачивают в структуру. Она хранит одно число, приведённое к диапазону $[0,\mathrm{mod})$, и перегружает операторы.

const long long MOD = 1'000'000'007;

struct mint {
    long long v;
    mint() : v(0) {}
    mint(long long x) : v(((x % MOD) + MOD) % MOD) {}   // приводим сразу

    mint& operator+=(const mint& o){ v += o.v; if (v>=MOD) v-=MOD; return *this; }
    mint& operator-=(const mint& o){ v -= o.v; if (v<0)   v+=MOD; return *this; }
    mint& operator*=(const mint& o){ v = v * o.v % MOD; return *this; }
};
mint operator+(mint a, const mint& b){ return a += b; }
mint operator-(mint a, const mint& b){ return a -= b; }
mint operator*(mint a, const mint& b){ return a *= b; }

mint pw(mint a, long long e){
    mint r(1);
    while (e > 0) { if (e & 1) r *= a; a *= a; e >>= 1; }
    return r;
}

Все обратные к массиву за $O(n)$

Пусть даны $a_1,\dots,a_n$ без нулей, и нужно найти $a_i^{-1}$ по простому модулю. Достаточно одного обращения: считаем префиксные произведения $P_i=a_1\cdots a_i$, обращаем только $P_n$, а дальше идём назад и получаем $$ a_i^{-1}=P_n^{-1}\cdot P_{i-1}\cdot S_{i+1}, $$ где $S_{i+1}=a_{i+1}\cdots a_n$ — суффиксное произведение, уже посчитанное. Тот же приём в частном случае факториалов даёт «обратные факториалы» за $O(n+\log p)$.

// обратные ко всем a[0..n-1] по простому модулю
vector<mint> pref(n+1, 1), inv(n);
for (int i = 0; i < n; ++i) pref[i+1] = pref[i] * a[i];
mint allInv = pw(pref[n], MOD-2);              // единственное обращение
for (int i = n-1; i >= 0; --i) {
    inv[i] = allInv * pref[i];                 // домножаем на суффикс
    allInv *= a[i];                            // переходим к суффиксу без a[i]
}

5. Китайская теорема об остатках 29:56

Задача Сунь-цзы: решить систему сравнений $$ x\equiv a_i \pmod{m_i},\qquad i=1,\dots,n. $$ С одним сравнением всё ясно: решения — это $a_1+m_1 k$, $k\in\mathbb Z$.

С двумя сравнениями приравниваем выражения $x=a_1+m_1k_1$ и $x=a_2+m_2k_2$: $$ m_1k_1-m_2k_2=a_2-a_1. $$ Это линейное диофантово уравнение относительно $k_1,k_2$; оно разрешимо тогда и только тогда, когда $$ \gcd(m_1,m_2)\mid (a_2-a_1). $$ Если решение есть, пара сравнений заменяется одним — по модулю наименьшего общего кратного. Далее действуем индуктивно: сливаем первые два сравнения, к результату присоединяем третье, и так далее. Если на каком-то шаге решения нет, то и вся система несовместна.

Явная формула: пусть $M=\prod m_i$, $M_i=M/m_i$ и $N_i\equiv M_i^{-1}\pmod{m_i}$. Тогда $$ x\equiv \sum_{i=1}^{n} a_i\,M_i\,N_i \pmod M. $$ На практике пишут инкрементальный алгоритм, который в каждый момент хранит одно сравнение $x\equiv a\pmod m$.

// слить x ≡ a (mod m) и x ≡ b (mod n) в одно сравнение
// x = a + m*t,  m*t ≡ b - a (mod n)
bool crt_merge(long long &a, long long &m, long long b, long long n) {
    long long x, y, g = extgcd(m, n, x, y);      // m*x + n*y = g
    if ((b - a) % g != 0) return false;          // нет решения
    long long t = ((b - a) / g % (n / g) * x) % (n / g);
    long long l = m / g * n;                      // lcm(m, n)
    a = ((a + m * t) % l + l) % l;
    m = l;
    return true;
}

Применение: хэши по нескольким модулям

Если хранить значение по нескольким маленьким простым модулям, то КТО говорит, что это эквивалентно хранению числа по модулю их произведения. Отсюда — «псевдообоснование» того, почему строковые хэши по нескольким модулям устойчивы к взлому: чем больше попарно простых модулей, тем больше период. Оговорка: произведение не простое, поэтому часть свойств поля теряется.

6. Гармонический ряд, делители и оценки 39:03

Ключевой факт, на котором держится почти вся теория делителей: $$ \sum_{i=1}^{n}\frac1i=\Theta(\log n). $$ Отсюда сразу получается суммарное число делителей всех чисел до $n$. Единица делит $n$ чисел, двойка — примерно $n/2$, тройка — $n/3$, и так далее: $$ \sum_{i=1}^{n}\tau(i)=\sum_{d=1}^{n}\left\lfloor\frac{n}{d}\right\rfloor=\Theta(n\log n), $$ где $\tau(i)$ — число делителей числа $i$. Именно поэтому «перебрать все делители всех чисел до $n$» работает за $n\log n$, а не за $n\sqrt n$.

Отдельно полезна оценка числа делителей одного числа. Теоретически $\tau(n)$ растёт медленнее любой степени — для любого $\varepsilon>0$ при больших $n$ имеем $\tau(n)\lt n^{\varepsilon}$, — но на практике для реальных диапазонов удобны такие ориентиры.

7. Дискретное логарифмирование (baby-step giant-step) 46:48

Задача: по данным $a,b,n$ найти $x$ с $$ a^{x}\equiv b \pmod n. $$ Сначала рассматриваем случай $\gcd(a,n)=1$ (тогда и $\gcd(b,n)=1$, иначе решений нет, так как $a^x$ взаимно прост с $n$). Из теоремы Эйлера показатель определён по модулю $\varphi(n)$: если $x$ — решение, то и $x+\varphi(n)$ — решение. Значит, достаточно перебрать $x$ от $0$ до $\varphi(n)-1$, но это медленно.

Применим метод встречи посередине. Положим $k=\lceil\sqrt n\rceil$ и запишем $$ x=x_1+k\,x_2,\qquad 0\le x_1\lt k,\quad 0\le x_2\le n/k. $$ Тогда $$ a^{x_1+kx_2}\equiv b \iff a^{x_1}\equiv b\cdot\bigl(a^{-1}\bigr)^{k x_2}\pmod n. $$ Поскольку $a$ обратим, $a^{-k}$ существует. Считаем все «маленькие шаги» $a^{x_1}$ и складываем в хэш-таблицу (значение → $x_1$). Затем «гигантскими шагами» перебираем $x_2=0,1,2,\dots$, вычисляя $b\cdot(a^{-k})^{x_2}$, и ищем совпадение. Первое совпадение даёт $x$.

// a^x = b (mod n), gcd(a, n) = 1
long long dlog(long long a, long long b, long long n) {
    long long k = 1;
    while (k * k < n) ++k;                        // k = ceil(sqrt(n))
    unordered_map<long long,long long> baby;
    long long cur = 1;
    for (long long x1 = 0; x1 < k; ++x1) {        // маленькие шаги
        baby.emplace(cur, x1);
        cur = mulmod(cur, a, n);                  // mulmod — см. §4, иначе переполнение при больших n
    }
    long long invA = inverse(a, n);
    long long giantStep = power(invA, k, n);      // a^{-k}
    long long g = b;                              // b * (a^{-k})^{x2}
    for (long long x2 = 0; x2 <= k; ++x2) {
        auto it = baby.find(g);
        if (it != baby.end()) return it->second + k * x2;
        g = mulmod(g, giantStep, n);
    }
    return -1;                                    // решения нет
}

8. Первообразный корень 53:20

Первообразный корень существует по любому простому модулю (а также по модулям $2$, $4$, $p^{\alpha}$, $2p^{\alpha}$ при нечётном простом $p$). Чтобы его найти, раскладываем $p-1$ на простые множители. Тогда

Почему этого достаточно: порядок $g$ делит $p-1$. Если порядок меньше $p-1$, то он делит $(p-1)/q$ для некоторого простого $q$, и тогда соответствующая проверка провалится. Если же все проверки пройдены, порядок ровно $p-1$.

Первообразных корней много — их ровно $\varphi(p-1)$ штук (в среднем примерно половина остатков). Поэтому можно брать случайное $g$ и проверять критерий; ожидание числа попыток невелико.

Зачем нужно. Если известно, что $a=g^{\alpha}$ и $b=g^{\beta}$, то уравнение $a^{x}\equiv b\pmod p$ превращается в линейное сравнение по модулю $p-1$: $$ g^{\alpha x}\equiv g^{\beta} \iff \alpha x\equiv \beta\pmod{p-1}, $$ которое решается раз и навсегда. Это позволяет обслуживать много запросов дискретного логарифма быстро, заранее один раз найдя первообразный корень.

9. Формула включений-исключений 58:57

Нужно посчитать размер объединения множеств $A_1,\dots,A_n$, зная размеры пересечений любых наборов. Формула включений-исключений утверждает: $$ \Bigl|\bigcup_{i=1}^{n} A_i\Bigr| \;=\; \sum_{\emptyset\ne S\subseteq\{1,\dots,n\}} (-1)^{|S|+1}\Bigl|\bigcap_{i\in S}A_i\Bigr|, $$ то есть суммируем размеры всех пересечений, чередуя знаки по размеру набора.

Пример: число беспорядков

Беспорядок (derangement) — перестановка $\pi$ длины $n$, где $\pi(i)\ne i$ для всех $i$. Пусть $A_i$ — множество перестановок с $\pi(i)=i$. Тогда $|A_i|=(n-1)!$, а размер пересечения набора из $s$ множеств равен $(n-s)!$. Беспорядки — это дополнение объединения, поэтому $$ !n = n!-\sum_{s=1}^{n}(-1)^{s+1}\binom{n}{s}(n-s)! = n!\sum_{k=0}^{n}\frac{(-1)^{k}}{k!}. $$

10. Проверка простоты: Миллер–Рабин и Поллард 1:18:29

Тест Ферма и числа Кармайкла

Если $n$ простое и $\gcd(x,n)=1$, то по малой теореме Ферма $x^{n-1}\equiv 1\pmod n$. Значит, если для какого-то $x$ это не выполняется — число составное. Но обратное неверно: существуют числа Кармайкла — составные $n$, для которых сравнение выполнено для всех $x$, взаимно простых с $n$ (эквивалентно: $n$ бесквадратно и для каждого простого $p\mid n$ выполнено $(p-1)\mid(n-1)$). Теста Ферма недостаточно.

Сильный тест (Миллер–Рабин)

Наложим более сильное условие, обязанное выполняться для простых. Запишем $n-1=2^{k}\cdot a$ с нечётным $a$. Тогда разность раскладывается в произведение: $$ x^{n-1}-1=(x^{a}-1)\prod_{i=0}^{k-1}\bigl(x^{a\cdot 2^{i}}+1\bigr). $$ Если $n$ простое, произведение по модулю $n$ равно нулю, а в поле $\mathbb Z/n\mathbb Z$ нет делителей нуля, значит, хотя бы один множитель равен нулю. Отсюда проверка:

bool isPrime(long long n) {
    if (n < 2) return false;
    for (long long p : {2,3,5,7,11,13,17,19,23,29,31,37,41})
        if (n % p == 0) return n == p;
    long long a = n - 1, k = 0;
    while (!(a & 1)) { a >>= 1; ++k; }             // n-1 = 2^k * a
    for (long long base : {2,3,5,7,11,13,17,19,23,29,31,37,41}) {
        long long y = power(base, a, n);
        if (y == 1 || y == n - 1) continue;
        bool ok = false;
        for (long long i = 0; i < k - 1; ++i) {
            y = mulmod(y, y, n);                    // не y*y%n — переполнение при n ~ 10^18
            if (y == n - 1) { ok = true; break; }
        }
        if (!ok) return false;
    }
    return true;
}

Факторизация: метод Полларда ($\rho$)

Пусть $n=pq$ с близкими $p\approx q$ — это самый трудный для разложения случай. Выбираем функцию по модулю $n$, обычно $f(x)=x^2+1\bmod n$, и строим последовательность $a_0$ (случайно), $a_{i+1}=f(a_i)$. Смотрим на $\gcd(a_i-a_j,\,n)$: он может быть $1$, $p$, $q$ или $n$. Угадали делитель — разложили. Используется алгоритм «черепахи и зайца» (Флойд): один указатель шагает на 1, другой на 2; их разность даёт нужное сравнение по модулю $p$.

long long pollard(long long n) {
    if (!(n & 1)) return 2;
    long long x = 2, y = 2, d = 1;
    while (d == 1) {
        x = (mulmod(x, x, n) + 1) % n;               // черепаха: 1 шаг
        y = (mulmod(y, y, n) + 1) % n;               // заяц: 2 шага
        y = (mulmod(y, y, n) + 1) % n;
        d = std::gcd(std::llabs(x - y), n);
    }
    return d;   // 1 < d < n — делитель; d == n — нужен перезапуск
}
// затем рекурсивно разложить найденный делитель и частное

11. Floor-сумма за $O(\log)$ 1:39:55

Требуется научиться быстро считать сумму $$ S(n,m,a,b)=\sum_{x=0}^{n-1}\left\lfloor\frac{a x+b}{m}\right\rfloor. $$ Полезно переписать её геометрически. Слагаемое $\lfloor(ax+b)/m\rfloor$ равно числу целых $y\ge 1$ с условием $my\le ax+b$. Значит, $S$ — это количество целых точек $(x,y)$ в полосе $0\le x\lt n$, лежащих не выше прямой $y=(ax+b)/m$.

В нормализованном виде ($0\le a,b\lt m$) прямая начинается в точке $(0,b)$, а к $x=n$ поднимается до уровня $Y=\bigl\lfloor\frac{an+b}{m}\bigr\rfloor$. Если $Y=0$, под прямой нет целых точек, и сумма равна нулю. Иначе поворачиваем картинку и суммируем по уровням $y$: точки на высоте $y$ существуют при $x\ge\lceil(my-b)/a\rceil$, а таких $x$ ровно $n-\lceil(my-b)/a\rceil$. Переходя от $\lceil\cdot\rceil$ к $\lfloor\cdot\rfloor$ и делая замену $y'=y-1$, получаем рекурсию $$ S(n,m,a,b)=(n-1)Y-S\bigl(Y,\ a,\ m,\ m-b-1\bigr),\qquad Y=\left\lfloor\frac{an+b}{m}\right\rfloor. $$ Модуль на каждом шаге уменьшается как в алгоритме Евклида, поэтому глубина рекурсии $O(\log(a+m))$.

// sum_{x=0}^{n-1} floor((a*x + b) / m), все аргументы неотрицательны
long long floor_sum(long long n, long long m, long long a, long long b) {
    long long ans = 0;
    if (a >= m) { ans += (long long)((__int128)(n - 1) * n / 2 * (a / m)); a %= m; }
    if (b >= m) { ans += n * (b / m);                b %= m; }
    long long Y = (long long)(((__int128)a * n + b) / m);
    if (Y == 0) return ans;
    // поворот осей: та же сумма, но модуль теперь a
    ans += (long long)((__int128)(n - 1) * Y) - floor_sum(Y, a, m, m - b - 1);
    return ans;
}

12. Минимум floor-суммы 1:57:44

Родственная задача — найти минимум линейной функции по модулю: $$ g(n,m,a,b)=\min_{0\le x\le n}\Bigl((a x+b)\bmod m\Bigr). $$ Она сводится к floor-сумме (через подсчёт числа $x$ с остатком в заданном отрезке), но у неё есть и красивое прямое рекурсивное решение.

Геометрический смысл второго (интуитивного) алгоритма из лекции: «прыгаем» с шагом $a$ по модулю $m$ и интересуемся только первыми точками, попавшими в каждый блок длины $a$. Остаток, отвечающий первой точке блока, каждый раз смещается на фиксированную величину по модулю $a$; задача снова переходит в себя с заменой модуля $m$ на $a$. Это тот же самый обмен ролями «число/модуль», что и в алгоритме Евклида.

13. Семинар: $\varphi$, суммы НОД и формула Мёбиуса 2:05:16

13.1. Количество чисел, взаимно простых с $n$ 2:05:16

Задача: посчитать количество чисел от $1$ до $n$, взаимно простых с $n$ (то есть $\gcd=1$). Ответ — это $\varphi(n)$, но полезно уметь выводить формулу через включения-исключения.

Разложим $n$ на различные простые делители $p_1,\dots,p_k$. Число не взаимно просто с $n$ тогда и только тогда, когда делится хотя бы на одно $p_i$. Для набора $S$ индексов чисел, делящихся на $\prod_{i\in S}p_i$, ровно $\lfloor n/\prod_{i\in S}p_i\rfloor$. Тогда $$ \varphi(n)=\sum_{S\subseteq\{1,\dots,k\}}(-1)^{|S|}\left\lfloor\frac{n}{\prod_{i\in S}p_i}\right\rfloor. $$ Работает за $O(2^{k})$, где $k=\omega(n)$ — число различных простых делителей; на практике это очень мало.

13.2. Сумма НОД по всем подотрезкам за $O(n\log C)$ 2:09:28

Дан массив, нужно просуммировать $\gcd$ по всем подотрезкам; $C$ — максимальное значение. Зафиксируем правую границу $r$ и будем смотреть на $\gcd(l\dots r)$ при уменьшении $l$. Ключевой факт: при движении $l$ налево $\gcd$ может только уменьшаться и каждый раз делится на предыдущее значение, поэтому различных значений $\gcd$ на отрезке $O(\log C)$. Храним список пар «значение $\gcd$ — количество левых границ с таким значением».

При сдвиге правой границы на $r+1$ добавляем новое число $a_{r+1}$: у всех старых значений берём $\gcd$ с $a_{r+1}$, приписываем новый отрезок из одной границы и склеиваем соседние равные значения. Суммарно по всем $r$ получается $O(n\log C)$.

// sum of gcd over all subsegments, O(n log C)
long long ans = 0;
vector<pair<int,long long>> cur;      // {gcd, number of left borders}
for (int r = 0; r < n; ++r) {
    vector<pair<int,long long>> nxt;
    nxt.push_back({a[r], 1});
    for (auto [g, cnt] : cur) {
        int ng = std::gcd(g, a[r]);
        if (ng == nxt.back().first) nxt.back().second += cnt;
        else                        nxt.push_back({ng, cnt});
    }
    cur = nxt;
    for (auto [g, cnt] : cur) ans += (long long)g * cnt;
}

13.3. Сумма НОД по всем подмножествам за $O(C\log C)$ 2:16:11

Дан набор чисел, нужно просуммировать $\gcd$ по всем непустым подмножествам.

1. Пусть $g[d]$ — сколько чисел набора делится на $d$. Считаем за $O(C\log C)$: для каждого $d$ перебираем его кратные и суммируем счётчики.
2. Число непустых подмножеств с $\gcd$, делящимся на $d$, равно $2^{g[d]}-1$.
3. Пусть $f[d]$ — число подмножеств с $\gcd$, равным ровно $d$. Тогда $$ 2^{g[d]}-1=\sum_{d\mid m} f[m]. $$ Преобразование (обращение) выполняется сверху вниз: $f[d]=\bigl(2^{g[d]}-1\bigr)-\sum_{m=2d,3d,\dots} f[m]$.
4. Ответ: $\sum_d d\cdot f[d]$.

Такую операцию «обращения по делителям» и называют формулой (функцией) Мёбиуса; по сути это то же включение-исключение по простым делителям.

14. Семинар: решёта, линейное решето и мультипликативные функции 2:22:49

14.1. Простые до $n$ при $O(\sqrt n)$ памяти 2:22:49

Задача: посчитать количество и сумму простых от $1$ до $n$ за $O(n\log\log n)$ времени и $O(\sqrt n)$ памяти. Прямое решето требует $O(n)$ памяти, поэтому используем блочное (сегментное) решето.

1. Решетом Эратосфена до $\sqrt n$ находим все простые $p\le\sqrt n$ — их $O(\sqrt n)$.
2. Любое составное $m\le n$ имеет простой делитель $\le\sqrt n$. Поэтому, идя по блокам длины $\sqrt n$, вычёркиваем в каждом блоке числа, кратные найденным простым.
3. Для каждого простого храним последнее посещённое кратное, чтобы не пересчитывать левую границу блока делением.

Суммарная работа $\sum_{p\le\sqrt n} n/p = \Theta(n\log\log n)$ (гармонический ряд по простым сходится к $\log\log n$).

14.2. Линейное решето (минимальный простой делитель) 2:30:31

Хотим за честное $O(n)$ найти для каждого числа его минимальный простой делитель $\mathrm{spf}(i)$. Идём по $i$ от $2$ до $n$; если $\mathrm{spf}(i)$ ещё не назначен, то $i$ простое. Затем для каждого простого $p\le\mathrm{spf}(i)$ помечаем $i\cdot p$ как имеющее $\mathrm{spf}=p$. Каждое число получает делитель ровно один раз, поэтому суммарно $O(n)$.

vector<int> spf(n + 1, 0), primes;
for (int i = 2; i <= n; ++i) {
    if (spf[i] == 0) { spf[i] = i; primes.push_back(i); }
    for (int p : primes) {
        if (p > spf[i] || (long long)i * p > n) break;
        spf[i * p] = p;                 // минимальный простой делитель
    }
}

14.3. Мультипликативные функции на линейном решете 2:34:32

Зная $\mathrm{spf}(i)$, легко считать мультипликативные функции для всех чисел.

• Число делителей $\tau(i)$. Храним $e(i)$ — степень, в которой $\mathrm{spf}(i)$ входит в $i$. Если $p^2\mid i$, то $e(i)=e(i/p)+1$ и $\tau(i)=\tau(i/p)\cdot\frac{e(i)+1}{e(i/p)+1}$; иначе $e(i)=1$ и $\tau(i)=2\tau(i/p)$.
• Сумма делителей $\sigma(i)$. Аналогично, домножаем на геометрическую прогрессию по степени $p$.
• Функция Эйлера $\varphi(i)$. Если $p^2\mid i$, то $\varphi(i)=p\,\varphi(i/p)$; иначе $\varphi(i)=(p-1)\,\varphi(i/p)$, где $p=\mathrm{spf}(i)$.

Простейшая реализация «в лоб» тоже быстра: для каждого $i$ делим на $\mathrm{spf}(i)$, пока делится, — суммарно это почти линейно, так как степень вхождения минимального простого делителя мала.

14.4. Сумма $b_k=\sum_{k\mid i} a_i$ 2:41:22

Нужно для каждого $k$ посчитать сумму элементов на позициях, кратных $k$. Наивно (для каждого $k$ пройти по его кратным) это $O(n\log n)$ — гармонический ряд. Есть более экономный вариант: прогонять преобразование отдельно для каждого простого $p$, идя по массиву с шагом $p$ и прибавляя $a[i]$ к $a[i/p]$; суммарно $n/2+n/3+n/5+\dots=\Theta(n\log\log n)$.

// b[k] = sum of a[i] for i divisible by k
vector<long long> b(a.begin(), a.end());
for (int p : primes)                       // все простые до n
    for (int i = n / p; i >= 1; --i)
        b[i] += b[i * p];                  // накапливаем кратные
// теперь b[k] содержит сумму по всем i, кратным k

15. Семинар: НОД на отрезке, GCD-преобразование и трюк с $\varphi$ 2:37:44

15.1. Прибавление на отрезке и НОД на отрезке 2:37:44

Требуется поддерживать прибавление числа ко всем элементам отрезка и запрос $\gcd$ на отрезке. Ключевой факт: $$ \gcd(a_l,\dots,a_r)=\gcd\bigl(a_l,\ a_{l+1}-a_l,\ \dots,\ a_r-a_{r-1}\bigr). $$ Перейдём к разностному массиву $d_i=a_i-a_{i-1}$. Прибавление на отрезке $[l,r]$ меняет только два элемента разностного массива: $d_l$ и $d_{r+1}$ (если $r\lt n$).

• На разностном массиве строим дерево отрезков, поддерживающее $\gcd$ и точечные изменения.
• Значение $a_l$ получаем из префиксной суммы разностей (Fenwick или дерево).
• Запрос: $\gcd(a_l,\ \gcd(d_{l+1},\dots,d_r))$.

Модификация — за $O(\log n)$, запрос — за $O(\log n)$. Более того, в цепочке пересчёта $\gcd$ каждый раз уменьшается, делясь на предыдущее значение, поэтому операции с НОД в дереве на практике очень быстрые.

// a_l через Fenwick по разностям, gcd по разностям — деревом отрезков
// range_add(l, r, x): add(l, x); if (r+1 < n) add(r+1, -x);
// range_gcd(l, r): return gcd(point_query(l), seg_gcd(l+1, r));

15.2. Сумма попарных НОД и трюк с $\varphi$ 2:47:22

Нужно для каждого $i$ посчитать сумму $\sum_{j\lt i}\gcd(a_i,a_j)$, а затем просуммировать по всем парам.

Используем тождество $$ \sum_{d\mid n}\varphi(d)=n, $$ то есть $\gcd(x,y)=\sum_{d\mid x,\ d\mid y}\varphi(d)$. Тогда суммирование по всем парам даёт $$ \sum_{\text{pairs}}\gcd(x,y)=\sum_{d}\varphi(d)\cdot\binom{c_d}{2}, $$ где $c_d$ — число элементов массива, делящихся на $d$. Можно считать и «на лету»: поддерживать $c_d$ и для нового $x$ добавлять $\sum_{d\mid x}\varphi(d)\cdot c_d$.

15.3. Преобразование по делителям в задаче о попарных НОД 2:49:07

Тот же аппарат используется, чтобы для каждого нового числа $x$ быстро узнать «сколько уже добавленных чисел взаимно просты с $x$»: поддерживаем $c_d$ — количество добавленных чисел, делящихся на $d$, и отвечаем включением-исключением по бесквадратным делителям $x$: $$ \bigl|\{y:\gcd(x,y)=1\}\bigr|=\sum_{d\mid x}(-1)^{\omega(d)}c_d, $$ где суммирование идёт по бесквадратным $d$. Коэффициенты $(-1)^{\omega(d)}$ — это значения функции Мёбиуса.

16. Семинар: включения-исключения в задачах 2:57:21

16.1. «Кубастые» числа 2:57:21

Натуральное число называется кубастым, если его можно представить в виде $a^3 b$ с $a\gt 1$, $b\ge 1$; нужно посчитать количество кубастых чисел до $n$ (для $n\le 10^{18}$).

Число кубасто тогда и только тогда, когда делится на $p^3$ для некоторого простого $p$. Обозначим $A_p$ — множество чисел, делящихся на $p^3$; нужно найти $|\bigcup_p A_p|$. Для набора различных простых $p_1,\dots,p_s$ размер пересечения равен $\lfloor n/(p_1\cdots p_s)^3\rfloor$. По формуле включений-исключений $$ \Bigl|\bigcup_p A_p\Bigr|=\sum_{d\ge 2,\ d\ \text{бесквадратно}}(-1)^{\omega(d)+1}\left\lfloor\frac{n}{d^3}\right\rfloor. $$ Поскольку $d^3\le n$, достаточно $d\le n^{1/3}\le 10^6$; все нужные $d$ можно получить решетом, накопив знак $(-1)^{\omega(d)+1}$.

16.2. Массивы с ограничением сверху 3:07:04

Сколько существует массивов длины $n$ с элементами из отрезка $[0,k]$, сумма которых равна $s$? Без ограничения сверху ответ даёт метод шаров и перегородок: $$ \binom{s+n-1}{n-1}. $$ Ограничение $a_i\le k$ учитываем включением-исключением по позициям, где $a_i\ge k+1$: выбираем $c$ «плохих» позиций, вычитаем из суммы по $(k+1)$ на каждой и распределяем остаток свободно. Итог $$ \sum_{c=0}^{\lfloor s/(k+1)\rfloor}(-1)^{c}\binom{n}{c}\binom{s-(k+1)c+n-1}{n-1}. $$ Это в точности коэффициент при $x^s$ у многочлена $(1+x+\dots+x^k)^n$.

16.3. Пара с максимальным НОК 3:16:48

Лёгкая версия: найти пару с минимальным (а затем с максимальным) НОК. Идея для минимума: перебрать кандидата в $\gcd$, взять два наименьших элемента массива, делящихся на него, и обновить ответ через их НОК; если найденный «НОК» не является настоящим $\gcd$, всегда существует пара не хуже, так что оптимум будет найден.

Сложная версия — максимум — сводится к одной задаче о двух взаимно простых числах. Построим множество $S$ из всевозможных делителей всех $a_i$. Утверждение: максимум произведения двух взаимно простых чисел из $S$ равен максимальному НОК пары исходных чисел.

Оставшуюся задачу решают двумя указателями. Массив сортируют по убыванию; поддерживая структуру, которая умеет удалять число и отвечать «есть ли во множестве число, взаимно простое с данным $x$». Ответ на запрос даётся включением-исключением по бесквадратным делителям $x$: $$ \bigl|\{y\in S:\gcd(x,y)=1\}\bigr|=\sum_{d\mid x,\ d\ \text{бескв.}}(-1)^{\omega(d)}\,c_d, $$ где $c_d$ — число элементов $S$, делящихся на $d$. При удалении элемента «пересчитываем» только счётчики его делителей. Так как каждый элемент имеет немного делителей, суммарно структура работает быстро.