Rubyで素数で遊ぶ(prime モジュール)

Ruby の [prime モジュール](https://docs.ruby-lang.org/ja/latest/library/prime.html)で遊びます。

## はじめに

`require 'prime'` を実行します。


```rb
[1] pry(main)> require 'prime'
=> true
[2] pry(main)> RUBY_VERSION
=> "2.5.1"
```

## N 以下の素数を求める

```rb
[5] pry(main)> Prime.each(20)
=> #<Prime::EratosthenesGenerator:0x00007ff862c50708 @last_prime_index=-1, @ubound=20>
[6] pry(main)> Prime.each(20).to_a
=> [2, 3, 5, 7, 11, 13, 17, 19]
```

## 小さな値から順に N 個の素数を求める

```rb
[7] pry(main)> Prime.take(10)
=> [2, 3, 5, 7, 11, 13, 17, 19, 23, 29]
```

## N 番目の素数を求める

```rb
[9] pry(main)> Prime.take(10).last
=> 29
```

`take` はリストを作ります。N 番目の要素のみが欲しいのでリストは不要ですね。

以下のように `next` を指定回数だけ呼び出すと、特定の要素のみ取得できます。

```rb
[11] pry(main)> p = Prime.each; (10-1).times { p.next }; p.next
=> 29
```

以下のように `each_with_index` と `find` を組み合わせて N 番目の要素を取得することもできます。

```rb
[12] pry(main)> Prime.each_with_index.find {|p, i| i == 10-1 }[0]
=> 29
```

## N を素因数分解する

[prime_division](http://docs.ruby-lang.org/ja/2.2.0/class/Prime.html#I_PRIME_DIVISION) メソッドを使います。

```rb
[13] pry(main)> Prime.prime_division(10) # 10 = 2 * 5
=> [[2, 1], [5, 1]]
[14] pry(main)> Prime.prime_division(20) # 20 = 2^2 * 5
=> [[2, 2], [5, 1]]
[15] pry(main)> Prime.prime_division(123456789)
=> [[3, 2], [3607, 1], [3803, 1]]
```

`prime_division` は、`[素数, 個数]` の形で素数がまとめられます。フラットにするには以下のようにします。

```rb
[16] pry(main)> Prime.prime_division(2600).map {|p, e| [p] * e }.flatten
=> [2, 2, 2, 5, 5, 13]
```

## N が素数かどうか調べる

[prime?](https://docs.ruby-lang.org/ja/latest/class/Prime.html#I_PRIME--3F) メソッドを使います。

```rb
[17] pry(main)> Prime.prime?(13)
=> true
[18] pry(main)> Prime.prime?(1024)
=> false
```

[Integer#prime?](https://docs.ruby-lang.org/ja/latest/method/Integer/i/prime=3f.html) を使うこともできます。

```rb
[3] pry(main)> 13.prime?
=> true
[4] pry(main)> 1024.prime?
=> false
```

## N 以上の最小の素数を求める

```rb
[19] pry(main)> Prime.find {|p| p >= 1_000_000 }
=> 1000003
```

## N 以下の最大の素数を求める

```rb
[22] pry(main)> Prime.take_while {|p| p < 1_000_000 }.last
=> 999983
```

## N 以上の素数を、小さな順から M 個求める

```rb
[23] pry(main)> Prime.lazy.drop_while {|p| p < 30 }.take(10).to_a
=> [31, 37, 41, 43, 47, 53, 59, 61, 67, 71]
```

`lazy` をはさまないと、`drop_while` が無限の長さの配列を作ろうとして結果が返ってこなくなります。

## 約数の個数を求める

$n$ を素因数分解して次のように表せるとき、

$$
n = p_1^{e_1} \cdot p_2^{e_2} \cdots p_m^{e_m}
$$

$n$ の約数の個数は、次式で求めることが出来ます。

$$
(e_1 + 1) \cdot (e_2 + 1) \cdots (e_m + 1)
$$

```rb
require 'prime'

def number_divisors(n)
  Prime.prime_division(n).map {|p, e| e + 1 }.inject(:*)
end
```

```rb
[3] pry(main)> number_divisors(2 ** 2 * 5)
=> 6
[4] pry(main)> number_divisors(1024)
=> 11
```

※`Prime.prime_division` は、[Integer#prime_division](https://docs.ruby-lang.org/ja/latest/method/Integer/i/prime_division.html)を代わりに使うこともできます。

## 約数の総和を求める

n を素因数分解して次のように表せるとき、

$$
n = p_1^{e_1} \cdot p_2^{e_2} \cdots p_m^{e_m}
$$

n の約数の総和は次式で表せます。

$$
(1 + p_1 + p_1^2 + \cdots + p_1^{e_1}) \cdot (1 + p_2 + p_2^2 + \cdots + p_2^{e_2}) \cdot \cdots \cdot (1 + p_m + p_m^2 + \cdots + p_m^{e_m})
$$

等比数列の和の公式より、上の式は

$$
\frac{p_1^{e_1 + 1} - 1}{p_1 - 1} \cdot \frac{p_2^{e_2 + 1} - 1}{p_2 - 1} \cdot \cdots \cdot \frac{p_m^{e_m + 1} - 1}{p_m - 1}
$$

となります。

```rb
require "prime"

def sum_divisors(n)
  Prime.prime_division(n).map {|p, e| (p ** (e + 1) - 1) / (p - 1) }.inject(:*)
end
```

```rb
[1] pry(main)> sum_divisors(12)
=> 28
[2] pry(main)> sum_divisors(1800)
=> 6045
```

## 約数を列挙する

```rb
require 'prime'

def _divisors(primes, k)
  Enumerator.new do |y|
    if primes.size == k
      y << 1
    else
      p, e = primes[k]
      _divisors(primes, k + 1).each {|d|
        (e + 1).times {|e1|
          y << p ** e1 * d
        }
      }
    end
  end
end

def divisors(n)
  _divisors(Prime.prime_division(n), 0).to_a.sort
end
```

```rb
[1] pry(main)> divisors(10)
=> [1, 2, 5, 10]
[2] pry(main)> divisors(1024)
=> [1, 2, 4, 8, 16, 32, 64, 128, 256, 512, 1024]
```

追記: コメント欄に @scivola さんによる[約数を列挙するコード](https://qiita.com/norioc/items/5c270ad6ff5d22d1efd0#comment-da13f6b749dba662ec63)があります。

## 参考

- [library prime](https://docs.ruby-lang.org/ja/latest/library/prime.html)
- [無限リストを map 可能にする Enumerable#lazy](http://magazine.rubyist.net/?0041-200Special-lazy)
- [約数](http://d.hatena.ne.jp/inamori/20091110/p1)
    - 約数の個数、総和、列挙の求め方を参考にしました。
- [約数の総和を求める二つの公式と証明](https://mathtrain.jp/yakusuwa)

## 更新履歴

- 2019/7/10 コメント欄で @scivola さんにご指摘いただいた点について修正。感謝です。