階乗と二重階乗
階乗とは次式で定義される演算です。
\[n!=n\,(n-1)\,(n-2) \cdots 2 \cdot 1, \quad 0!=1\]
たとえば、$5$ の階乗は
\[5!=5\cdot 4\cdot 3\cdot 2\cdot 1=120\]
のように計算します。階乗計算は単純な再帰アルゴリズムで実装できます:
\[f(n)=nf(n-1)\]
たとえば $4!$ を計算するときには、$f=4$ からスタートして、$f$ と自然数 $3,\ 2,\ 1$ の積を順次 $f$ 自身に代入していきます。
\[\begin{align*}&f=4\\[6pt]&f=f\times 3\\[6pt]&f=f\times 2\\[6pt]&f=f\times 1\end{align*}\]
Python で階乗関数を実装してみましょう。
# PYTHON_FACTORIAL # In[1] # 階乗関数 def fact(x): if x == 0: return 1 else: return x * fact(x - 1) print(fact(30)) # 265252859812191058636308480000000
math モジュールには階乗を計算する factorial() が用意されています。
# In[2] import math # 30の階乗を計算 val = math.factorial(30) print(val) # 265252859812191058636308480000000
二重階乗は
\[\begin{align*}n!!&=n\,(n-2)\,(n-4) \cdots 3 \cdot 1\quad (n=\mathrm{odd\:number})\\[6pt]n!!&=n\,(n-2)\,(n-4) \cdots 4 \cdot 2\quad (n=\mathrm{even\:number})\\[6pt]0!!&=(-1)!!=1\\[6pt]\end{align*}\]
によって定義されます。たとえば、$5$ の二重階乗を計算すると
\[5!!=5\cdot 3\cdot 1=15\]
となります。
math.factorial()
math.factorial(x) は x の階乗を返します。
引数 x に非整数や負数を渡すと ValueError を返します。
# MATH_FACTORIAL # In[1] # mathモジュールをインポート import math # 10の階乗を計算 val = math.factorial(10) print(val) # 3628800
scipy.special.factorial()
scipy.special.factorial() は整数 n を受け取ってn の階乗を返します。デフォルト (exact = False) では処理速度を優先して近似計算を実行しますが、exact に True を渡すと正確な値を返します。
# SCIPY_FACTORIAL # In[1] from scipy import special n = 50 # n!を近似計算 val = special.factorial(n) print(val) # 3.041409320171338e+64 # In[2] # n!の正確な値を計算 val = special.factorial(n, True) print(val) # 30414093201713378043612608166064768844377641568960512000000000000
mpmath.factorial()
mpmath.factorial() あるいは mpmath.fac() を使うと、任意精度で階乗計算を実行できます。
# MPMATH_FACTORIAL # In[1] from mpmath import * # 精度を10桁に設定 mp.dps = 10 # 10の階乗 x = factorial(10) print(x) # 3628800.0
引数に巨大数を受け取った場合は、スターリングの近似公式を使って階乗の近似値を計算します。
# In[2] # スターリングの公式を使って巨大数の階乗を計算 x = factorial(10**15) print(x) # 1.178796412e+14565705518096756
mpmath.fac2()
mpmath.fac2(x) は $x$ の二重階乗 $x!!$ を返します。
# PYTHON_FACTORIAL_2 # In[1] from mpmath import * # 精度を10桁に設定 mp.dps = 10 # 15!!を計算 val = fac2(15) print(val) # 2027025.0
sympy.factorial()
sympy.factorial(x) は x の階乗を返します。
# SYMPY_FACTORIAL # In[1] import sympy # 10!を計算 val = sympy.factorial(10) print(val) # 362880
SymPy は記号計算に対応するパッケージなので、引数に文字式を渡すと factorial オブジェクトが生成されます。
# In[2] # 記号nを定義 sympy.var('n') # (n+2)!を計算 s = sympy.factorial(n + 2) print(s) # factorial(n + 2)
$(n+2)!$ を $n!$ で割ると $(n+2)(n+1)$ になることを確認してみましょう。
# In[3] # n!を計算 t = sympy.factorial(n) # (n+2)!をn!で割って数式を簡略化する expr = sympy.simplify(s/t) print(expr) # (n + 1)*(n + 2)
sympy.factorial2()
sympy.factorial2(x) は x の 二重階乗 を返します。
# SYMPY_FACTORIAL_2 # In[1] import sympy # 10の二重階乗を計算 val = sympy.factorial2(10) print(val) # 3840
コメント
【GPT対話】「Pythonで階乗を計算する方法」
生徒:先生、Pythonで階乗を計算する方法を教えていただけますか?
先生:もちろんです。階乗を計算するには、いくつの階乗を求めたいかを指定し、それを計算するための関数を作成します。Pythonでは、再帰的なアプローチやループを使って階乗を計算することができます。どちらの方法を使って説明しましょうか?
生徒:両方の方法について教えていただけますか?
先生:もちろんです。まず、再帰的なアプローチから始めましょう。再帰関数を使用すると、関数自体の中で自身を呼び出すことができます。階乗の再帰的な計算は以下のようになります。
この関数では、nが0になったときに再帰が終了します。それ以外の場合、nをn-1で再帰的に掛け続けることで階乗を計算します。
生徒: 分かりました。では、具体的な例を見せていただけますか?
先生: もちろんです。例えば、5の階乗を計算する場合は次のようにします。
結果として、120が表示されます。
生徒: 理解しました。次にループを使った方法について教えていただけますか?
先生: はい、次はループを使った方法です。ループを使って階乗を計算するには、1からnまでの数を順番に掛け合わせていく方法です。
の関数では、1からnまでの数を順番に掛け合わせていくためにループを使用しています。
生徒: なるほど、再帰的なアプローチとループを使ったアプローチがありますね。どちらを使えばいいですか?
先生: どちらの方法を使っても階乗を計算することができますが、一般的には小さい数の階乗ではどちらの方法も効率的です。ただし、非常に大きな数の階乗を計算する場合は、再帰的なアプローチはスタックの深さによる制約がありますので、ループを使ったアプローチが好まれます。
【AI連載小説】科学とコードの交差点(30)「階乗計算の高速アルゴリズム」
Pythonサークルの中で、明信、瞳、開誠の3人が最速の階乗計算アルゴリズムについて熱心に討論していました。部屋には数式やコードのノートが広がり、情熱的な議論が交わされています。
瞳: さて、最速の階乗計算アルゴリズムについて考えてみましょう。それぞれどんな方法が思い浮かびますか?
開誠: まずはメモ化再帰が挙げられますね。計算結果をメモしておいて、再帰で同じ計算が出てきたらその結果を再利用する方法です。
明信: そうだな、でも再帰だと大きな数だとスタックオーバーフローが心配だし、関数呼び出しにもオーバーヘッドがある。
瞳: それに加えて、数学的な手法も視野に入れましょう。例えば、二項係数を用いて高速に計算する方法があります。
開誠: それは確かに効率的だけど、実装が複雑になりがちだよね。
明信: そうだな。なるべくシンプルで高速な方法が良いな。
瞳: Pythonのmathモジュールにあるfactorial関数も考慮すべきかもしれませんが、これも大きな数になると限界がありますね。
開誠: そうだね。でも、この辺りの関数もCで実装されていることが多いから、そのアルゴリズムを調査してみるのもいいかもしれない。
明信: なるほど。Cの実装を調べつつ、可能であればそれをPythonに組み込むというのは一つの手かもしれない。
瞳: それにしても、最速を求めるって難しいよね。計算量やオーバーヘッド、メモリ使用量など、いろいろな要素を考慮しないと。
開誠: そうだね。でも、ちょっとずつ試行錯誤していくしかないかもしれない。
3人はmathモジュールのfactorial関数が内部でどのようなアルゴリズムを使用しているかについて調査を始めました。彼らはPythonのソースコードやドキュメントを参照しながら、factorial関数がどのように実装されているのかを把握しようとしていました。
開誠: まずはmathモジュールのソースコードを見てみよう。どのファイルに入っているかな。
明信: よし、mathモジュールはCで書かれているから、Pythonのソースではなく、Cの方を見る必要があるな。
瞳: そうだね。PythonのmathモジュールのC実装は、CPythonと呼ばれるCで書かれたPythonの実装に含まれているはずだ。
三人はオンラインでCPythonのソースコードを検索し、mathモジュールのC実装を見つけました。
開誠: ここにmathモジュールのC実装があるようだ。では、factorial関数を探してみよう。
三人はソースコードをスクロールし、factorial関数の実装箇所を見つけました。
瞳: ここだね。factorial関数がどのように実装されているか見てみよう。
開誠: これは再帰を使った実装みたいだね。でも、再帰だと大きな数になるとスタックオーバーフローの可能性がある。
明信: そうだな。でも、この関数はPyLong_AsUnsignedLongLongを使って、Cのunsigned long long型に変換しているみたいだ。これならかなり大きな数も扱えそうだ。
瞳: なるほど、Cの型に変換しているから効率的に計算できるんだ。でも、再帰は避けた方がいいのかもしれないね。
開誠: そうだね。もっと効率的な実装方法はないか、他の関数も調べてみよう。