ゼータ変換について勉強する機会になった。
ゼータ変換とメビウス変換
- 下位集合に関するゼータ変換: F(S)=∑T⊆Sf(T)
- 上位集合に関するゼータ変換: F(S)=∑T⊇Sf(T)
これをO(2NN)でやるのが高速ゼータ変換。普通の総和以外にも、×,max,minとかでもまったく同じ方法が使える。結合法則と交換法則が成り立っていればよいので、可換半群なら大丈夫なはず。
- 下位集合に関するメビウス変換: f(S)=∑T⊆S(−1)∣S∖T∣F(T)
- 上位集合に関するメビウス変換: f(S)=∑T⊇S(−1)∣T∖S∣F(T)
ゼータ変換の逆変換。つまり、上のFからfを導出する変換。
ゼータ変換と違って引き算を使うので、アーベル群でないといけない。演算がminとかだとそもそもfを一意に復元できないので、当然ではある。
(用語法はよくわからないけど、競プロ界ではだいたいみんなnaoya_tさんの記事にある使いわけをしているっぽいのでそれに従っている。文献(これとか)によっては、ここでいう下位集合に対するゼータ変換がMöbius transform、その逆変換がMöbius inversionと呼ばれていたりする。)
∪畳み込みと∩畳み込み
f, gを
F(S):=T⊆S∑f(T)
G(S):=T⊆S∑g(T)
と下位集合に関してゼータ変換して、点毎に掛け算したものをHとする:
H(S):=F(S)×G(S)=T1⊆S∑f(T1)×T2⊆S∑g(T2)=T⊆S∑T1∪T2=T∑f(T1)×g(T2)
Hをメビウス変換すると
h(S):=T1∪T2=S∑f(T1)×g(T2)
という関数が得られて、∪に関する畳み込みができたことになる。
同様に、上位集合に関してゼータ変換して掛け合わせて逆変換すると、∩に関する畳み込みができる。
分配法則を使っているので、この掛け算は足し算と合わせて環になっていないといけない。
約数、倍数に関するゼータ変換とメビウス変換
上と同じ話を約数、倍数について繰り返す。
- 約数に関するゼータ変換: F(x)=∑d∣xf(d)
- 倍数に関するゼータ変換: F(x)=∑x∣mf(m)
ただし、定義域は{0,1,...,N}とかに決まっているとして、その中だけで計算する。この変換はO(NloglogN)でできる。
上のFからfを導出する変換。同じオーダーでできる。歴史的には約数・倍数に関するもののほうが最初に考えられたらしい。
LCM畳み込みとGCD畳み込み
やはり上と同じ話をLCMとGCDについて書き下してみる。
f, gを
F(x):=d∣x∑f(d)
G(x):=d∣x∑g(d)
と約数に関してゼータ変換して、掛け合わせたものをHとする:
H(x):=F(x)×G(x)=d1∣x∑f(d1)×d2∣x∑g(d2)=d∣x∑lcm(d1,d2)=d∑f(d1)×g(d2)
Hを逆変換すると
h(x):=lcm(d1,d2)=x∑f(d1)×g(d2)
という関数が得られて、lcmに関する畳み込みができたことになる。
同様に、倍数に関してゼータ変換して掛け合わせて逆変換すると、gcdに関する畳み込みができる。
問題
lcm(Ai,Aj)=AiAj/gcd(Ai,Aj)で、gcd(Ai,Aj)は高々1000000なので、同じgcd毎にまとめて和を計算できれば良さそう。
最終的には、dp[x]がgcd(Ai,Aj)=xとなるような組み合わせについてのAiAjの総和になってほしいので、とりあえずdpを次のように初期化していく。
dp[Ai]:=dp[Ai]+Ai
dpを倍数に関してゼータ変換してdp[x]:=dp[x]⋅dp[x]と掛け合わせて逆変換すると、
dp[x]=gcd(Ai,Aj)=x∑AiAj
になっている。このままでは不必要な組み合わせについても加算されているので、
dp[Ai]:=dp[Ai]−Ai2
として添え字が同じ組み合わせを除いたあと、AiAjとAjAiのダブルカウントを除くためにすべて半分にする: dp[x]:=dp[x]/2。(ここではなくて一番最後に半分にしてもよい。)
これでdp[x]はgcd(Ai,Aj)=xとなるような組み合わせについてのAiAj(i<j)の総和になったので、あとはgcdで割りながら総和を取れば答えになる: ∑x=11000000dp[x]/x