Lucy’s DP

包除原理より \(n\) 以下の数のうち、\(p_1, p_2, \ldots, p_j\) と互いに素なものの個数 \(\pi_j(n)\) は
\begin{align}
\pi_j(n)=\sum \mu(d) \lfloor n/d \rfloor \\
\end{align}と書ける。ただし、\(\mu\) は \(p_1, p_2, \ldots, p_j\) を素因数に持つ数からなる整除束の部分束のメビウス関数とする。

\(n\) を固定する。\(\lfloor n / i\rfloor\) を基底とする整数係数の加群 \(V\) を考えて \(\boldsymbol{v} = \sum a_i \lfloor n / i\rfloor\) に対する作用素 \(D_m\) を \(D_m \boldsymbol{v} = \sum a_i \lfloor n / (im)\rfloor\) と定義する。\(\mathbb{N}\) と \(V\) の元を区別して、\(\sum a_i \lfloor n / i\rfloor\) は \(\sum a_i\boldsymbol{e}_{i}\) のように書くべきかもしれないが、分かりやすいので混同する。
\(i\) 番目の素数 \(p_i\) を用いて \(f_i = 1-D_{p_i}\) とすると $$\pi_j(n) = \left(\prod_{i=1}^j f_i \right) \left\lfloor \frac{n}{1} \right\rfloor $$ になる（左辺と右辺で型が違いますが、自然に解釈してください）。合成数 \(m\) の素因数には \(\sqrt{m}\) 以下のものが必ず存在するから、\(p_i^2 \leq n\) となる最大の \(i\) を \(k\) と置くと$$\pi(n) = k+\pi_k(n)$$である。

基底が小さい場合に帰着したいので \(i\) の降順に \(f_i\) を作用させて
$$\pi_k(n) = f_1\left(f_2\left(\cdots \left(f_{k-1}\left(f_k\left\lfloor \frac{n}{1} \right\rfloor\right)\right)\right)\right) $$の順に計算する。\(\left\lfloor \frac{n}{d} \right\rfloor\) の形の基底は \(O(\sqrt n)\) 個しか存在しない。また、\(\left\lfloor \frac{n}{d}\right\rfloor < p_m^2\) ならば \(f_m\) は恒等作用素として扱ってよい。つまり、dp の更新は in-place に行い \(\sqrt {\left\lfloor \frac{n}{d}\right\rfloor} < p_m\) となる基底には何もしない。

• \(d \leq \sqrt{n}\) となる基底 \(\left\lfloor \frac{n}{d}\right\rfloor\) に対して、恒等でない作用が \(O(\sqrt{n/d})\) 回起きる。作用の回数の和は \(O\left(\int_{x=1}^{\sqrt{n}} \sqrt{n/x}dx \right) = O(n^{3/4})\)。
• \(d \geq \sqrt{n}\) となる基底 \(\left\lfloor \frac{n}{d}\right\rfloor\) に対して、恒等でない作用が \(O(\sqrt{n/d}) = O(n^{1/4})\) 回起きる。このような基底は \(O(\sqrt{n})\) 個しかないので、作用の回数の和は \(O(n^{3/4})\)。

Meissel-Lehmer algorithm

基本方針は \(f_k, f_{k-1}, \ldots, \) を順に作用させて、\(\pi_k(n)\) を \(\pi_i(j)\) \((1 \leq j \leq n^{2/3})\) と \(\pi_0(j)\) \((1 \leq j \leq n)\) の線形結合で表すことである。
\begin{align}
\pi(n) &= k+f_1f_2\cdots f_{k-1}\left(f_k\left\lfloor \frac{n}{1} \right\rfloor\right) \\
&= k+f_1f_2\cdots f_{k-1}\left(\left\lfloor \frac{n}{1} \right\rfloor – \left\lfloor \frac{n}{p_k} \right\rfloor\right) \\
&= k+f_1f_2\cdots f_{k-1}\left\lfloor \frac{n}{1} \right\rfloor – \pi\left(\left\lfloor \frac{n}{p_k} \right\rfloor\right)\\
&\vdots\\
&= k+f_1f_2\cdots f_{l}\left\lfloor \frac{n}{1} \right\rfloor – \sum_{i=l+1}^k \pi\left(\left\lfloor \frac{n}{p_i} \right\rfloor\right)\\
\end{align}を得る。ただし、\(p_i^2 \leq n\) となる最大の \(i\) を \(k\)、\(p_i^3 \leq n\) となる最大の \(i\) を \(l\) と置く。
\(\sum_{i=l+1}^k \pi\left(\left\lfloor \frac{n}{p_i} \right\rfloor\right)\) は \(n^{2/3}\) 以下の数に対してエラトステネスの篩を適用すれば計算できる。\(\left\lfloor \frac{n}{d} \right\rfloor \leq n^{2/3}\) となる基底に対しては、それ以上作用素を掛けずそのままにして、計算を進めると、$$f_1f_2\cdots f_{l}\left\lfloor \frac{n}{1} \right\rfloor = \sum_i \sum_{\left\lfloor \frac{n}{d} \right\rfloor \leq n^{2/3}} a_{i, d} f_1f_2 \cdots f_i \left\lfloor \frac{n}{d}\right\rfloor + \sum_{d < n^{1/3}} b_d \left\lfloor \frac{n}{d}\right\rfloor $$と展開できる。\(\left\lfloor \frac{n}{d} \right\rfloor > n^{2/3}\) となる基底は \(O(n^{1/3})\) 個しかないので、\(O(n^{2/3})\)で展開できる。\(f_1f_2 \cdots f_i \left\lfloor \frac{n}{d}\right\rfloor\)は\(\left\lfloor \frac{n}{d}\right\rfloor\)以下の数のうち、\(p_1, p_2, \ldots, p_i\) で割り切れないものの個数だから、累積和を高速に得られるデータ構造（BITなど）でエラトステネスの篩を管理すれば、\(l = O(n^{1/3} / \log(n))\) より \(O(n^{2/3}\log(n))\) で列挙できる。エラトステネスの篩を \(O(n^{1/3})\) 個ずつに区切って、それぞれ別々に更新すれば、空間計算量 \(O(n^{1/3})\) にできる。

一般化

完全乗法的関数 \(g\) に対して、\(G(n) = \sum_{1 \leq i \leq n} g(i), G'(n) = \sum_{1 \leq i \leq n\\ p \text{ is prime}} g(i)\) と置く。\(f_i = 1-g(p_i)\left\lfloor \frac{\cdot}{p_i}\right\rfloor\)と置くと$$G'(n)=\sum_{i=1}^k g(p_i)+ f_1f_2\cdots f_k G(n)$$が成り立つ。\(g(p_i)=1\) と置くと \(G'(n)=\pi(n)\) である。

\(G(\lfloor n/1\rfloor), G(\lfloor n/2\rfloor), \ldots, G(\lfloor n/n\rfloor)\) が \(O(1)\) で求まるならば、Meissel-Lehmer algorithm により \(O(n^{2/3}\log(n))\) で \(G'(n)\) が求まる。アルゴリズムを少し変更すると \(G'(\lfloor n/1\rfloor), G'(\lfloor n/2\rfloor), \ldots, G'(\lfloor n/n\rfloor)\) が同様の計算量で列挙できる。\(1 \leq n/j \leq n^{2/3}\) に対して、$\pi_i(\lfloor n/j \rfloor)$ の取得および列挙は BIT によって \(O(\log(n))\) でできる。\(i\) の降順に必要になるので、破壊的に変更して構わない。\(n^{2/3} \leq n/j \leq n\) に対する \(\pi_i(\lfloor n/j \rfloor)\) の列挙は、次の様にする。\(1 \leq i \leq \pi(n^{1/3}), 1\leq j \leq n^{1/3}\) なので、全て保持しても \(O(n^{3/2} / \log(n))\) で済む。上記の説明では、\(\lfloor n/i\rfloor\) の線形結合を保持したが、今回は \(\boldsymbol{v}_{\lfloor n/d \rfloor}^{(i)}=f_if_{i+1}f_{i+2} \cdots f_l \lfloor n/d\rfloor\) を（ベクトルではなく値として）各 \(d\) について列挙する。\(\boldsymbol{v}_{\lfloor n/d \rfloor}^{(i-1)} = f_{i-1}\boldsymbol{v}_{\lfloor n/d \rfloor}^{(i)} = \boldsymbol{v}_{\lfloor n/d \rfloor}^{(i)} – f(p_{i-1})\boldsymbol{v}_{\lfloor n/(p_{i-1}d) \rfloor}^{(i)}\) となる。ゼータ変換に対して同様のアルゴリズムを適用することで (G’) から (G) を求めることもできる。すなわち、\(f_i^{-1}=1+g(p_i) \left\lfloor \frac{\cdot}{p_i} \right\rfloor\) とすると \(G(n) = f_1^{-1}f_2^{-1}..f_k^{-1}(G'(n)-\sum_{i=1}^k g(p_i))\) となる。

floorに関する問題

• \(a_{n+1} = \left\lfloor \frac{3}{2}a_n\right\rfloor\)について、\(a_0=L,L+1,\ldots,R\)としたときの\(a_n\)の和 (ARC135F)
• \(a_{n+1} = a_n+1-\left\lceil \frac{1}{n}a_n\right\rceil\)について、\(\sum_{i=L}^Ra_i\) (ARC135E)
• \(\left\lfloor \frac{x}{a}\right\rfloor – \left\lfloor \frac{y}{b}\right\rfloor + c = 0\) の整数解 \((x, y)\) (ARC123E)