無平方数（Squarefree number）と メビウス関数（その１）

### メビウス関数

[メビウス関数(Wikipedia)](https://ja.wikipedia.org/wiki/%E3%83%A1%E3%83%93%E3%82%A6%E3%82%B9%E9%96%A2%E6%95%B0)は以下のように定義されます。

```math 
\begin{equation}
  \mu(n)=
  \begin{cases}
    0       & \text{$nは平方因子を持つとき$,} \\
    1       & \text{$nは偶数個の相異なる素数の積のとき$,} \\
    -1      & \text{$nは奇数個の相異なる素数の積のとき$。}
  \end{cases}
\end{equation}
```
$n \le 12$では

|n | 1| 2| 3| 4| 5| 6| 7| 8| 9| 10|11 |12 | 
|:---:|:---:|:---:|:---:|:---:|:---:|:---:|:---:|:---:|:---:|:---:|:---:|:---:|
|$\mu(n)$|1|-1|-1|0|-1|1|-1|0|0|1|-1|0|



名前は聞いたことがあるけど何に使うの？と思っている人は多いのではないでしょうか？　今回はなるほどと思えるような問題を解いていきたいと思います。

### 無平方数（Squarefree number）を数える

[無平方数（Squarefree number）(Wikipedia)](https://ja.wikipedia.org/wiki/%E5%B9%B3%E6%96%B9%E5%9B%A0%E5%AD%90%E3%82%92%E3%82%82%E3%81%9F%E3%81%AA%E3%81%84%E6%95%B4%E6%95%B0)とは、$2^2=4$や$3^2=9$などの平方数で割ることのできない正の整数のことで、例えば$2 \cdot 3 \cdot 5 = 30$が無平方数となります。

**【例題】N以下の無平方数の数をsf(N)と表すときsf($10^8$)を求めよ**
 
もちろんメビウス関数を使って$\mu(n) \neq 0$のnをすべて数えても良いのですが、もう少し賢い方法があります。それをまずsf(40)の場合で考えてみます

1. まず平方因子を持つ数を数えて、それをnから引く
1.  $n\le 40$以下の平方因子を持つ数は$2^2=4、3^2=9　または  5^2=25$の倍数
1.  この**4か9か25**の倍数の数はまさに「[包除原理（Inclusion-exclusion principle）でProject Euler Problem 1を解く（その2）](https://qiita.com/masa0599/items/32923278dc6087b3473d)」と同様に考えることができるので、30以下のkの倍数の数を$m(k)$と表すと

```math 
\begin{align}
&sf(40) = 40- (m(4)+m(9)+m(25))\\
&+ (m(4 \times 9)+m(9 \times 25)+m(25 \times 4))\\
&- m(4 \times 9 \times 25)\\
\end{align}
```
と表されます。要は掛け合わせる個数が奇数なら$-$, 偶数なら$+$となるということですね。40はm(1)と表し、40を超えるものを消し、さらに$m(k^2)=m2(k)$とすると、

```math 
\begin{align}
&sf(40) = m(1)- (m(2^2)+m(3^2)+m(5^2))\\
&+ (m((2 \times 3)^2))\\
&= m2(1) - m2(2) - m2(3) - m2(5) + m2(6) \cdots (1)\\
\end{align}
```

この$(1)$と上のメビウス関数の値を見比べると**係数が一致している**のが分かります。すなわち無平方数だけを抜き出しその素因数の個数で$\pm$が決まるわけですね。これを$\mu(n)$を使って書き、$m2(k) = \lfloor n/k^2 \rfloor$ なので

```math 
\begin{align}
&sf(40) = \mu(1)m(1) + \mu(2) m2(2) + \mu(3)  m2(3) + \mu(5)  m2(5) + \mu(6)  m2(6)\\
&= \sum_{k=1}^{\lfloor \sqrt{40} \rfloor} \mu(k) m2(k) \\
&= \sum_{k=1}^{\lfloor \sqrt{40} \rfloor} \mu(k)  \lfloor 40/k^2 \rfloor \\
\end{align}
```
これを一般化すると以下の式が得られます。

```math 
\begin{align}
&\large \mathrm{sf(n)} = \sum_{k=1}^{\lfloor \sqrt{n} \rfloor} \mu(k)  \lfloor n/k^2 \rfloor \\
\end{align}
```
これをpythonのコードにするとsympyのmobius関数を使えば1行で書くことができます。

```python 
import sympy as sp
N = 10**8
print(f"SF({N})={sum(sp.mobius(k)*(N//k**2) for k in range(1,int(N**(0.5))+1))}")
# SF(100000000)=60792694
```

今回は平方因子ない0の無平方数の数を考えましたが[次回（その２）](https://qiita.com/masa0599/items/ebeaa36d908299146c7f)は、平方因子が$r$個ある場合の数を考えます。

### 【別解】メビウス関数を使わない解法

メビウス関数を使わない解法も見つけましたので紹介します。平方因子を含む数の個数を求めてnから引くという考え方です

```math 
\begin{align}
&sq(n,i) : n以下のiが最大平方因子の数の個数\\
&例えばsq(16,2)=|\{4,8,12\}|=3, (16は最大平方因子が4)\\ \\
&この最大平方因子iを使って平方因子を含む数は\\
&i^2j　(jは無平方数) と表せる\\
&するとn以下の最大平方因子がiの個数は\\
&jの個数すなわち\lfloor n/i^2 \rfloor 以下の無平方数の個数になる\\
&式にすると\\
&sq(n,i) = sf(\lfloor n/i^2 \rfloor)\\
&sf(n)はすべての最大平方因子iで合計してnから引けばよいので\\
&sf(n) = n - \sum_{i=2}^{\sqrt{n}}sq(n,i) = n - \sum_{i=2}^{\sqrt{n}}sf(\lfloor n/i^2 \rfloor)\\
&という漸化式が得られます\\
\end{align}
```
これは素直にpythonのコードにすると以下のようなシンプルなものになります。

```python 
from functools import lru_cache
@lru_cache(maxsize = None)
def SqFree(n):
  return n - sum(SqFree(n//i**2) for i in range(2,int(n**(0.5))+1))
print(SqFree(10**8))
# 60792694
```


（開発環境：Google Colab）

この考え方はProject Euler、Problem 193: Squarefree Numbersを解くのに役に立ちます。