Python で素因数分解と約数を求める。


def usage():
    n = 2**2 * 3**2 * 5**2
    
    print('# 素因数')
    fct = factorize(n)
    print(fct, num(fct))
    
    print('# 約数')
    for div in divisorize(fct):
        print(div, num(div))


def factorize(n):
    fct = []  # prime factor
    b, e = 2, 0  # base, exponent
    while b * b <= n:
        while n % b == 0:
            n = n // b
            e = e + 1
        if e > 0:
            fct.append((b, e))
        b, e = b + 1, 0
    if n > 1:
        fct.append((n, 1))
    return fct


def divisorize(fct):
    b, e = fct.pop()  # base, exponent
    pre_div = divisorize(fct) if fct else [[]]
    suf_div = [[(b, k)] for k in range(e + 1)]
    return [pre + suf for pre in pre_div for suf in suf_div]


def num(fct):
    a = 1
    for base, exponent in fct:
        a = a * base**exponent
    return a


if __name__ == '__main__':
    usage()


実装

素因数分解

ひたすら割り続けています。実際にコードを見ると単純に割り続けているだけであることがわかるかと思います。

素数の求め方は、ここを参考にしています。
素数と仲良しになる3つのプログラム - ぴよぴよ.py
ttp://cocodrips.hateblo.jp/entry/2014/02/02/011119

◯ 約数

素因数分解の結果をもとに、組み合わせを、全て列挙しています。これはコードよりも先に、下の実行結果を見た方が、何をしているか伝わりやすいかもしれません。

下は sample.py の実行結果です。900 を素因数分解して、約数を全て求めています。この実行結果から、約数の実装が、素因数分解の結果をもとに、組み合わせを、全て列挙しているだけであることが、もし伝われば幸いです。

$ sample.py
# 素因数
[(2, 2), (3, 2), (5, 2)]

# 約数
[(2, 0), (3, 0), (5, 0)] 1
[(2, 0), (3, 0), (5, 1)] 5
[(2, 0), (3, 0), (5, 2)] 25
[(2, 0), (3, 1), (5, 0)] 3
[(2, 0), (3, 1), (5, 1)] 15
[(2, 0), (3, 1), (5, 2)] 75
[(2, 0), (3, 2), (5, 0)] 9
[(2, 0), (3, 2), (5, 1)] 45
[(2, 0), (3, 2), (5, 2)] 225
[(2, 1), (3, 0), (5, 0)] 2
[(2, 1), (3, 0), (5, 1)] 10
[(2, 1), (3, 0), (5, 2)] 50
[(2, 1), (3, 1), (5, 0)] 6
[(2, 1), (3, 1), (5, 1)] 30
[(2, 1), (3, 1), (5, 2)] 150
[(2, 1), (3, 2), (5, 0)] 18
[(2, 1), (3, 2), (5, 1)] 90
[(2, 1), (3, 2), (5, 2)] 450
[(2, 2), (3, 0), (5, 0)] 4
[(2, 2), (3, 0), (5, 1)] 20
[(2, 2), (3, 0), (5, 2)] 100
[(2, 2), (3, 1), (5, 0)] 12
[(2, 2), (3, 1), (5, 1)] 60
[(2, 2), (3, 1), (5, 2)] 300
[(2, 2), (3, 2), (5, 0)] 36
[(2, 2), (3, 2), (5, 1)] 180
[(2, 2), (3, 2), (5, 2)] 900
$


書き換え 1 - 文法

もし divisorize で2重リスト内包表記や条件式(三項演算子)を使わないならこんな感じになります。

def divisorize(fct):
    div = []
    if not fct:
        div.append([])
        return div
    b, e = fct.pop()  # base, exponent
    pre_div = divisorize(fct)
    suf_div = [[(b, k)] for k in range(e + 1)]
    for pre in pre_div:
        for suf in suf_div:
            div.append(pre + suf)
    return div
1. 二重のリスト内包表記

Python では可読性の低い書き方を嫌う文化にあるような気がします。Effective Python でも「項目8:リスト内包表記には、3つ以上の式を避ける」と述べられていました。でも、これくらい簡単なら使ってもいいかなと思って使いました。確かに読みづらいかも笑。

def divisorize(fct):
    ...
    return [pre + suf for pre in pre_div for suf in suf_div]
def divisorize(fct):
    div = []
    ...
    for pre in pre_div:
        for suf in suf_div:
            div.append(pre + suf)
    return div


二重のリスト内包表記は、二次元リストを1次元にするときにも使ったりします。

lst2d = [[i + 3*j for i in range(3)] for j in range(3)]
lst2d == [[0, 1, 2], [3, 4, 5], [6, 7, 8]]
lst1d = [elm for lst1d in lst2d for elm in lst1d]
lst1d == [0, 1, 2, 3, 4, 5, 6, 7, 8]
2. 条件式

よく三項演算子と言われるものです。これも、可読性を下げるものの一つであるような気がします。ただ、再帰呼び出しの k = 0 の時の動作を定義する時には、これは結構使えるかなと思ったりもします。k = 0 のときの動作という、特殊な条件で全体のコードを切ってしまうことがなくなります。

def divisorize(fct):
    b, e = fct.pop()  # base, exponent
    pre_div = divisorize(fct) if fct else [[]]
    ...
def divisorize(fct):
    # いきなり not fct のときなんていう
    # 特殊な条件を示されてもピンとこない。
    if not fct:
        return [[]]
    b, e = fct.pop()  # base, exponent
    pre_div = divisorize(fct)
    ...

条件式 (しばしば "三項演算子" と呼ばれます) は最も優先度が低いPython の演算です。
x if C else y という式は最初に条件 x ではなく C を評価します; C が true の場合 x が評価され値が返されます; それ以外の場合には y が評価され返されます。
6.12. 条件式 (Conditional Expressions)


例えば階乗や

factorial = lambda n: n * factorial(n - 1) if n else 1
def factorial(n):
    if n:
        return n * factorial(n -1)
    else:
        return 1


ユークリッド互除法なんかも、こんな風に書けます。

gcd = lambda a, b: gcd(b, (a % b)) if a % b else b
def gcd(a, b):
    if a % d:
        return gcd(b, (a % b))
    else:
        return b


ただ Guido は lambda 式は嫌いらしいです。

私は lambda が、いいと思ったことがない。

  • 不自由(たった1行しか書けない)
  • 紛らわしい(引数リストの括弧がない)
  • 普通の関数定義で代用できる。

I've never liked lambda

  • crippled (only one expression)
  • confusing (no argument list parentheses)
  • can use a local function instead

Python Regrets


PEP 8 でも lambda は、関数の定義では使用しないように指示されています。lambda を使用するのは、原則、関数オブジェクトを引数にとる関数に代入するときだけ。また変数などに代入せず、そのまま引数として与えてください。

# OK
list(map(lambda x: x**2, range(10)))

# NG
f = lambda x: x**2
list(map(f, range(10)))

lambda 式を識別子に直接束縛する代入文ではなく、常に def 文を使ってっください。
Always use a def statement instead of an assignment statement that binds a lambda expression directly to an identifier.

Yes:
def f(x): return 2*x

No:
f = lambda x: 2*x

最初の形式は、結果としてえられる関数オブジェクトの名前が、一般的な <lambda> ではなく 'f' という名前がつけられていることを意味しています。関数オブジェクトに文字列の名前が与えられていることは、一般に例外が発生した時にそれをトレースバックさせたり、関数名を文字列で出力させる際に役立ちます。代入文を使うことは(代入文で lambda 式を変数に束縛してから map, filiter などの高階関数に引数として与えることは)、 def 文にはなく lambda 式にある、たった1つの利点(すなわち、lambda 式は、より大きな式の中に埋め込められるということ)を無意味なものにしてしまいます。
The first form means that the name of the resulting function object is specifically 'f' instead of the generic '<lambda >'. This is more useful for tracebacks and string representations in general. The use of the assignment statement eliminates the sole benefit a lambda expression can offer over an explicit def statement (i.e. that it can be embedded inside a larger expression)
Programming Recommendations - PEP8

書き換え 2 - データ構造

素因数分解 factorize が返すデータ構造を [(2, 2), (3, 2), (5, 2)] ではなく [2, 2, 3, 3, 5, 5] にして実装して見ました。しかし約数を求める速度が 10 倍、遅くなってしまいました笑
divisor.py

約数を求めるときは、[(2, 2), (3, 2), (5, 2)] として取り扱った方が良さそうです。しかし、それ以外の場合は、[2, 2, 3, 3, 5, 5] として取り扱っても、いいかなという気もします。シンプルですしね。

def factorize_naive(n):
    b = 2
    fct = []
    while b * b <= n:
        while n % b == 0:
            n //= b
            fct.append(b)
        b = b + 1
    if n > 1:
        fct.append(n)
    return fct


def divisorize_naive(fct):
    # bit list
    bit_len = len(fct)
    r = range
    bit_list = [[(i >> j) % 2 for j in r(bit_len)] for i in r(2**bit_len)]
    # divisor
    div = [[f for f, b in zip(fct, bit) if b] for bit in bit_list]
    div[0] = [1]
    # exclude duplication, such as  [2, 2, 0], [2, 0, 2], [0, 2, 2].
    div = [list(t) for t in set([tuple(d) for d in div])]
    return div


def num_naive(fct):
    a = 1
    for b in fct:
        a = a * b
    return a

for 文 それとも while 文?

素因数分解を求めるときは while 文の方が for 文よりも速い気配があるので、while 文で書きました。while 文の方が綺麗に書けますしね。約数を求めるときも while でかけなくもないですが、あまりに汚くなってしまうので、for 文で書きました。速度の比較はしていません。

ちなみに、エラトステネスの篩をもちいて素数のリストを得ようとするときは for 文の方が速かったです。なぜかは、わかりません。while 文と for 文の比較は、下のページでやりました。