はじめに

本記事は下記の@Altaka4128氏の三次元有限要素法をSTLファイルを活用して簡単に利用できるようにしたツールの使用方法や実装方法についての解説です。

モチベーションとして、マウスクリックと簡単な数値入力だけで気軽に三次元の有限要素法を試してみたいという思いがありました。
基本は上記参考記事のプログラムをベースにしています。参考記事のものは、メッシュ情報を手書きしていたので、STLファイルから読み込めるように拡張しました。勉強がてら作っているので、誤作動や誤動作を起こす可能性があります。自己責任でご使用ください。

本ツールの使用方法

ツールのダウンロード
下記のGitHubのページからZipでダウンロードして展開するか、git cloneでダウンロードしてください。
https://github.com/hijimasa/python-fem-tool

必要なパッケージのインストール

git clone https://github.com/hijimasa/python-fem-tool.git
cd python-fem-tool
pip install -r requirements.txt

プログラムの起動
```
python main.py
```
以下のようなウインドウが表示されれば成功です。
読み込みファイルの指定
読み込むSTLファイルパスを記載の上、「読み込み」ボタンからファイルの読み込みを実行してください。
読み込むファイルの長さ単位はメートル指定にしたほうが良いと思います。
今回はサンプルとして添付しているtest.stlを読み込んでみます。
固定端の指定
x, y, z方向の力を未入力の状態（あるいは全て0）で表示されたメッシュのノードをクリックすると、ノードが赤色になり固定端になります。右クリックで固定端ではなくなります。
今回は奥の面を構成する４点を固定端として選択しました。固定端に選択したノードが赤くなっているのが解ると思います。
力点の指定
x, y, z方向の力を入力した状態で表示されたメッシュのノードをクリックすると、ノードが緑色になり力の大きさに応じた緑のベクトルが生成されます。右クリックで力点ではなくなります。
今回は手前のノードにz軸下向きに1kNの力を加えてみます。図のように力点に矢印が表示されます。矢印が見切れてしまう場合には右クリック長押しでマウスカーソルを移動させると拡大縮小が出来ます。
各種物性の入力
ヤング率やポアソン比、密度などの項目を適宜設定してください。
また、部材全体にかかる重力の有無を選択できます。
今回はデフォルトのまま実行してみます。
解析開始ボタンを押す
「解析開始」ボタンを押すことで、有限要素法が実行され、結果表示のスケールに合わせた倍率で拡大された変位が表示されます。
デフォルトでは変形の見え方が小さかったので、今回はスケールを大きめで設定しています。
詳細な結果は、pythonを実行した場所に生成されるtest.txtに記載されます。

本ツールの実装

本ツールの作成にあたって重要なポイントを３つ解説します。

STLファイルのメッシュの切り直し

STLファイルからのメッシュ情報の読み出し自体はnumpy-stlライブラリを利用することで簡単に実装できます。しかし、今回参考にしたプログラムはメッシュが四面体により構成されたテトラヘドロンメッシュである必要があります。使用するツールによってSTLに設定されるメッシュの切り方が異なるので、場合によってはテトラヘドロンメッシュにはなっていない場合があります。そのため、STLファイルのメッシュからテトラヘドロンメッシュにメッシュを再構成する必要があります。今回はtetGenというライブラリを使用してテトラヘドロンメッシュへ変換しています。
当該箇所のソースコードは以下のとおりです。

    # STLファイルの読み込み
    stl_mesh = mesh.Mesh.from_file(file_path)

    # 頂点情報の取得
    points = np.unique(stl_mesh.vectors.reshape(-1, 3), axis=0)

    # 面（ファセット）の定義
    faces = []
    for triangle in stl_mesh.vectors:
        face = []
        for vertex in triangle:
            index = np.where((points == vertex).all(axis=1))[0][0]
            face.append(index)
        faces.append(face)

    # TetGenオブジェクトの作成
    tet = tetgen.TetGen(points, np.array(faces))

    # メッシュの生成
    nodes, elems = tet.tetrahedralize(order=1)

ノードの選択機能の実装

今回メッシュを三次元表示するために、matplotlibライブラリを使用しています。
クリックされたノードを特定するためには、ノードを描画するArtist(matplotlibの描画オブジェクト)を保管する必要があります。今回は各ノードをそれぞれ1点だけ描画するscatter（散布図用Artist）で表現しています。クリック時にクリックされたArtistが取得できるので、保管されたArtistと比較することでクリックされたノードを特定できます。
当該箇所のソースコードは以下のとおりです。

# STLメッシュの描画
def drawStl(ax, nodes, elems):
   ...
   # ノードを描画するArtistをnode_scatterリストに保管
    for i in range(len(nodes)):
        node_scatter.append(ax.scatter(nodes[i, 0], nodes[i, 1], nodes[i, 2], color="blue", picker=True))

# ノードがクリックされたときの処理
def on_node_click(event):
    artist = event.artist
    if event.mouseevent.button == 1:  # 左クリック
        ...
        for i in range(len(nodes)):
            if node_scatter[i] == artist:
                # 設定された力の値に応じて固定端または力点を設定
                ...
fig = Figure(figsize=(6, 6))
ax = fig.add_subplot(111, projection="3d")
...
# イベントリスナーを登録
fig.canvas.mpl_connect("pick_event", on_node_click)

変形後メッシュの描画

各ノードの変形量が解析結果から得られるので、それをスケールに設定した値を掛けて各ノードの位置に加算しています。単純な処理ですが、実際に表示してみるとそれらしい形状を得ることが出来たので良かったです。
当該箇所のソースコードは以下のとおりです。

    # 有限要素法を実行
    fem = FEM(fem_nodes, fem_elems, bound)
    fem.analysis()
    fem.outputTxt("test")
    # 各ノードの変形量を取得
    displacement = fem.outputDisplacement()
    ...
    # スケール値を取得
    try:
        scale = float(entry_scale.get())
    except ValueError:
        print("wrong scale is setted")
        scale = 10000.0

    # 変形したノードを生成
    result_nodes = np.zeros((len(nodes), 3), dtype=np.float32)
    for i in range(len(nodes)):
        for j in range(3):
            result_nodes[i][j] = nodes[i][j] + scale * displacement[i][j]

    #result_nodesを用いて四面体ごとにメッシュを描画
    for elem in elems:
        tetra = result_nodes[elem]
        verts = [
            [tetra[0], tetra[1], tetra[2]],
            [tetra[0], tetra[1], tetra[3]],
            [tetra[0], tetra[2], tetra[3]],
            [tetra[1], tetra[2], tetra[3]]
        ]
        # 再描画時に消せるようにリストに保管
        draw_result.append(ax.add_collection3d(Poly3DCollection(verts, edgecolor="b", alpha=0.0)))
    canvas.draw()

おわりに

本記事では、STLファイルを入力とする３次元有限要素法ツールの使用方法と実装方法について解説しました。狙い通りマウスクリックと簡単な数値入力で有限要素法を試せるツールが出来たので、個人的には満足しています。
当初はこんなツールが作れるライブラリなんてあるのかと疑問でしたが、探せば色々見つかるものですね。今後もこんな形であったら良いなというツールを作っていけたら良いなと思っています。