数値計算の基礎
本稿ではPythonで行う数値計算を具体的に書こうと思います。Part1では数値計算の思考プロセスについて書いたので、それを実現するためのプログラムを示します。まず数値計算の手法で解くことのできる数式ですが、今回はcos関数を可視化します(Part1の1.に対応)。
y = \cos{\theta}
cos関数は数値計算の手法で解くことのできる数式であり、複雑な式の変形は必要ないです。式の変形が必要となるのはもう少し難しい方程式に対してですので、その辺に関しては後にまわします。まず数値計算に慣れるところから始めようと思います。
数値計算の手法で解くことのできる数式を用意したら、定義域の設定をしますが(Part1の2.に対応)、今回は0[rad]から20[rad]までの区間でcos関数を書き出すことにします。Pythonのmathモジュールに含まれるcos(x)関数はx[rad]の余弦を計算します。単位がラジアンであることに注意してください。次に、cos関数を計算するためにこれを細かく分割します(Part1の3.に対応)。この分割が細ければ細かいほど多くの点をプロットすることになり、可視化されたときに点の集合がなめらかな曲線を引きますが、計算量が多くなることで時間がかかることに注意してください。今回は0.1[rad]ずつ角度を増やして、その角度におけるcos関数を出力しますと。Part1で紹介した積分の計算では、積分区間を細かく分割して、分割した区間とその値域を平均したものをかけて求めた台形の面積を全て足しあわせることで積分を求めるという説明しました。これは数値計算の手法として区分求積法を使うことに等価です。cos関数も同様に区間の分割を行います。分割した各点に対応するcos関数の値を求め、その各点をグラフにプロットすることで可視化が実現します。可視化にはmatplotlibというライブラリーのpylabというモジュールを利用します(Part1の5.に対応)。今回はPartの4.にあるようなデータの保存はしません。データを保存する必要があるような問題が出てこない限り、データの保存はフォルダを汚すだけです。
長々と説明しましたが、いよいよプログラムをみていくことにしましょう。
import math, pylab #①
steps = 200 #②
dtheta = 0.1
theta = [0.0]*steps #③
cos = [0.0]*steps
for i in range(steps): #④
theta[i] = i*dtheta
cos[i] = math.cos(theta[i])
print(theta) #⑤
print(cos)
pylab.plot(theta, cos) #⑥
pylab.show()
このプログラムではまず①でmathモジュールとpylabモジュールをインポートします。②で定義域と分割数を設定しています。0.1[rad]ずつ角度を増やして、最終的には200回繰り返すことで定義域の最大が20[rad]になるようにしています。③では定義域としての角度と値域としてのcos関数についてリストを繰り返す回数分だけ生成します。④では角度リストのi番目の値を求めて、それをmathモジュールのcos関数に渡し、結果をcos関数リストのi番目に格納します。このfor文が終了した時点で、角度リストとcos関数リストに値が格納されていることを確認するために⑤のprint文を書いてみました。これは数値計算には本質的に関係してこないのですが、僕は簡単な計算が正しく行われているのかを確認するためによく使っています。ことごとくprint関数でデータを出力することでどこで間違えたのか、あるいは正しいのかを簡単に知ることができます。また、print関数ではなくwrite関数を使うとデータの書き出しができたりもします。⑥でデータの可視化を行います。pylabモジュールのplot関数に角度リストとcos関数リストを渡すことで二次元のグラフが出来上がります。
最後に
本稿では数値計算の初歩をみました。こんな感じで数値計算の入門的な内容を扱っていきます。Pythonに慣れている人は考え方さえわかってしまえばプログラムの部分はそこまで難しくなかったと思います。慣れていない人はコピー&ペーストではなく自分の手で打ち込んでみてください。自分の手で打ち込むことで高い学習効果が得られると誰かが言っていたような気がしたりしなかったりします。プログラミングで数学を勉強したいという高校生や大学生にとってもPythonは学習が容易だと思っているので、難しく考えずにやってみて欲しいです。cos関数は高校生でも習いますよね。中学生は放物線を描いてみるといいと思います。
余談ですが、僕は開発環境としてPyCharmを使っています。