Views
Document
Event
File
Folder
Image
Link
News Item
Photo
Photo Album
Topic
Weblog
Wiki Page
数値計算
Document Actions
~/python/14 で、やるのよ！
はじめに
- 計算機の使い道として、数値的な計算を行うことがあります。
- 計算の方法を「アルゴリズム」と呼び、
- その計算方法を具体化したものをプログラムといいます。
- 最終回は、数値計算のいくつかの例を取り上げます。
階乗の計算
- 階乗 n! の計算には、
-
n! = 1 * 2 * ... * (n-1) * n
-
n! = 1 (n = 0), n * (n-1)!
-
- の２つの方法が考えられる。
階乗の計算プログラム
- ２つの方法をそれぞれ、fact1(), fact2() の関数にしてみます。
# -*- coding: utf-8 -*- # a.py # 階乗の計算 # def fact1(n): f = 1 for i in range(1,n+1): f = f * i return f def fact2(n): if n == 0: return 1 return n * fact2(n-1)
実行例：(２つの方法を比較)
>>> from a import * >>> fact1(10) 3628800 >>> fact2(10) 3628800
- あたりまえですが、同じ結果になります。
説明
- ２番目のように，fact2 の定義の中に fact2 を使用しているようなアルゴリズムを（つまり，関数の定義が終わってないのにそれをその関数の中で使ってしまう）「再帰的」アルゴリズムと呼びます。
- 再帰的アルゴリズムはプログラムの構文解析やゲームの手の探索などにも用いられることがあり，「プログラムの記述を簡単にして，複雑なことは計算機に任せる」手法と考えることができます。
Newton 法による平方根の計算
平方根の計算
- ここでは、
- xk+1 = xk - f(xk)/f'(xk)
- において
- f(x) = x2 - a, f'(x) = 2*x
- であることから、次のようなプログラムを書くことができます。
Newton法のプログラム
# -*- coding: utf-8 -*- # b.py # Newton 法による平方根の計算 # def Sqrt(x, e=0.0001): p = 0.0 q = x while abs(q - p) > e: print q p = q q = p - (p*p - x)/(2.*p) return q
補足
- 但し、このプログラムでは上記の式における a を入力 x とし、
- xk を p, xk+1 を q とし、
- p と q の差が e 以下になったら、計算を終了するようにしています。
- また、q の初期値を x としています。
実行例：(math モジュールによる計算結果と比較)
>>> from b import * >>> Sqrt(2.0) 2.0 1.5 1.41666666667 1.41421568627 1.4142135623746899 >>> import math >>> math.sqrt(2.0) 1.4142135623730951
補足
- e の値はオプションなので、省略すると関数定義に書かれた 0.0001 が使われます。この値を調節してみてください (Sqrt(2.0, 0.01)のように呼び出す)
- この例のように，いい加減な値から良い解に近づいていくような計算方法を一般に「漸近法（ぜんきんほう）」と呼びます。（Newton法は漸近法の１種です）
問題1
- 三乗根を計算する関数 Cbrt(x[, e]) を q1.py に作成し、x ** (1./3) の結果と比較しなさい。但し、収束判断値 e はオプションである。
問題2
- 2つの正の整数の最大公約数を求めよう。ユークリッドの互除法
は次のように再帰的に表現することができる
- もし m % n (m を n で割った剰余) が 0 であれば n が最大公約数である
- でなければ m, n の最大公約数は n, (m % n) の最大公約数として求める
- これを q2.py に gcd という名前の関数で作成しなさい
- 実行例
>>> from q2 import gcd >>> gcd(1071, 1029) 21
解答は、"メール"：mailto:tkikuchi＋ci2006＠is.kochi-u.ac.jp?subject=14_Numeric で送ってください。
また，今回の演習の結果できたファイルはそのまま ~/python/14 に残しておいてください。評価の対象とします。
