【科学技術計算講座1-6】Pythonでグラフを表示する

前回はPythonで一分子反応をオイラー法で計算するプログラムを作りました。

Pythonでオイラー法をプログラミングします。for文による繰り返し、数値誤差、出力フォーマットなどPythonプログラミングの基本を解説します。科学技術計算講座1「一分子反応をオイラー法でシミュレーション」の第5回目です。

AとBの濃度を数値データとして出力するところまで出来ました。本日は、そのデータをグラフにしてみましょう。

Pythonでグラフを表示する

Pythonでグラフを表示するには、Matplolibというライブラリを使います。ライブラリとはPythonで様々なことをできるようにするツールのことで、MatplotlibはPythonで色々なグラフを描くことができるようにしてくれます。

早速ですが、以下のプログラムをテキストエディタで書いて、Pythonで実行してみてください。実行方法は、今までと同じように、Anaconda Promptを使ってください。

import matplotlib.pyplot as plt

x = [0, 1, 2, 3, 4, 5]
y = [10, 12, 15, 20, 17, 16]

plt.plot(x, y)
plt.show()

すると、次のようなグラフが出力されます。

Pythonのグラフ

プログラムの1行目は、Matplotlibを使うよという宣言です。これは、Matplotlibをpltという名前(オブジェクトといいます)で読み込んでいます。

次の x と y のところですが、これは数値データを与えているところです。

x = [0, 1, 2, 3, 4, 5]

まず、x ですが、前回とは違って [ ] の中に複数の数値が書かれています。これは配列(Pythonではリストといいます)といって、箱のようなものです。箱のなかに、順番に数値が並んで入っていると思ってください。この箱の名前を x としています。 y についても同様に、6つの数値が入っています。

それぞれ6つデータがあるので、x と y をセットで考えると、(x, y) = (0, 10), (1, 12), (2, 15), (3, 20), (4, 17), (5, 16) の6点が定義されています。

次に、plt.plot(x, y) ですが、これはグラフの横軸と縦軸に相当するデータを指定します。この場合、x が横軸、yが縦軸のデータです。x、y は上で定義している配列のデータです。

最後の plt.show() でグラフを表示しています。このように、データを与えてやれば、Pythonで簡単にグラフを描くことができるのです。

計算しながらグラフにする

先程の例は与えられたデータをグラフにしていましたが、計算しながらデータを定義することもできます。

例えば、$y=(x-50)^2$ のグラフを書いてみましょう。プログラムは次のようになります。

import matplotlib.pyplot as plt

x = []
y = []

for i in range(100):
    x.append(i)
    y.append((i - 50)**2)

plt.plot(x, y)
plt.show()

今度は x = [] と、[ ] の中に数値が入っていません。これは、空の配列を作るという意味です。とりあえず、中身の無い空の箱だけ作っています。

そして、前回出てきたfor文で、100回ループさせています。繰り返しているのは、まず x.append(i) です。これは、先程の x という箱に、i を入れていくという命令です。i は0から 99 まで1ずつ増えていきますので、x には、0から99まで順番に数値が入っていきます。

次の y.append((i - 50)**2) ですが、これは y の箱に、(i - 50)**2 を入れていくという命令になります。**2 は2乗を表します。したがって、(i - 50) の2乗を y の配列に入れていきますが、i は0から99まで増えるので、それぞれの i で計算された値が追加されていきます。

グラフを描画する部分は同じです。このプログラムを実行してみてください。次のような $y=(x-50)^2$ グラフが出力されるはずです。

Pythonのグラフ2

一分子反応のシミュレーションをグラフにする

では、前回作った一分子反応の計算プログラムをグラフが描けるように変更してみましょう。AとBの濃度を時系列にプロットしたいので、グラフの横軸が時間、縦軸がAとBの濃度になるようにします。

上の例と同じように、まずは空の配列を作りましょう。配列の名前は、時間が gt、Aの濃度が ga、Bの濃度が gb としておきます。前回最後に作ったプログラムを開いて、for文の直前に以下を記述します。

gt = []
ga = []
gb = []

次に計算中に時間 t、A濃度 a、B濃度 b をこの配列に追加していきます。forループの中のprint文の下に、次の3行を書きましょう。字下げを合わせておいてくださいね。

    gt.append(t)
    ga.append(a)
    gb.append(b)

最後に、このデータを描画する命令を書きます。この時、字下げをしないで、元に戻しておいてください。字下げするとforループの中にあると見なされて、繰り返し描画されてしまいます。

plt.plot(gt, ga)
plt.plot(gt, gb)
plt.xlabel('t')
plt.ylabel('A,B')
plt.legend(['A', 'B'])
plt.show()

最初の plt.plot の2行は、横軸に時間 gt、縦軸にA、Bそれぞれの濃度(ga と gb)の配列データを指定しています。

次の plt.xlabel('t') ですが、これはX軸のラベルを「t」と指定しています。同様にY軸は、「A,B」のラベルをつけます。

次の plt.legend(['A', 'B']) は、凡例のラベルを「A」と「B」に指定して表示する命令です。今回のように複数のグラフを表示する場合は、凡例をつけるとわかりやすいです。凡例は配列と同じで、[ ] の中に書いてください。最後は、plt.show() でグラフを表示します。

では、このプログラムを保存して、Pythonで実行してみてください。次のようなグラフが出力されましたか?

プログラムを表示
import matplotlib.pyplot as plt

# parameter
a = 1.0
b = 0.0
t = 0.0
k = 0.1
dt = 0.1

# graph data array
gt = []
ga = []
gb = []

# time loop
for i in range(1000):
    fa = -k * a
    fb = k * a
    
    a = a + dt * fa
    b = b + dt * fb
    t = t + dt

    print('i: {0:4d}, t: {1:6.2f}, A: {2:9.6f}, B: {3:9.6f}'.format(i, t, a, b))
    
    gt.append(t)
    ga.append(a)
    gb.append(b)

# graph plot
plt.plot(gt, ga)
plt.plot(gt, gb)
plt.xlabel('t')
plt.ylabel('A,B')
plt.legend(['A', 'B'])
plt.show()
一分子反応のシミュレーション結果

グラフを見れば一目瞭然です。Aは減少して60秒くらいまでには、ほぼゼロになっています。BはAが減ると同時に増加して行きます。

まとめ

今回は、Pythonでグラフを表示する所をプログラミングしました。数値データだけでなく、グラフと合わせて見ることで、解きたい現象がどうなっているのかよくわかるようになります。Pythonでは、このようにグラフを簡単に描くことができるので、ぜひ活用してみてください。

さてグラフは出来て、数値データも画面に出力されているのですが、レポートなどにまとめる時に、数値データをファイルに出力しておいた方が便利という場合があります。そこで、次回はデータをファイルに出力する方法について学んでいきたいと思います。


←前回 Pythonでオイラー法をプログラミングする

→次回 Pythonでデータをファイルに出力する

全体の目次

スポンサーリンク
科学技術計算のご相談は「キャットテックラボ」へ

科学技術計算やCAEに関するご相談、計算用プログラムの開発などお困りのことは「株式会社キャットテックラボ」へお問い合わせください。

お問い合わせはこちら

フォローする