はじめに
有限要素法というものがあります。
企業の設計の紹介とかで見るカラフルな色がついてどこに力がかかっているかわかるやつです。
それを見るツールはあるのですが、これもしかしたらBlenderでも見れるのでは,
と思い、プログラムを作りました。
幸いBlenderはPythonでの操作に対応しており、以下のようにできました。
今回はそのプログラムについて解説していきたいと思います。
入力ファイル群とプログラムは以下に置いておきます。
Blenderの準備
Blenderは無料で使えオープンソースで開発されている高機能3Dモデリングソフトです。
使い方に関しては、解説書を公開している方がいましたので、そちらで勉強しました。
私はダウンロードのしやすさからSteamからBlenderをダウンロードしています。
ライブラリのインストール
Blenderには独自のPythonが入っており、普段使用しているAnaconda環境度は別に、Blender内で使用するライブラリをインストールする必要があります。
初期のBlenderにはライブラリがあまりインストールされていません。
そこで以下のサイトを参考に今回必要なmatplotlibをインストールします。
今回、私はGドライブ直下にSteamフォルダがあり、そこにソフトをダウンロードしています。
上記サイトを参考にすると、以下の位置にpython.exeが見つかりました。
（3.1の部分はBlenderのバージョンなので読み替えてください。）
G:\Steam\steamapps\common\Blender\3.1\python\bin\python.exe
そして順番に、コマンドプロンプトで以下を実行します。
上から、
- pipの有効化
- pipのアップグレード
- インストール済みライブラリ一覧の表示
となっています。
G:\Steam\steamapps\common\Blender\3.2\python\bin\python.exe -m ensurepip
G:\Steam\steamapps\common\Blender\3.2\python\bin\python.exe -m pip install --user --upgrade pip
G:\Steam\steamapps\common\Blender\3.2\python\bin\python.exe -m pip list
pip listで調べると、以下のように最低限しか入っていないことがわかります。
Package Version
------------------ ---------
autopep8 1.6.0
certifi 2021.10.8
charset-normalizer 2.0.10
Cython 0.29.26
idna 3.3
numpy 1.22.0
pip 22.1.2
pycodestyle 2.8.0
requests 2.27.1
setuptools 58.1.0
toml 0.10.2
urllib3 1.26.8
zstandard 0.16.0
そして、今回のプログラムで必要なmatplotlibをダウンロードします。
G:\Steam\steamapps\common\Blender\3.2\python\bin\python.exe -m pip install --user -U matplotlib
バージョンを変えてしまった時の注意
SteamでBlenderをアップデートすると、インストールしたライブラリが使えなくなってしまいます。
あたらしいバージョンのほうでライブラリを再インストールしようとして、なぜか古いバージョンのフォルダにダウンロードされ、ライブラリが新しいほうに入らなくなるので注意です。
解決法としては、無理やりですが、
バージョンが3.1から3.2になったときは、
G:\Steam\steamapps\common\Blender\3.1\python\lib\site-packages
の中身を
G:\Steam\steamapps\common\Blender\3.2\python\lib\site-packages
にコピーすることで解決します。
BlenderのPythonでの操作
Blender Pythonと検索すると色々やり方は出てきますが、一応説明します。
まずBlenderを開いたら、上のほうにあるタブから、Scriptingを選びます。
出てきた画面で左下にあるのがIPythonシェルで、対話モードでPythonを実行することにできます。
右上のところからは、普段の.pyファイルから実行することができます。
今回は、右上の部分から新規作成し、プログラムをコピペして三角矢印を押して実行します。
画像で出てくるタブが見つからないときは、隠れてしまっている可能性があるので、タブをマウスのホイールでクリックしながら左に動かして出しましょう。
データの用意
今回は、3次元8節点六面体要素 の可視化をしていきたいと思います。
2次元や3次元四面体要素、2次要素などいろいろありますが、今回紹介する方法を使用すればご自分でも実装できるかもしれません。
https://monoist.itmedia.co.jp/mn/articles/0909/02/news100_2.html
今回必要なファイルは、
- 節点番号、座標ファイル
- 要素番号、コネクティビティファイル
- 要素値ファイル
の3つです。
順番に解説していきます。
節点番号、座標ファイル
文字通り、節点の番号、座標（xyz）を記入したファイルです。
以下のように、各行に順番に節点の情報を書いていきます。
| 節点番号 | x座標 | y座標 | z座標 |
|---|---|---|---|
| 1 | 0 | 0 | 0 |
| 2 | 0.1 | 4 | 0 |
| 3 | 2 | 0. 05 | 2 |
ファイル形式は、Pythonで読み込んで配列に変換出来たらなんでもいいですが、
今回はスペース区切りのファイルにします。
要素番号、コネクティビティファイル
各要素を構成する節点番号をまとめたファイルです。
わかりやすく説明するため、2次元3節点要素で解説します。
上記のようなメッシュを切った場合、
要素0は0,1,4
要素1は1,2,3
要素2は1,3,4
要素3は0,4,5
という節点番号から構成されています。
これをまとめます。
今回の3次元8節点六面体要素 では、要素の形は6面体で、節点8個から成り立っています。
そしてその要素を構成する節点をまとめたファイルです。
以下のように、要素番号、構成節点番号を順番に書いたものを1行として書いていっています。
| 要素番号 | 構成節点番号1 | 構成節点番号2 | 構成節点番号3 | 構成節点番号4 | 構成節点番号5 | 構成節点番号6 | 構成節点番号7 | 構成節点番号8 |
|---|---|---|---|---|---|---|---|---|
| 1 | 1 | 4 | 5 | 2 | 8 | 11 | 25 | 78 |
| 2 | 2 | 8 | 1 | 3 | 72 | 52 | 44 | 25 |
| 3 | 10 | 5 | 7 | 50 | 43 | 20 | 1 | 4 |
ファイル形式は、Pythonで読み込んで配列に変換出来たらなんでもいいですが、
今回はカンマ区切りのファイルにします。
要素値ファイル
メッシュのみを可視化するなら、上記2つの節点、要素だけでできます。
しかし、有限要素法では、よく要素における値（応力やひずみなど）を持っています。
今回は、この要素値によって色を付けたいと思います。
節点値を持っていることもありますが、今回は対応していないです。
ファイルは以下のように、要素番号、その要素における値、をまとめたものです。
| 要素番号 | 値 |
|---|---|
| 1 | 2.35 |
| 2 | 5.2 |
| 3 | 2.2 |
ファイル形式は、Pythonで読み込んで配列に変換出来たらなんでもいいですが、
今回はスペース区切りのファイルにします。
Blenderでの表示プログラム
いよいよ、blenderを使用して可視化をしていきます。
まず完成したプログラムが以下で、ここについて解説していきます。
全プログラム
import numpy as np
import matplotlib.pyplot as plt
import csv
#条件 ファイルから読み込んでもよい
#節点数
num_node = 1331
#要素数
num_eleme = 1000
#表示倍率
amp = 10
node = np.empty((num_node,3), dtype=np.float64) #節点座標
eleme = np.empty((num_eleme,8),dtype=np.int32) #各要素のコネクティビティ #つまりある六面体elementを構成する接点node番号(1スタートに注意)
eleme_value = np.empty((num_node), dtype=np.float64) #要素における値 これで色を付ける
#節点番号、座標ファイルの読み込み
with open('H:/programing/blender/output_disp2.dat') as f:
reader = csv.reader(f, delimiter=' ')
l = [row for row in reader]
for i in range(num_node):
node[i,0] = l[i][1].replace('d','e')
node[i,1] = l[i][2].replace('d','e')
node[i,2] = l[i][3].replace('d','e')
#要素番号、コネクティビティファイルの読み込み
with open('H:/programing/blender/input_eleme.txt') as f:
reader = csv.reader(f)
l = [row for row in reader]
for i in range(num_eleme):
eleme[i] = l[i][1:9]
#要素値ファイルの読み込み
with open('H:/programing/blender/output_ave_strain2.dat') as f:
reader = csv.reader(f, delimiter=' ')
#reader = csv.reader(f)
l = [row for row in reader]
for i in range(num_eleme):
eleme_value[i] = l[i][1].replace('d','e')
#正規化する関数
#https://deepage.net/features/numpy-normalize.html
def min_max(x, axis=None):
min = x.min(axis=axis, keepdims=True)
max = x.max(axis=axis, keepdims=True)
result = (x-min)/(max-min)
return result
#0から1に正規化した要素値を0から255に変更
eleme_value = min_max(eleme_value, axis=None) * 255
#要素の値からmatplotlibのカラーマップ機能を使用して(R, G, B, A)に変換
cm_color = [plt.cm.jet(int(value)) for value in eleme_value]
import bpy
for item in bpy.data.meshes:
bpy.data.meshes.remove(item)
for item in bpy.data.materials:
bpy.data.materials.remove(item)
verts = np.zeros((8,3), dtype=np.float64) #ある六面体要素を構成する6接点のxyz座標
for i in range(num_eleme):
for j in range(8): #節点の8
#eleme[i,j]つまり要素値が1始まりなので0始まりのPython配列に合わせて-1している
verts[j,0] = node[eleme[i,j]-1,0]
verts[j,1] = node[eleme[i,j]-1,1]
verts[j,2] = node[eleme[i,j]-1,2]
#表示倍率を適用
verts *= amp
#面を構成する節点の場所
faces = [[0,1,2,3], [0,4,5,1], [1,5,6,2], [2,6,7,3], [0,3,7,4], [4,7,6,5]]
#メッシュの作成
msh = bpy.data.meshes.new(name=f"cubemesh{i}")
msh.from_pydata(verts, [], faces)
msh.update()
#オブジェクトの作成
obj = bpy.data.objects.new(name=f"cube{i}", object_data=msh)
scene = bpy.context.scene
scene.collection.objects.link(obj)
#マテリアルの作成
material = bpy.data.materials.new(f"cubematerial{i}")
obj.data.materials.append(material)
material.use_nodes = True
#マテリアルの設定
color = cm_color[i]
p_BSDF = material.node_tree.nodes["Principled BSDF"]
p_BSDF.inputs[0].default_value = color
p_BSDF.inputs[9].default_value = 0
p_BSDF.inputs[17].default_value = 1
#環境光で全体を明るく
#環境色を(R, G, B, A)で指定、今回は白
bpy.data.worlds["World"].node_tree.nodes["Background"].inputs[0].default_value = (1, 1, 1, 1)
#強さ
bpy.data.worlds["World"].node_tree.nodes["Background"].inputs[1].default_value = 1
条件
プログラムを動かすのにあたリ、節点数、要素数、表示倍率
が必要です。
これがなくても動くようなプログラムを作ることもできそうですが、有限要素法を使っているならあるはずなので、今回は使用します。
外部から読み込むプログラムを作ることも可能です。
#条件 ファイルから読み込んでもよい
#節点数
num_node = 1331
#要素数
num_eleme = 1000
#表示倍率
amp = 10
データの読み込み
先ほど用意した3ファイルを読み込んでいきます。
blenderのカレントディレクトリは<built-in function dir>と出て場所がわからないので、ファイルの場所を絶対パスでスラッシュを使ってパスを書きます。
今回は、カンマ区切り/スペース区切りのファイルですので、Python標準ライブラリのcsvライブラリを使用します。
読み込んだものは、numpyライブラリのndarrayオブジェクトに入れます。
ちなみにこれはFEMプログラムの名残だったりします。
また、節点座標に関しては、Fortranの倍精度浮動小数点を表した際のdを読み込めるようにeに変えていますが、元からeでも、ただの小数でも、整数でも読み込めます。
また、区切り文字に関しては、カンマ区切りだけでなく、スペース区切りや任意の区切り文字に対応しています。
変更したいときは、reader = csv.reader(f, delimiter=' ')というようにdelimiter変数で変更します。
今回は試しとして、節点ファイル、要素値ファイルはスペース区切り、要素番号コネクティビティファイルはカンマ区切りで用紙して読み込んでみましょう。
node = np.empty((num_node,3), dtype=np.float64) #節点座標
eleme = np.empty((num_eleme,8),dtype=np.int32) #各要素のコネクティビティ #つまりある六面体elementを構成する接点node番号(1スタートに注意)
eleme_value = np.empty((num_node), dtype=np.float64) #要素における値 これで色を付ける
#節点番号、座標ファイルの読み込み
with open('H:/programing/blender/output_disp2.dat') as f:
#reader = csv.reader(f) #カンマ区切り用
reader = csv.reader(f, delimiter=' ')
l = [row for row in reader]
for i in range(num_node):
node[i,0] = l[i][1].replace('d','e')
node[i,1] = l[i][2].replace('d','e')
node[i,2] = l[i][3].replace('d','e')
#要素番号、コネクティビティファイルの読み込み
with open('H:/programing/blender/input_eleme.txt') as f:
reader = csv.reader(f)
l = [row for row in reader]
for i in range(num_eleme):
eleme[i] = l[i][1:9]
#要素値ファイルの読み込み
with open('H:/programing/blender/output_ave_strain2.dat') as f:
reader = csv.reader(f, delimiter=' ')
#reader = csv.reader(f)
l = [row for row in reader]
for i in range(num_eleme):
eleme_value[i] = l[i][1].replace('d','e')
#正規化する関数
#https://deepage.net/features/numpy-normalize.html
def min_max(x, axis=None):
min = x.min(axis=axis, keepdims=True)
max = x.max(axis=axis, keepdims=True)
result = (x-min)/(max-min)
return result
#0から1に正規化した要素値を0から255に変更
eleme_value = min_max(eleme_value, axis=None) * 255
#要素の値からmatplotlibのカラーバー機能を使用して色に変換
cm_color = [plt.cm.jet(int(value)) for value in eleme_value]
色の作成
要素の値からグラデーションの色を作成します。
今回はmatplotlibのカラーマップを使用して要素の値から(R, G, B, A)に変換します。
カラーマップでは0から255の整数に応じて色が返ってくるため、0から1に正規化したのちに0から255の整数に変換します。
#正規化する関数
def min_max(x, axis=None):
min = x.min(axis=axis, keepdims=True)
max = x.max(axis=axis, keepdims=True)
result = (x-min)/(max-min)
return result
#0から1に正規化した要素値を0から255に変更
eleme_value = min_max(eleme_value, axis=None) * 255
#要素の値からmatplotlibのカラーマップ機能を使用して色に変換
cm_color = [plt.cm.jet(int(value)) for value in eleme_value]
もとからあるオブジェクトを削除
Blenderは最初には立方体が置いてあります。
なのでこれを削除しましょう。
for item in bpy.data.meshes:
bpy.data.meshes.remove(item)
for item in bpy.data.materials:
bpy.data.materials.remove(item)
メッシュ・オブジェクト・マテリアルの作成
要素ごとに順番にメッシュ・オブジェクト・マテリアルを作成していきます。
先にプログラムを出すと、このようになっています。
verts = np.zeros((8,3), dtype=np.float64) #ある六面体要素を構成する6接点のxyz座標
for i in range(num_eleme):
for j in range(8): #節点の8
#eleme[i,j]つまり要素値が1始まりなので0始まりのPython配列に合わせて-1している
verts[j,0] = node[eleme[i,j]-1,0]
verts[j,1] = node[eleme[i,j]-1,1]
verts[j,2] = node[eleme[i,j]-1,2]
#表示倍率を適用
verts *= amp
#面を構成する節点の場所
faces = [[0,1,2,3], [0,4,5,1], [1,5,6,2], [2,6,7,3], [0,3,7,4], [4,7,6,5]]
#メッシュの作成
msh = bpy.data.meshes.new(name=f"cubemesh{i}")
msh.from_pydata(verts, [], faces)
msh.update()
#オブジェクトの作成
obj = bpy.data.objects.new(name=f"cube{i}", object_data=msh)
scene = bpy.context.scene
scene.collection.objects.link(obj)
#マテリアルの作成
material = bpy.data.materials.new(f"cubematerial{i}")
obj.data.materials.append(material)
material.use_nodes = True
#マテリアルの設定
color = cm_color[i]
p_BSDF = material.node_tree.nodes["Principled BSDF"]
p_BSDF.inputs[0].default_value = color
p_BSDF.inputs[9].default_value = 0
p_BSDF.inputs[17].default_value = 1
メッシュの作成
まず6面体のメッシュを作成しましょう。
メッシュを作成するためには、メッシュを構成する頂点の座標vertsと、面を構成する頂点をまとめたfacesが必要です。6面体の節点の順番によっては、立方体が構成されなかったり、面が裏返っていることがあります。その時はうまくやってください。（丸投げですみません。）
オブジェクトの作成
オブジェクトを作成します。
オブジェクトの中にメッシュ、マテリアルなどの情報を含むようになります。
マテリアルの作成
マテリアルを作成し、色、透明度、粗さなどを設定していきます。
p_BSDF.inputs[i]というのは、画像のようにBlender内のShadingにあるプリシンプルBSDFの左に丸がある部分の上からの0からの番号です。
この番号ですが、どうもBlenderのバージョンによって変わってしまうので、ほかのサイトの記事と違っていたりすることもあるので注意してください。
すなわち、今回は
p_BSDF.inputs[0]で色を決め、
p_BSDF.inputs[9]で粗さを0、つまりつるつるにし、
p_BSDF.inputs[17]で伝搬を1に。つまり透明にしています。
ライティング
光源については触らず、環境光で全体を明るくします。
#環境光で全体を明るく
#環境色を(R, G, B, A)で指定、今回は白
bpy.data.worlds["World"].node_tree.nodes["Background"].inputs[0].default_value = (1, 1, 1, 1)
#強さ
bpy.data.worlds["World"].node_tree.nodes["Background"].inputs[1].default_value = 1
観察
ちゃんと色がついているか観察してみましょう。
Shadingから見てもいいですが、今回は、Layoutタブ内で見ていきます。
右上の丸いマークのうち、右から2番目を押します。
さらに右の下矢印をクリックすれば環境を変えることもできます。
結果がこちらです。
ちゃんと色がついていますね
まとめ
以上のプログラムで、Blender内で有限要素法の結果を表示することができるようになりました。Blender内で観察するもよし、エクスポートするもよしです。
BlenderをPythonで操作するというのがそもそもニッチな世界で、解説記事を書いてくださった先人には感謝しかありません。
今後としては、アニメーションもできるようにしたり、任意の位置でスライスして中身を見るようにできるようにしてみたいです。
参考
有限要素法関連の記事