sum of multiplicative function

これはなに？

Min_25さんの以下の記事を自分が読みやすいように書き換えたものです。
min-25.hatenablog.com
Min_25さんによる実装はこちら
https://ideone.com/9F9amv

要求知識

・平方根の上下で分けるテク
・fenwick tree
・線形篩の"気持ち"を知っていると理解の助けになる箇所が一部ある

定義

定義域が正整数全体である複素数値関数 f が次を満たすとき、f は乗法的(multiplicative)であるという。
$\forall a,b\ ({\rm gcd} (a,b)=1 \Rightarrow f(ab)=f(a)f(b))$

定義から、f≠0 が multiplicativeなら f(1)=1 であることに注意。以下単に f とかいたとき、それは multiplicative であるとする。

以下では値域が複素数であることは使っていないので、演算がO(1)でできる環なら同様に計算できるはず。

記号の定義

・ $\displaystyle{ F_{{\rm prime}}(x):=\sum_{p\leq x,p:{\rm prime}}f(p)}$ 　素数のとこだけの和
・ $F(x):=\sum_{i=1}^{x} f(i)$ 　全部の和
・ ${\rm lpf}(x)$ ： x の最小素因子。 ${\rm lpf}(1)=\infty$
・ $p_k$ ： k番目の素数。1-based index
・ $\pi(x)$ ： x以下の素数の個数
・ $f_k(x):= \begin{cases} f(x) & ({\rm lpf}(x) \geq p_k) \\ 0 & (\text{otherwise}) \end{cases}$ 　最小素因子がでかいとこだけ生き残るやつ
・ $F_k(x):=\sum_{i=1}^{x} f_k(i)$ 　最小素因子がでかいとこだけの和。 $F_1=F$ に注意。
・ $H_N:=\{\lfloor \frac{N}{1}\rfloor, \lfloor \frac{N}{2}\rfloor,\ldots,\lfloor \frac{N}{\lfloor\sqrt{N}\rfloor}\rfloor\}$
・ $L_N:=\{1,2,\ldots,\sqrt{N}\}$
$H_N,L_N$ については N が明らかな時には単にH,Lと表す。

何ができるか

任意の素数べき $q=p^e$ に対し $f(q)$ は $q,p,e$ から $O(1)$ で計算可能であるとする。
全ての $x\in H_N\cup L_N$ に対して $F_{{\rm prime}}(x)$ が求まっているとき、全ての $x\in H_N\cup L_N$ に対する $F(x)$ を時間計算量 $O(N^{2/3})$ 、空間計算量 $O(N^{1/2})$ で求めることができる。

Q. 元記事に比べて仮定が強すぎない？（元記事では $F_{{\rm prime}}(x)$ の計算からやっているはず）
A. $F_{{\rm prime}}(x)$ の計算には本編と同じようなテクを使うので、説明の順番的に後回し

$F_k$ という概念を用意したことからお察しのとおり、kが大きいところから順に $F_k$ を計算していって、最終的に $F_1$ にたどり着きます。
前処理として $\sqrt{N}$ 以下の素数を列挙しておきます。
本編は3つのパートからなります。

お断り：この記事は気持ちを理解することを優先しているため、境界条件については適当です。正確な部分は元記事を見るか自分でよく考えて実装しましょう。

$x\in H\cup L$ に対して $F_{\pi(N^{1/3})+1}(x)$ を求める

要するに、 $N^{1/3}$ より大きな素因子しか持たない部分について計算しましょうという話です。 $k:=\pi(N^{1/3})+1$ 、 $p:=p_k$ とします。
${\rm lpf}(i) \geq p$ を満たすような i は「1」「素数」「素数の2乗」「相異なる2つの素数の積」のいずれかです。ということで、 $F_k$ は $F_{{\rm prime}}$ を使うことで簡単に計算できます。
具体的には x の値に応じて次の通りです。
・ $1\leq x < p$ のとき
$i\leq x$ の範囲には ${\rm lpf}(i)\geq p$ となる i は「1」しか存在しません。
よって $F_k(x)=1$ です。
・ $p\leq x < p^2$ のとき
$i\leq x$ の範囲には ${\rm lpf}(i)\geq p$ となる i は「1」「素数」しか存在しません。
よって $F_k(x)=1+(F_{{\rm prime}}(x)-F_{{\rm prime}}(p-1))$ です。
・ $p^2\leq x \leq N$ のとき
$i\leq x$ の範囲には ${\rm lpf}(i)\geq p$ となる i は「1」「素数」「素数の2乗」「相異なる2つの素数の積」しか存在しません。
よって $F_k(x)=1+(F_{{\rm prime}}(x)-F_{{\rm prime}}(p-1)) + \sum_{p\leq p'\leq \sqrt{x}}f(p'^2) + \sum_{p\leq p'\leq \sqrt{x}}f(p')(F_{{\rm prime}}(x/p')-F_{{\rm prime}}(p'))$ です。
第4項は「p'×(p'よりでかい素数p'')」の形で表される数について、p''で和を取ったものです。

計算量解析
上2つのケースは各xに対しO(1)です。そもそもxの候補が全体で $O(\sqrt{N})$ 個しかないので、この部分は合わせて $O(\sqrt{N})$ です。
3つ目のケースは各xに対し $O(\pi(\sqrt{x}))$ です。xとして $\lfloor \frac{N}{1}\rfloor, \lfloor \frac{N}{2}\rfloor,\ldots,\lfloor \frac{N}{\lfloor\sqrt[3]{N}\rfloor}\rfloor$ を計算するので、
$O(\sum_{i=1}^{N^{1/3}}\pi(\sqrt{N/i})) = O(\int_1^{N^{1/3}} (N/i)^{1/2} / \log (N^{2/3}) di) = O(N^{2/3}/\log N)$ です。
空間計算量は、 $F_k$ の値の保持以外に領域を使っていないため、明らかに $O(N^{1/2})$ です。

$x\in H\cup L$ に対して $F_{\pi(N^{1/6})+1}(x)$ を求める

ここが一番難しい。
さっきは $N^{1/3}$ まででしたが今度は $N^{1/6}$ までやります。 $F_k$ から $F_{k-1}$ を作るというステップを順番に刻んでいて $N^{1/6}$ まで行きます。更新はin-placeにやることにします。
説明のためだけに新たな記号を用意します（えぇ……）。
$x\in H\cup L$ に対し、 $x':=\max\{y \mid y\in H\cup L, y< x\}$ とします。 $H\cup L$ に入ってるやつで x の次に小さい値ですね。ただし、 $x=1$ のとき $x'=0$ とします。
$x$ が $N^{2/3}$ 以上か以下かで更新の仕方を変えるんですが、そもそも $N^{2/3}$ 以上か以下かで値の持ち方も変えます。
つまり、今までは、各 $x\in H\cup L$ に対して $F_k(x)$ を持っていましたがそれをやめて、 $x< N^{2/3}$ のときには、一旦 $F_k(x)-F_k(x')$ を持つことにします。 $x\geq N^{2/3}$ のときは今まで通り $F_k(x)$ を持ちます。
そうして、 $x< N^{2/3}$ の部分をfenwick treeで管理します。これで今まで通り $F_k(x)$ を取ってこれます。
要するに、「点更新クエリが飛んでくるなら累積和じゃなくてfenwick treeかなーやっぱりｗ」です。

更新の仕方を考えましょう
・ $N^{2/3}\leq x \leq N$ のとき
「もらうDP」の気持ちで更新します。
$F_k(x)=\sum_{e=0}^{\infty} f(p_k^e) F_{k+1}(x/p_k^e)$ を使って愚直に計算するだけです。
$x/p_k^e< N^{2/3}$ のときの $F_{k+1}(x/p_k^e)$ の値はfenwick treeにrange sumクエリを投げて取ってくることになるのでlogがつくことに気をつけましょう。
・ $1\leq x < N^{2/3}$ のとき
「配るDP」の気持ちで更新します。
まずは x=1から初めて $xp_k,xp_k^2,\ldots$ に配ります。そうしたら、配った先の地点 $y$ で、今度は更に $yp_{k+1},yp_{k+1}^2,\ldots,yp_{k+2},yp_{k+2}^2,\ldots,$ に配ります。線形篩と同じですね。
配る時にはfenwick treeにpoint addクエリを投げることになるのでlogがつくことに気をつけましょう。
$N^{2/3}$ 以上の範囲は上述の通り別に計算しているので、配る範囲は $N^{2/3}$ 未満までであることにも気をつけましょう。
配り元としては $N^{2/3}/ p_k < N^{1/2}$ までを考えればいいので、 $f_{k+1}(x)=F_{k+1}(x)-F_{k+1}(x')$ の値が問題なくとれることにも注意しましょう。
初手で配り先は $p_k$ 以上になるので、2個目以降の素数としては、 $N^{2/3}/p_k\leq N^{1/2}$ 以下のみを考慮すれば良く、事前に列挙した範囲で収まることにも気をつけましょう。
気をつけることがいっぱい。

計算量解析
ここも難しい
・ $N^{2/3}\leq x \leq N$ の範囲について
面倒なのでfenwick treeにクエリを投げなくていいところも投げてしまうことにすると、クエリを投げる回数の総和は
$O(\sum_{N^{1/6}< p\leq N^{1/3}, p:{\rm prime}}\sum_{i=1}^{N^{1/3}}\log_p(N/i))$ となります。
$\sum_{i=1}^{N^{1/3}}\log_p(N/i) =O(N^{1/3}\log N/\log p)=O(N^{1/3})$ に気をつけて変形すると $O(N^{1/3}/\log N \cdot N^{1/3})$ となるので、クエリ1回に $O(\log N)$ かかることから、全体の計算量は $O(N^{2/3})$ です。
・ $1\leq x < N^{2/3}$ の範囲について
配る先がいくつあるかがクエリを投げる回数になります。
配る先は、lpf が $N^{1/6}$ 以上であるような $N^{2/3}$ 未満の数全てです。そのような数は $O(N^{2/3}/\log N)$ 個らしいので*1、クエリ1回に $O(\log N)$ かかることから、全体の計算量は $O(N^{2/3})$ です。
$F_k$ の計算時には $F_k,F_{k+1}$ の値の一部をfenwick treeで持っていますが、定数倍の違いしかないので空間計算量は $O(N^{1/2})$ です。

$x\in H\cup L$ に対して $F(x)$ を求める

$F_k(x)=\sum_{e=0}^{\infty} f(p_k^e) F_{k+1}(x/p_k^e)$ を使って愚直に計算するだけです。やるだけ。

計算量解析
k を1つ進めるために、各 xに対し、 $\log_{p_k}(x)$ 回の計算が必要なので、全体で
$\sum_{p\leq N^{1/6}, p:{\rm prime}} \sum_{i=1}^{N^{1/2}}(\log_p(N/i)+\log_p(i))$ くらいの更新が必要で
$O(\sum_{2\leq p\leq N^{1/6}} \sum_{i=1}^{N^{1/2}}\log_p(N))=O(N^{1/2}\log N \sum_{2 \leq p\leq N^{1/6}} 1/\log p) =O(N^{2/3})$ になる。
$F_k$ の計算時には $F_k,F_{k+1}$ 以外の値を保持していないので空間計算量は $O(N^{1/2})$ です。

結局 $F_{{\rm prime}}$ は？

TODO：書く

2<=p<=N^{1/6}は (Lucy) DP で O(N^{2/3}/log(N))。
N^{1/6}<p<=N^{1/3} は、 N^{2/3} 未満は Fenwick で愚直更新で O(N^{2/3}/log(N)*log(N))、N^{2/3} 以上は DP + Fenwick で O(N^{1/3}/log(N)*N^{1/3}*log(N))。
N^{1/3}<p<=N^{1/2} は、Lucy DP で O(N^{2/3}) (N^{2/3}/log(N)) (?) 想定でした。
— Min_25 (@min_25_) 2020年4月7日

Lucy DPというのは $O(N^{3/4})$ prime counting でおなじみのやつです（多分）

*1:Buchstab function - Wikipedia

メモ

yukicoderでゆるふわgolf