はじめに
久しぶりの投稿である。
まずは下の図を見ていただきたい。
これは、節点数 10000、要素数 19410 の三角形要素を用いた2次元FEM計算結果の出力である。FEMの計算は3秒そこそこで完了するのであるが、この図1枚を描くのに55秒かかっている。
これをなんとか高速化したいと昔から考えていたのだが、今回ある程度早くできる方法が見つかったので、これを紹介したい。
使っている人はすでに使っていると思うが、自分にとっては発見であったので、記録として残しておきたいと思った次第である。
紹介する方法は、matplotlib polygon
で検索をかけていて見つけた以下のページの内容を応用したののである。
これまでの方法
まずは、これまで使っていた三角形の塗りつぶし方法を紹介する。
下のような図を出力するコードを考える。
import matplotlib.pyplot as plt
import numpy as np
import time
def main():
start = time.time()
npoin,nele,nod=12,12,3
x=np.array([1,2,3,4,1,2,3,4,1,2,3,4])
y=np.array([1,1,1,1,2,2,2,2,3,3,3,3])
node=np.zeros((nod+1,nele),dtype=np.int32)
node[:,0]=np.array([1,2,5,1])
node[:,1]=np.array([5,2,6,2])
node[:,2]=np.array([2,3,6,3])
node[:,3]=np.array([6,3,7,1])
node[:,4]=np.array([3,4,7,2])
node[:,5]=np.array([7,4,8,3])
node[:,6]=np.array([5,6,9,1])
node[:,7]=np.array([9,6,10,2])
node[:,8]=np.array([6,7,10,3])
node[:,9]=np.array([10,7,11,1])
node[:,10]=np.array([7,8,11,2])
node[:,11]=np.array([11,8,12,3])
fsz=10
xmin,xmax=0,5
ymin,ymax=0,4
plt.figure(figsize=(5,4),facecolor='w')
plt.xlim([xmin,xmax])
plt.ylim([ymin,ymax])
col=['#ffff00', '#00ffff', '#ff0000']
for ne in range(0,nele):
icol=node[-1,ne]-1
nn=node[0:nod,ne]-1
xa=x[nn]
ya=y[nn]
plt.fill(xa,ya,facecolor=col[icol],edgecolor='#aaaaaa',lw=0.5)
ds=0.1
for i in range(0,npoin):
xs,ys,ss=x[i]-ds,y[i]-ds,str(i+1)
plt.text(xs,ys,ss,va='center',ha='center',fontsize=fsz)
fnameF='fig_test1.jpg'
plt.savefig(fnameF, dpi=100, bbox_inches="tight", pad_inches=0.1)
#plt.show()
print(time.time() - start)
#---------------
# Execute
#---------------
if __name__ == '__main__': main()
上記コードに示す、これまでの方法の概要は以下の通り。
点の座標、三角形の定義は、少しFEMチックにしている。
- 座標点数 npoin=12、これらの座標点で構成される三角形の数をnele=12とする。
- 配列 xおよびyで、(x,y)座標を定義する。
- 配列nodeに、三角形を構成する店の番号および色を指定する番号を格納する。例えば、
node[:,0]=np.array([1,2,5,1])
という場合、はじめの3つの数値 1,2,3が三角形を構成する点の番号、4つめの1は色を指定する番号である。 - 使う色は、colというリストに
col=['#ffff00', '#00ffff', '#ff0000']
という形で格納する。 - 三角形の塗りつぶしは、
for
文を回しながら、三角形を構成する点 (xa, ya)を定義して、plt.fill
により一つずつ塗りつぶしていく。これまでの方法はこの処理に相当な時間を費やしていた。
時短のための方法
時短のためのコードは以下に示す通り。
import matplotlib.pyplot as plt
import numpy as np
import time
from matplotlib.patches import Polygon
from matplotlib.collections import PatchCollection
def main():
start = time.time()
npoin,nele,nod=12,12,3
x=np.array([1,2,3,4,1,2,3,4,1,2,3,4])
y=np.array([1,1,1,1,2,2,2,2,3,3,3,3])
node=np.zeros((nod+1,nele),dtype=np.int32)
node[:,0]=np.array([1,2,5,1])
node[:,1]=np.array([5,2,6,2])
node[:,2]=np.array([2,3,6,3])
node[:,3]=np.array([6,3,7,1])
node[:,4]=np.array([3,4,7,2])
node[:,5]=np.array([7,4,8,3])
node[:,6]=np.array([5,6,9,1])
node[:,7]=np.array([9,6,10,2])
node[:,8]=np.array([6,7,10,3])
node[:,9]=np.array([10,7,11,1])
node[:,10]=np.array([7,8,11,2])
node[:,11]=np.array([11,8,12,3])
fsz=10
xmin,xmax=0,5
ymin,ymax=0,4
plt.figure(figsize=(5,4),facecolor='w')
plt.xlim([xmin,xmax])
plt.ylim([ymin,ymax])
patches=[]
lcol=[]
col=['#ffff00', '#00ffff', '#ff0000']
for ne in range(0,nele):
icol=node[-1,ne]-1
nn=node[0:nod,ne]-1
xa=x[nn]
ya=y[nn]
points=[[xa[0],ya[0]],[xa[1],ya[1]],[xa[2],ya[2]],[xa[0],ya[0]]]
polygon=Polygon(xy=points,closed=True)
patches.append(polygon)
lcol.append(col[icol])
p=PatchCollection(patches,facecolor=lcol,edgecolor='#aaaaaa',lw=0.5)
plt.gca().add_collection(p)
ds=0.1
for i in range(0,npoin):
xs,ys,ss=x[i]-ds,y[i]-ds,str(i+1)
plt.text(xs,ys,ss,va='center',ha='center',fontsize=fsz)
fnameF='fig_test2.jpg'
plt.savefig(fnameF, dpi=100, bbox_inches="tight", pad_inches=0.1)
#plt.show()
print(time.time() - start)
#---------------
# Execute
#---------------
if __name__ == '__main__': main()
時短案の概要は以下の通り。
- まず、matplotlib.pyplotとnumpyに加え、モジュールを追加する。
from matplotlib.patches import Polygon
from matplotlib.collections import PatchCollection
- 三角形を構成する座標を格納するリスト
patches
と色情報を格納するリストlcol
を準備する。 -
for
文でPolygon
を作成してリストpatches
に格納していくと同時に、色情報をリストlcol
に格納していく。 -
Polygon
は閉じるように定義しないと線をうまく引いてくれないので注意。 -
PatchCollection
の引数に、plt.fill
で指定したのと同じ要領で、指定したのと同様に、三角形の座標情報patches
、塗りつぶし色facecolor=lcol
、境界職edgecolor='@aaaaaa'
、線の太さlw=0.5
を指定する。 - あとは
add_collection
で描画する。 -
for
文の中で作図はせず、PatchCollection
に保存してあとで一気に作図するのが決め手のようである。
効果
上記2種類のコードの実行結果は、
これまでの方法・時短案ともに0.1秒程度で差はない。
しかし、最初に示した約2万個の三角形を塗りつぶす作業を複数回行う図となると、時短案採用により、作図時間は55秒からでは15秒まで短縮できており、作業時間を27%に短縮できる。
これで、作図による待ち時間を大幅に短縮できたが、計算時間3秒に対して作図時間15秒は直感的には長い気がする。機会をみてさらなる時短を考えていきたい。
最後になるが、計算環境に触れておく。
- M1 Pro MacBook Pro 14''
- Memory 16GB
- macOS Ventura 13.4.1 (c)
- Python 3.11.2 (環境はpyenvで構築)
(2023.08.05 追記)
最初に紹介した図は、ベクトル図も含まれているが、ベクトルのfor文で1要素1要素書くのではなく、LineCollectionを使うことにより、1枚あたりの作図時間は7秒にまで短縮できた。
LineCollectionの活用は、以下の記事を参考にした。
https://salad-bowl-of-knowledge.github.io/hp/python/2019/10/04/multilines.html
以 上