【FreeCAD使い方講座】FreeCADのmatplotlibを使ってGyroid関数を3次元プロットしてみる

※ 当ページには【広告/PR】を含む場合があります。

2024/03/25

【FreeCAD初心者ガイド】FreeCAD付属のmatplotlibを使って3次元モデリングに便利なプロトタイプ関数を思索する

蛸壺の中の工作室｜FreeCADのmatplotlibを使ってGyroid関数を3次元プロットしてみる

以前の記事で、

を簡単に使ってみる方法を紹介したことがあります。

これにより、FreeCADユーザーはわざわざmatlibplot専用のPython環境を整えてあげなくても、FreeCADを立ち上げる感覚そのままにnumpyやscipyから提供される高度な関数をある程度利用することが可能です。

【FreeCAD初心者ガイド】FreeCAD付属のmatplotlibを使って3次元モデリングに便利なプロトタイプ関数を思索する

FreeCADに標準で組込まれているPython製高機能プロッター・matplotlibの基本的な使い方を紹介

今回の記事では、matplotlibを使ったPythonコードの活用の一例として、以下の記事で紹介していた

「Gyroid関数」

の解析解を導いて3次プロットしてみるまでを紹介します。

【pythonで動くG-codeビルダー】G-coordinator/gcoordinatorをLinuxにインストール＆動作確認してみる

Pythonソースコードで3Dプリンター向けのgcodeが生成できる・『G-coordinator』をLinuxへ導入します。

基礎からのFreeCAD 三訂版 (I/O BOOKS)

FreeCAD入門

FreeCADの使い方いろんな物を作りながら3DCADを学ぼう!: 3Dプリンターでの印刷も掲載

図面の描き方がやさしくわかる本

Gyroid関数の解析解

まずGyroid関数は以下のような陰関数で与えられます。

\displaystyle { sin(x) cos(y) + sin(y) cos(z) + sin(z) cos(x) = 0 }

Eq. (gyroid-1)

Gyroidは陰関数で書くとサイクリックな項の和として表せるため、綺麗な見た目をしています。

ただし、このまま３次元グラフを描く際に、Mathematicaのような高度な数学的な処理が可能なソフトウェアが必要になります。

今回はmatplotlibから3次元プロットを試みたいのですが、数値解的な手法を使うと、ソルバをどう設定すれば精度が良いか...など他に考慮しなければならない点が多くなります。

こういったことから、matplotlibで３次元グラフ化するためには、

解析解(陽関数)

を得ることが出来ないかを最初に考えると良いでしょう。

式変形した結果だけを示しますが、Gyroid関数は以下のような2変数関数として解析解を持ちます。

\displaystyle { z = \arcsin \biggl\{ \frac{-\cos{x}\sin{y}\cos{y} \pm \sin{x}\sqrt{\sin^2{x} + \cos^2{y} - \cos^2{x} \sin^2{y}}}{\sin^2{x} + \cos^2{y}} \biggr\} }

Eq. (gyroid-2)

陰関数形とは違い割と綺麗ではないですが、計算としてはより単純です。

注意が必要なのは、分子の中にある平方根の取り扱いで、ルートの中身が負の値をとるとき、

z

は実数解を持ちません。

また、二次方程式由来の

\pm

からも分かる通り、Gyroid関数の解析解は「2枚張り」の対で存在しています。

matplotlibで3次元描画する際にはこの点を注意してスクリプトを書いていきましょう。

基礎からのFreeCAD 三訂版 (I/O BOOKS)

FreeCAD入門

FreeCADの使い方いろんな物を作りながら3DCADを学ぼう!: 3Dプリンターでの印刷も掲載

図面の描き方がやさしくわかる本

matplotlibでの描画

では先程の解析解の考え方を基にして、Pythonコードに落とし込んでいきます。

開く/閉じる

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

#--- メッシュグリッドの生成(1回目) ---

N = 1000
R1 = -np.pi/2
R2 = np.pi/2
x = np.linspace(R1, R2, N)
y = np.linspace(R1, R2, N)
X, Y = np.meshgrid(x, y)

#--- メッシュグリッドの生成(2回目) ---

A1 = np.sin(X)**2 + np.cos(Y)**2
A2 = np.cos(X)*np.sin(Y)
Z0 = A1 - A2**2
X,Y = np.where(Z0 >= 0, X, np.nan), np.where(Z0 >= 0, Y, np.nan)

#--- Gyroid関数の計算 ---

A1 = np.sin(X)**2 + np.cos(Y)**2
A2 = np.cos(X)*np.sin(Y)
A3 = (-A2*np.cos(Y) - np.sin(X)*np.sqrt(A1 - A2**2))/A1
A4 = (-A2*np.cos(Y) + np.sin(X)*np.sqrt(A1 - A2**2))/A1
Z1 = np.arcsin(A3)
Z2 = np.arcsin(A4)

#--- ワイヤーフレーム図として描画 ---

fig = plt.figure()

ax = fig.add_subplot(111, projection='3d')
ax.plot_wireframe(X,Y,Z1)
ax.plot_wireframe(X,Y,Z2)

plt.show()