pythonで等高線(コンター図)を描画する

エクセルで書こうとすると手間がかかる等高線ですが、pythonであれば簡単に描くことができます。
pythonで等高線(コンター図)を描画するプログラムはいろいろな方が紹介されていますが、その中でも私が一番良いと思うものを紹介します。

参考にさせていただいたページ

https://hk29.hatenablog.jp/entry/2020/01/31/161508

等高線図の作図

万有引力のポテンシャルを例に等高線図を作成します。
等高線図と3D図を作成します。

どちらもmatplotlibに含まれている関数を使います。
等高線図はcontourを、3D図はplot_surfaceを使います。

先ほどのページで紹介されている自作の関数を定義して、描画します。

万有引力のポテンシャルとは?

万有引力の法則とは、「2つの物体の間には、物体の質量に比例し、2物体間の距離の2条に反比例する引力が作用する」という法則です。
万有引力の大きさFは、物体の質量をM, m、物体間の距離をrとして、

\displaystyle F = G\frac{Mm}{r^{2}}

とあらわされます。

万有引力のポテンシャルは、
\displaystyle U = G\frac{Mm}{r}

となります。

万有引力のポテンシャルを計算

1つの物体が作る万有引力のポテンシャルを描画するのでは面白くないので、複数の物体が作るポテンシャルを描画してみます。
とはいえ、簡単のために物体同士の相互作用は考えず、ポテンシャルの単純な和で表現します。
物体の座標はランダム関数を使って生成します。

ソースコード


import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits import mplot3d
import random

### 等高線図と3D図の作成例1:ジャグ配列のndarryをグラフセットする場合の関数
def my_contur_and_3Dgraph(X, Y, Z, function_name):
 print('X -> ' + str(X) + '\n' + str(type(X)))
 print('Y -> ' + str(Y) + '\n' + str(type(Y)))
 print('Z -> ' + str(Z) + '\n' + str(type(Z)))

 fig = plt.figure(figsize=(9, 4))

 # 等高線
 ax1 = fig.add_subplot(121, facecolor="w")
 contour = ax1.contour(X, Y, Z, levels=[-100, -10, -5, -2, -1, -0.5, -0.1], cmap="bwr_r")
 contour.clabel(fmt='%1.1f', fontsize=11)
 # オプション
 ax1.set_title(function_name)
 ax1.set_xlabel(my_xlabel, fontsize=12)
 ax1.set_ylabel(my_ylabel, fontsize=12)
 plt.grid()
 print(type(contour)) # <class 'matplotlib.contour.QuadContourSet'>

 # 3D図
 ax2 = fig.add_subplot(122, projection="3d", facecolor="w")
 surf = ax2.plot_surface(X, Y, Z, cmap="bwr_r", alpha=0.9, cstride=1, rstride=1, lw=0.1)
 # オプション
 #fig.colorbar(surf, shrink=0.6, aspect=10) # カラーバー追加。surfでcmapを指定必要
 ax2.set_title(function_name)
 ax2.set_xlabel(my_xlabel, fontsize=12)
 ax2.set_ylabel(my_ylabel, fontsize=12)
 ax2.set_zlim(-10, 0)

 # 作図の実行
 #fig.savefig(now + '_' + function_name + "_contour.png")
 plt.show()
 plt.close()

random.seed(1)

x = np.linspace(-10, 10, 1000) #x軸の描画範囲の生成。0から10まで0.05刻み。
y = np.linspace(-10, 10, 1000) #y軸の描画範囲の生成。0から10まで0.05刻み。

X, Y = np.meshgrid(x, y)

#点電荷の位置座標と電荷
xi = []
yi = []
p = 30

for n in range(p):
 xi.append(random.uniform(-10, 10))
 yi.append(random.uniform(-10, 10))

<br data-mce-bogus="1">

#ポテンシャル
Z = 0
for n in range(p):
 Z = Z - 1/np.sqrt[I](X-xi[n])**2 + (Y-yi[n])**2

my_xlabel = 'X'
my_ylabel = 'Y'

my_contur_and_3Dgraph(X, Y, Z, 'gravity well')


描画結果

余談

今回、万有引力のポテンシャルで等高線図を描いてみたのは、スターウォーズの影響です。
映画ハン・ソロで、ブラックホール群を通る時に、巨大生物に襲われ、わざとブラックホールに近づかせることで化け物を倒す、というシーンがあります。
これをみて、安全にブラックホールの間を抜けるルートというのは、ブラックホールの作るポテンシャルをもとに決まるのではないか、と考えたのがきっかけです。
そう思って上の図をみると、なんとなく航空ルートを示しているように見えてきませんか?

脚注

脚注
I (X-xi[n])**2 + (Y-yi[n])**2

sciencompass34 has written 145 articles

はじめまして!”あおやぎ”と言います。
メーカーで研究開発の仕事をしています。このブログでは、私の専門分野である半導体やそれに関連する内容を紹介していきます。
半導体関連の知識をまとめたデータベースのようにしたいなと思っています。

Leave a Reply

Your email address will not be published. Required fields are marked *

You may use these HTML tags and attributes: <a href="" title=""> <abbr title=""> <acronym title=""> <b> <blockquote cite=""> <cite> <code> <del datetime=""> <em> <i> <q cite=""> <s> <strike> <strong>

このサイトはスパムを低減するために Akismet を使っています。コメントデータの処理方法の詳細はこちらをご覧ください