前回（その１）では無平方数（Squarefree number）の数がメビウス関数を使うと簡単に求まるということを紹介しました。今回はその応用編です。
平方因子を１個だけ持つ数を数える
例としては$2^2 \times 3$は平方因子が１個で$2^2 \times 3^2 \times 5$は２個になります。
【例題】N以下の平方因子を１個だけ持つ数をsp(1,N)と表すときsp(1,$10^8$)を求めよ
まず前回のように平方因子が$4,9,25$の3つの場合で考えてみます、平方因子が$4$のみの数の個数をSS(4)とすると以下のように表されます。（ベン図を書いてみると分かりやすいかもしれません）
\begin{align}
&ss(4) = m(4)- (m(4 \times 9)+m(25 \times 4))\\
&+ m(4 \times 9 \times 25)\\
\end{align}
同様にss(9)とss(25)を求めて足し合わせれば良いので、
\begin{align}
&sp(1,40) = ss(4)+ss(9)+ss(25) \\
&= m(4)+m(9)+m(25)\\
&- \color{red} {2} \times (m(4 \times 9)+m(9 \times 25)+m(25 \times 4))\\
&+ \color{red} {3} \times m(4 \times 9 \times 25)\\
\end{align}
この式は前回と似ていますが赤で表される係数がつきました。この係数はそれぞれ２個と３個から１個に与える影響ということで、$\large \binom{2}{1}, \binom{3}{1}$で求められます。
これを前回と同様に変形します。ここではnの素因数の数をnf(n)で表します。
\begin{align}
&sp(1,40) = -\binom{1}{1}(\mu(2) m2(2) + \mu(3) m2(3) + \mu(5) m2(5))\\
&+ \binom{2}{1}((\mu(6)m2(6))\\
&= \sum_{k=1}^{\lfloor \sqrt{40} \rfloor} (-1)\mu(k) m2(k) \binom{\mathrm{nf(k)}}{1}\\
&= \sum_{k=1}^{\lfloor \sqrt{40} \rfloor} (-1)\mu(k) \lfloor 40/k^2 \rfloor \binom{\mathrm{nf(k)}}{1}\\
\end{align}
これを一般化すると以下の式が得られます
\begin{align}
& \mathrm{sp(1,N)} = \sum_{k=1}^{\lfloor \sqrt{n} \rfloor} (-1)\mu(k) \lfloor N/k^2 \rfloor \binom{\mathrm{nf(k)}}{1} \cdots (1)\\
\end{align}
さらに平方因子を増やして考えます
平方因子をr個だけ持つ数を数える
【例題】N以下の平方因子をr個だけ持つ数をsp(r,N)と表すとき$\mathrm{sp(r,10^$)}$を求めよ
式(1)の符号はrが奇数と偶数で変わるのと、係数が、$\large \binom{\mathrm{nf(k)}}{r}$になるだけなので以下の式になります。
\begin{align}
& \mathrm{sp(r,N)} = \sum_{k=1}^{\lfloor \sqrt{n} \rfloor} (-1)^r\mu(k) \lfloor N/k^2 \rfloor \binom{\mathrm{nf(k)}}{r} \cdots (2)\\
\end{align}
これを素直にpythonのコードにして$\mathrm{sp(r,10^$) \ (r=0, \cdots ,6)} $を求めました。kの素因数の数$\mathrm{nf(k)}$はsympyのprimefactorsを使っています。
import sympy as sp
def SP(r, N):
ret = 0
for k in range(1, int(N**(0.5))+1):
ret += sp.mobius(k)*(N//k**2)*sp.binomial(len(sp.primefactors(k)),r)
return (-1)**r*ret
N = 10**8
print(f"N = {N}")
for r in range(0,6+1):
print(f"SP({r},N) = {SP(r, N)}")
# N = 100000000
# sp(0,N) = 60792694
# sp(1,N) = 33539196
# sp(2,N) = 5329747
# sp(3,N) = 329028
# sp(4,N) = 9257
# sp(5,N) = 78
# sp(6,N) = 0
今回は平方因子を$r$個持つ数を正確にかぞえる方法を考えましたが、次回（その３）では$N \rightarrow \infty$としたときに平方因子を$r$個持つ数の割合を、確率の考え方を用いて求めます。
（開発環境：Google Colab）
この考え方はProject Euler、Problem 632: Square Prime Factorsを解くのに役に立ちます。